MOAB: Mesh Oriented datABase  (version 5.2.1)
MOAB_iMeshP_unit_tests.cpp
Go to the documentation of this file.
00001 #include "iMeshP.h"
00002 #include "moab_mpi.h"
00003 #include <iostream>
00004 #include <algorithm>
00005 #include <vector>
00006 #include <sstream>
00007 #include <assert.h>
00008 #include <math.h>
00009 #include <map>
00010 #include <string.h>
00011 #include <stdio.h>  // remove()
00012 
00013 #if !defined( _MSC_VER ) && !defined( __MINGW32__ )
00014 #include <unistd.h>
00015 #endif
00016 
00017 #define STRINGIFY_( X ) #X
00018 #define STRINGIFY( X )  STRINGIFY_( X )
00019 const char* const FILENAME = "iMeshP_test_file";
00020 
00021 /**************************************************************************
00022                               Error Checking
00023  **************************************************************************/
00024 
00025 #define CHKERR                                                                                         \
00026     do                                                                                                 \
00027     {                                                                                                  \
00028         if( ierr )                                                                                     \
00029         {                                                                                              \
00030             std::cerr << "Error code  " << ierr << " at " << __FILE__ << ":" << __LINE__ << std::endl; \
00031             return ierr;                                                                               \
00032         }                                                                                              \
00033     } while( false )
00034 
00035 #define PCHECK                            \
00036     do                                    \
00037     {                                     \
00038         ierr = is_any_proc_error( ierr ); \
00039         CHKERR;                           \
00040     } while( false )
00041 
00042 // Use my_rank instead of rank to avoid shadowed declaration
00043 #define ASSERT( A )                                                                                                 \
00044     do                                                                                                              \
00045     {                                                                                                               \
00046         if( is_any_proc_error( !( A ) ) )                                                                           \
00047         {                                                                                                           \
00048             int my_rank = 0;                                                                                        \
00049             MPI_Comm_rank( MPI_COMM_WORLD, &my_rank );                                                              \
00050             if( 0 == my_rank )                                                                                      \
00051                 std::cerr << "Failed assertion: " #A << std::endl << "  at " __FILE__ ":" << __LINE__ << std::endl; \
00052             return 1;                                                                                               \
00053         }                                                                                                           \
00054     } while( false )
00055 
00056 // Test if is_my_error is non-zero on any processor in MPI_COMM_WORLD
00057 int is_any_proc_error( int is_my_error )
00058 {
00059     int result;
00060     int err = MPI_Allreduce( &is_my_error, &result, 1, MPI_INT, MPI_MAX, MPI_COMM_WORLD );
00061     return err || result;
00062 }
00063 
00064 /**************************************************************************
00065                            Test  Declarations
00066  **************************************************************************/
00067 
00068 class PartMap;
00069 
00070 /**\brief Consistency check for parallel load
00071  *
00072  * All other tests depend on this one.
00073  */
00074 int test_load( iMesh_Instance, iMeshP_PartitionHandle prtn, PartMap& map, int comm_size );
00075 
00076 /**\brief Test partition query methods
00077  *
00078  * Test:
00079  * - iMeshP_getPartitionComm
00080  * - iMeshP_getNumPartitions
00081  * - iMeshP_getPartitions
00082  */
00083 int test_get_partitions( iMesh_Instance, iMeshP_PartitionHandle prtn, const PartMap& );
00084 
00085 /**\brief Test part quyery methods
00086  *
00087  * Test:
00088  * - iMeshP_getNumGlobalParts
00089  * - iMeshP_getNumLocalParts
00090  * - iMeshP_getLocalParts
00091  */
00092 int test_get_parts( iMesh_Instance, iMeshP_PartitionHandle prtn, const PartMap& );
00093 
00094 /**\brief Test query by entity type
00095  *
00096  * Test:
00097  * - iMeshP_getNumOfTypeAll
00098  * - iMeshP_getNumOfType
00099  * - iMeshP_getEntities
00100  * -
00101  */
00102 int test_get_by_type( iMesh_Instance, iMeshP_PartitionHandle prtn, const PartMap& );
00103 
00104 /**\brief Test query by entity topology
00105  *
00106  * Test:
00107  * - iMeshP_getNumOfTopoAll
00108  * - iMeshP_getNumOfTopo
00109  * - iMeshP_getEntities
00110  * -
00111  */
00112 int test_get_by_topo( iMesh_Instance, iMeshP_PartitionHandle prtn, const PartMap& );
00113 
00114 /**\brief Test mapping from part id to part handle
00115  *
00116  * Test:
00117  * - iMeshP_getPartIdFromPartHandle
00118  * - iMeshP_getPartIdsFromPartHandlesArr
00119  * - iMeshP_getPartHandleFromPartId
00120  * - iMeshP_getPartHandlesFromPartsIdsArr
00121  * - iMeshP_getRankOfPart
00122  * - iMeshP_getRankOfPartArr
00123  */
00124 int test_part_id_handle( iMesh_Instance, iMeshP_PartitionHandle prtn, const PartMap& );
00125 
00126 /**\brief Test get part rank
00127  *
00128  * Tests:
00129  * - iMeshP_getRankOfPart
00130  * - iMeshP_getRankOfPartArr
00131  */
00132 int test_part_rank( iMesh_Instance, iMeshP_PartitionHandle prtn, const PartMap& );
00133 
00134 /**\brief Test querying of part neighbors
00135  *
00136  * Test:
00137  * - iMeshP_getNumPartNbors
00138  * - iMeshP_getNumPartNborsArr
00139  * - iMeshP_getPartNbors
00140  * - iMeshP_getPartNborsArr
00141  */
00142 int test_get_neighbors( iMesh_Instance, iMeshP_PartitionHandle prtn, const PartMap& );
00143 
00144 /**\brief Test querying of part boundary entities
00145  *
00146  * Test:
00147  * - iMeshP_getNumPartBdryEnts
00148  * - iMeshP_getPartBdryEnts
00149  */
00150 int test_get_part_boundary( iMesh_Instance, iMeshP_PartitionHandle prtn, const PartMap& );
00151 
00152 /**\brief Test querying of part boundary entities
00153  *
00154  * Test:
00155  * - iMeshP_initPartBdryEntIter
00156  * - iMeshP_initPartBdryEntArrIter
00157  */
00158 int test_part_boundary_iter( iMesh_Instance, iMeshP_PartitionHandle prtn, const PartMap& );
00159 
00160 /**\brief Test adjacent entity query
00161  *
00162  * Test:
00163  * - iMeshP_getAdjEntities
00164  */
00165 int test_get_adjacencies( iMesh_Instance, iMeshP_PartitionHandle prtn, const PartMap& );
00166 
00167 /**\brief Test entity iterators
00168  *
00169  * Test:
00170  * - iMeshP_initEntIter
00171  * - iMeshP_initEntArrIter
00172  */
00173 int test_entity_iterator( iMesh_Instance, iMeshP_PartitionHandle prtn, const PartMap& );
00174 
00175 /**\brief Test entity owner queries
00176  *
00177  * Test:
00178  * - iMeshP_getEntOwnerPart
00179  * - iMeshP_getEntOwnerPartArr
00180  * - iMeshP_isEntOwner
00181  * - iMeshP_isEntOwnerArr
00182  */
00183 int test_entity_owner( iMesh_Instance, iMeshP_PartitionHandle prtn, const PartMap& );
00184 
00185 /**\brief Test entity status
00186  *
00187  * Test:
00188  * - iMeshP_getEntStatus
00189  * - iMeshP_getEntStatusArr
00190  */
00191 int test_entity_status( iMesh_Instance, iMeshP_PartitionHandle prtn, const PartMap& );
00192 
00193 /**\brief Test information about entity copies for interface entities
00194  *
00195  * Test:
00196  * - iMeshP_getNumCopies
00197  * - iMeshP_getCopyParts
00198  */
00199 int test_entity_copy_parts( iMesh_Instance, iMeshP_PartitionHandle prtn, const PartMap& );
00200 
00201 /**\brief Test information about entity copies for interface entities
00202  *
00203  * Test:
00204  * - iMeshP_getCopies
00205  * - iMeshP_getCopyOnPart
00206  * - iMeshP_getOwnerCopy
00207  */
00208 int test_entity_copies( iMesh_Instance, iMeshP_PartitionHandle prtn, const PartMap& );
00209 
00210 /**\brief Test creation of ghost entities
00211  *
00212  * Test:
00213  * - iMeshP_createGhostEntsAll
00214  */
00215 int test_create_ghost_ents( iMesh_Instance, iMeshP_PartitionHandle prtn, const PartMap& );
00216 
00217 /**\brief Test commuinication of tag data
00218  *
00219  * Test:
00220  * - iMeshP_pushTags
00221  * - iMeshP_pushTagsEnt
00222  */
00223 int test_push_tag_data_iface( iMesh_Instance, iMeshP_PartitionHandle prtn, const PartMap& );
00224 int test_push_tag_data_ghost( iMesh_Instance, iMeshP_PartitionHandle prtn, const PartMap& );
00225 
00226 int test_exchange_ents( iMesh_Instance imesh, iMeshP_PartitionHandle prtn, const PartMap& map );
00227 
00228 /**************************************************************************
00229                               Helper Funcions
00230  **************************************************************************/
00231 
00232 class PartMap
00233 {
00234   public:
00235     int num_parts() const
00236     {
00237         return sortedPartList.size();
00238     }
00239 
00240     iMeshP_Part part_id_from_local_id( int local_id ) const
00241     {
00242         return sortedPartList[idx_from_local_id( local_id )];
00243     }
00244 
00245     int local_id_from_part_id( iMeshP_Part part ) const
00246     {
00247         return partLocalIds[idx_from_part_id( part )];
00248     }
00249 
00250     int rank_from_part_id( iMeshP_Part part ) const
00251     {
00252         return partRanks[idx_from_part_id( part )];
00253     }
00254 
00255     int rank_from_local_id( int id ) const
00256     {
00257         return partRanks[idx_from_local_id( id )];
00258     }
00259 
00260     int count_from_rank( int rank ) const
00261     {
00262         return std::count( partRanks.begin(), partRanks.end(), rank );
00263     }
00264 
00265     void part_id_from_rank( int rank, std::vector< iMeshP_Part >& parts ) const;
00266 
00267     void local_id_from_rank( int rank, std::vector< int >& ids ) const;
00268 
00269     const std::vector< iMeshP_Part >& get_parts() const
00270     {
00271         return sortedPartList;
00272     }
00273 
00274     const std::vector< int >& get_ranks() const
00275     {
00276         return partRanks;
00277     }
00278 
00279     int build_map( iMesh_Instance imesh, iMeshP_PartitionHandle partition, int num_expected_parts );
00280 
00281     static int part_from_coords( iMesh_Instance imesh, iMeshP_PartHandle part, int& id_out );
00282 
00283   private:
00284     inline int idx_from_part_id( iMeshP_Part id ) const
00285     {
00286         return std::lower_bound( sortedPartList.begin(), sortedPartList.end(), id ) - sortedPartList.begin();
00287     }
00288     inline int idx_from_local_id( int id ) const
00289     {
00290         return localIdReverseMap[id];
00291     }
00292 
00293     std::vector< iMeshP_Part > sortedPartList;
00294     std::vector< int > partRanks;
00295     std::vector< int > partLocalIds;
00296     std::vector< int > localIdReverseMap;
00297 };
00298 
00299 /**\brief Create mesh for use in parallel tests */
00300 int create_mesh( const char* filename, int num_parts );
00301 
00302 int create_mesh_in_memory( int rank, int size, iMesh_Instance imesh, iMeshP_PartitionHandle& prtn, PartMap& map );
00303 
00304 /**\brief get unique identifier for each vertex */
00305 int vertex_tag( iMesh_Instance imesh, iBase_EntityHandle vertex, int& tag );
00306 
00307 int get_local_parts( iMesh_Instance instance, iMeshP_PartitionHandle prtn, std::vector< iMeshP_PartHandle >& handles,
00308                      std::vector< iMeshP_Part >* ids = 0 )
00309 {
00310     iMeshP_PartHandle* arr = 0;
00311     int ierr, alloc = 0, size;
00312     iMeshP_getLocalParts( instance, prtn, &arr, &alloc, &size, &ierr );CHKERR;
00313     handles.resize( size );
00314     std::copy( arr, arr + size, handles.begin() );
00315     free( arr );
00316     if( !ids ) return iBase_SUCCESS;
00317 
00318     ids->resize( size );
00319     alloc            = size;
00320     iMeshP_Part* ptr = &( *ids )[0];
00321     iMeshP_getPartIdsFromPartHandlesArr( instance, prtn, &handles[0], handles.size(), &ptr, &alloc, &size, &ierr );CHKERR;
00322     assert( size == (int)ids->size() );
00323     assert( ptr == &( *ids )[0] );
00324     return iBase_SUCCESS;
00325 }
00326 
00327 static int get_entities( iMesh_Instance imesh, iBase_EntitySetHandle set, iBase_EntityType type,
00328                          iMesh_EntityTopology topo, std::vector< iBase_EntityHandle >& entities )
00329 {
00330     iBase_EntityHandle* array = 0;
00331     int junk = 0, size = 0, err;
00332     iMesh_getEntities( imesh, set, type, topo, &array, &junk, &size, &err );
00333     if( !err )
00334     {
00335         entities.clear();
00336         entities.resize( size );
00337         std::copy( array, array + size, entities.begin() );
00338         free( array );
00339     }
00340     return err;
00341 }
00342 
00343 static int get_part_quads_and_verts( iMesh_Instance imesh, iMeshP_PartHandle part,
00344                                      std::vector< iBase_EntityHandle >& elems,
00345                                      std::vector< iBase_EntityHandle >& verts )
00346 {
00347     int ierr = get_entities( imesh, part, iBase_FACE, iMesh_QUADRILATERAL, elems );CHKERR;
00348 
00349     verts.resize( 4 * elems.size() );
00350     std::vector< int > junk( elems.size() + 1 );
00351     int junk1 = verts.size(), count, junk2 = junk.size(), junk3;
00352     iBase_EntityHandle* junk4 = &verts[0];
00353     int* junk5                = &junk[0];
00354     iMesh_getEntArrAdj( imesh, &elems[0], elems.size(), iBase_VERTEX, &junk4, &junk1, &count, &junk5, &junk2, &junk3,
00355                         &ierr );CHKERR;
00356     assert( junk1 == (int)verts.size() );
00357     assert( count == (int)( 4 * elems.size() ) );
00358     assert( junk2 == (int)junk.size() );
00359     assert( junk4 == &verts[0] );
00360     assert( junk5 == &junk[0] );
00361     std::sort( verts.begin(), verts.end() );
00362     verts.erase( std::unique( verts.begin(), verts.end() ), verts.end() );
00363     return iBase_SUCCESS;
00364 }
00365 
00366 static int get_coords( iMesh_Instance imesh, const iBase_EntityHandle* verts, int num_verts, double* coords )
00367 {
00368     double* junk1 = coords;
00369     int junk2     = 3 * num_verts;
00370     int junk3;
00371     int ierr;
00372     iMesh_getVtxArrCoords( imesh, verts, num_verts, iBase_INTERLEAVED, &junk1, &junk2, &junk3, &ierr );
00373     if( iBase_SUCCESS != ierr ) return ierr;
00374     assert( junk1 == coords );
00375     assert( junk2 == 3 * num_verts );
00376     assert( junk3 == 3 * num_verts );
00377     return iBase_SUCCESS;
00378 }
00379 
00380 /**************************************************************************
00381                               Main Method
00382  **************************************************************************/
00383 
00384 #define RUN_TEST( A ) run_test( &A, #A )
00385 
00386 int run_test( int ( *func )( iMesh_Instance, iMeshP_PartitionHandle, const PartMap& ), const char* func_name )
00387 {
00388     int rank, size, ierr;
00389     MPI_Comm_rank( MPI_COMM_WORLD, &rank );
00390     MPI_Comm_size( MPI_COMM_WORLD, &size );
00391     iMesh_Instance imesh;
00392     iMesh_newMesh( 0, &imesh, &ierr, 0 );
00393     PCHECK;
00394 
00395     iMeshP_PartitionHandle prtn;
00396     iMeshP_createPartitionAll( imesh, MPI_COMM_WORLD, &prtn, &ierr );
00397     PCHECK;
00398 
00399     PartMap map;
00400 
00401 #ifdef MOAB_HAVE_HDF5
00402     if( rank == 0 ) { ierr = create_mesh( FILENAME, size ); }
00403     MPI_Bcast( &ierr, 1, MPI_INT, 0, MPI_COMM_WORLD );
00404     if( ierr )
00405     {
00406         if( rank == 0 ) { std::cerr << "Failed to create input test file on root processor.  Aborting." << std::endl; }
00407         abort();
00408     }
00409 
00410     ierr = test_load( imesh, prtn, map, size );
00411     if( ierr )
00412     {
00413         if( rank == 0 )
00414         {
00415             std::cerr << "Failed to load input mesh." << std::endl
00416                       << "Cannot run further tests." << std::endl
00417                       << "ABORTING" << std::endl;
00418         }
00419         abort();
00420     }
00421 #else
00422     // so we have MPI and no HDF5; in order to run the test we need to create the
00423     // model in memory, and then call sync to resolve shared ents, as if it was read
00424     ierr = create_mesh_in_memory( rank, size, imesh, prtn, map );
00425     MPI_Bcast( &ierr, 1, MPI_INT, 0, MPI_COMM_WORLD );
00426     if( ierr )
00427     {
00428         if( rank == 0 ) { std::cerr << "Failed to create mesh.  Aborting." << std::endl; }
00429         abort();
00430     }
00431 
00432 #endif
00433     int result = ( *func )( imesh, prtn, map );
00434     int is_err = is_any_proc_error( result );
00435     if( rank == 0 )
00436     {
00437         if( is_err )
00438             std::cout << func_name << " : FAILED!!" << std::endl;
00439         else
00440             std::cout << func_name << " : success" << std::endl;
00441     }
00442 
00443     iMesh_dtor( imesh, &ierr );CHKERR;
00444     return is_err;
00445 }
00446 
00447 int main( int argc, char* argv[] )
00448 {
00449     MPI_Init( &argc, &argv );
00450     int size, rank;
00451     MPI_Comm_rank( MPI_COMM_WORLD, &rank );
00452     MPI_Comm_size( MPI_COMM_WORLD, &size );
00453 
00454     if( argc > 2 && !strcmp( argv[1], "-p" ) )
00455     {
00456 #if !defined( _MSC_VER ) && !defined( __MINGW32__ )
00457         std::cout << "Processor " << rank << " of " << size << " with PID " << getpid() << std::endl;
00458         std::cout.flush();
00459 #endif
00460         // loop forever on requested processor, giving the user time
00461         // to attach a debugger.  Once the debugger in attached, user
00462         // can change 'pause'.  E.g. on gdb do "set var pause = 0"
00463         if( atoi( argv[2] ) == rank )
00464         {
00465             volatile int pause = 1;
00466             while( pause )
00467                 ;
00468         }
00469         MPI_Barrier( MPI_COMM_WORLD );
00470     }
00471 
00472     int num_errors = 0;
00473     num_errors += RUN_TEST( test_get_partitions );
00474     num_errors += RUN_TEST( test_get_parts );
00475     num_errors += RUN_TEST( test_get_by_type );
00476     num_errors += RUN_TEST( test_get_by_topo );
00477     num_errors += RUN_TEST( test_part_id_handle );
00478     num_errors += RUN_TEST( test_part_rank );
00479     num_errors += RUN_TEST( test_get_neighbors );
00480     num_errors += RUN_TEST( test_get_part_boundary );
00481     num_errors += RUN_TEST( test_part_boundary_iter );
00482     num_errors += RUN_TEST( test_get_adjacencies );
00483     num_errors += RUN_TEST( test_entity_iterator );
00484     num_errors += RUN_TEST( test_entity_owner );
00485     num_errors += RUN_TEST( test_entity_status );
00486     num_errors += RUN_TEST( test_entity_copy_parts );
00487     num_errors += RUN_TEST( test_entity_copies );
00488     num_errors += RUN_TEST( test_push_tag_data_iface );
00489     num_errors += RUN_TEST( test_push_tag_data_ghost );
00490     num_errors += RUN_TEST( test_create_ghost_ents );
00491     num_errors += RUN_TEST( test_exchange_ents );
00492 
00493     // wait until all procs are done before writing summary data
00494     std::cout.flush();
00495     MPI_Barrier( MPI_COMM_WORLD );
00496 
00497 #ifdef MOAB_HAVE_HDF5
00498     // clean up output file
00499     if( rank == 0 ) remove( FILENAME );
00500 #endif
00501 
00502     if( rank == 0 )
00503     {
00504         if( !num_errors )
00505             std::cout << "All tests passed" << std::endl;
00506         else
00507             std::cout << num_errors << " TESTS FAILED!" << std::endl;
00508     }
00509 
00510     MPI_Finalize();
00511     return num_errors;
00512 }
00513 
00514 // Create a mesh
00515 //
00516 //
00517 // Groups of four quads will be arranged into parts as follows:
00518 // +------+------+------+------+------+-----
00519 // |             |             |
00520 // |             |             |
00521 // +    Part 0   +    Part 2   +    Part 4
00522 // |             |             |
00523 // |             |             |
00524 // +------+------+------+------+------+-----
00525 // |             |             |
00526 // |             |             |
00527 // +    Part 1   +    Part 3   +    Part 5
00528 // |             |             |
00529 // |             |             |
00530 // +------+------+------+------+------+-----
00531 //
00532 // Vertices will be enumerated as follows:
00533 // 1------6-----11-----16-----21-----26-----
00534 // |             |             |
00535 // |             |             |
00536 // 2      7     12     17     22     27
00537 // |             |             |
00538 // |             |             |
00539 // 3------8-----13-----18-----23-----28-----
00540 // |             |             |
00541 // |             |             |
00542 // 4      9     14     19     24     29
00543 // |             |             |
00544 // |             |             |
00545 // 5-----10-----15-----20-----25-----30-----
00546 //
00547 // Element IDs will be [4*rank+1,4*rank+5]
00548 template < int size >
00549 struct EHARR
00550 {
00551     iBase_EntityHandle h[size];
00552     iBase_EntityHandle& operator[]( int i )
00553     {
00554         return h[i];
00555     }
00556     operator iBase_EntityHandle*()
00557     {
00558         return h;
00559     }
00560 };
00561 int create_mesh( const char* filename, int num_parts )
00562 {
00563     const char* tagname = "GLOBAL_ID";
00564     int ierr;
00565 
00566     iMesh_Instance imesh;
00567     iMesh_newMesh( 0, &imesh, &ierr, 0 );CHKERR;
00568 
00569     const int num_full_cols  = 2 * ( num_parts / 2 );
00570     const int need_half_cols = num_parts % 2;
00571     const int num_cols       = num_full_cols + 2 * need_half_cols;
00572     const int num_vtx        = 5 + 5 * num_cols - 4 * ( num_parts % 2 );
00573     std::vector< EHARR< 5 > > vertices( num_cols + 1 );
00574     std::vector< EHARR< 4 > > elements( num_cols );
00575     std::vector< int > vertex_ids( num_vtx );
00576     std::vector< iBase_EntityHandle > vertex_list( num_vtx );
00577     for( int i = 0; i < num_vtx; ++i )
00578         vertex_ids[i] = i + 1;
00579 
00580     // create vertices
00581     int vl_pos = 0;
00582     for( int i = 0; i <= num_cols; ++i )
00583     {
00584         double coords[15]       = { static_cast< double >( i ), 0, 0, static_cast< double >( i ), 1, 0,
00585                               static_cast< double >( i ), 2, 0, static_cast< double >( i ), 3, 0,
00586                               static_cast< double >( i ), 4, 0 };
00587         iBase_EntityHandle* ptr = vertices[i];
00588         const int n             = ( num_full_cols == num_cols || i <= num_full_cols ) ? 5 : 3;
00589         int junk1 = n, junk2 = n;
00590         iMesh_createVtxArr( imesh, n, iBase_INTERLEAVED, coords, 3 * n, &ptr, &junk1, &junk2, &ierr );CHKERR;
00591         assert( ptr == vertices[i] );
00592         assert( junk1 == n );
00593         assert( junk2 == n );
00594         for( int j = 0; j < n; ++j )
00595             vertex_list[vl_pos++] = vertices[i][j];
00596     }
00597 
00598     // create elements
00599     for( int i = 0; i < num_cols; ++i )
00600     {
00601         iBase_EntityHandle conn[16];
00602         for( int j = 0; j < 4; ++j )
00603         {
00604             conn[4 * j]     = vertices[i][j];
00605             conn[4 * j + 1] = vertices[i][j + 1];
00606             conn[4 * j + 2] = vertices[i + 1][j + 1];
00607             conn[4 * j + 3] = vertices[i + 1][j];
00608         }
00609         iBase_EntityHandle* ptr = elements[i];
00610         const int n             = ( i < num_full_cols ) ? 4 : 2;
00611         int junk1 = n, junk2 = n, junk3 = n, junk4 = n;
00612         int stat[4];
00613         int* ptr2 = stat;
00614         iMesh_createEntArr( imesh, iMesh_QUADRILATERAL, conn, 4 * n, &ptr, &junk1, &junk2, &ptr2, &junk3, &junk4,
00615                             &ierr );CHKERR;
00616         assert( ptr == elements[i] );
00617         assert( junk1 == n );
00618         assert( junk2 == n );
00619         assert( ptr2 == stat );
00620         assert( junk3 == n );
00621         assert( junk4 == n );
00622     }
00623 
00624     // create partition
00625     iMeshP_PartitionHandle partition;
00626     iMeshP_createPartitionAll( imesh, MPI_COMM_SELF, &partition, &ierr );CHKERR;
00627     for( int i = 0; i < num_parts; ++i )
00628     {
00629         iMeshP_PartHandle part;
00630         iMeshP_createPart( imesh, partition, &part, &ierr );CHKERR;
00631         iBase_EntityHandle quads[] = { elements[2 * ( i / 2 )][2 * ( i % 2 )],
00632                                        elements[2 * ( i / 2 ) + 1][2 * ( i % 2 )],
00633                                        elements[2 * ( i / 2 )][2 * ( i % 2 ) + 1],
00634                                        elements[2 * ( i / 2 ) + 1][2 * ( i % 2 ) + 1] };
00635         iMesh_addEntArrToSet( imesh, quads, 4, part, &ierr );CHKERR;
00636     }
00637 
00638     // assign global ids to vertices
00639     iBase_TagHandle id_tag = 0;
00640     iMesh_getTagHandle( imesh, tagname, &id_tag, &ierr, strlen( tagname ) );
00641     if( iBase_SUCCESS == ierr )
00642     {
00643         int tag_size, tag_type;
00644         iMesh_getTagSizeValues( imesh, id_tag, &tag_size, &ierr );CHKERR;
00645         if( tag_size != 1 ) return iBase_TAG_ALREADY_EXISTS;
00646         iMesh_getTagType( imesh, id_tag, &tag_type, &ierr );CHKERR;
00647         if( tag_type != iBase_INTEGER ) return iBase_TAG_ALREADY_EXISTS;
00648     }
00649     else
00650     {
00651         iMesh_createTag( imesh, tagname, 1, iBase_INTEGER, &id_tag, &ierr, strlen( tagname ) );CHKERR;
00652     }
00653     iMesh_setIntArrData( imesh, &vertex_list[0], num_vtx, id_tag, &vertex_ids[0], num_vtx, &ierr );CHKERR;
00654 
00655     // write file
00656     iBase_EntitySetHandle root_set;
00657     iMesh_getRootSet( imesh, &root_set, &ierr );
00658     iMeshP_saveAll( imesh, partition, root_set, filename, 0, &ierr, strlen( filename ), 0 );CHKERR;
00659 
00660     iMesh_dtor( imesh, &ierr );CHKERR;
00661 
00662     return 0;
00663 }
00664 int create_mesh_in_memory( int rank, int num_parts, iMesh_Instance imesh, iMeshP_PartitionHandle& partition,
00665                            PartMap& map )
00666 {
00667     const char* tagname = "GLOBAL_ID";
00668     int ierr;
00669 
00670     const int num_cols = 2;
00671     const int num_vtx  = 9;
00672     // we are on the top or botton row
00673     int bottom = rank % 2;  // 0  2
00674                             // 1  3
00675     std::vector< EHARR< 3 > > vertices( 3 );
00676     std::vector< EHARR< 2 > > elements( 2 );  // 4 elements per process
00677     std::vector< int > vertex_ids( num_vtx );
00678     std::vector< iBase_EntityHandle > vertex_list( num_vtx );
00679     int start = 1 + 2 * bottom + 10 * ( rank / 2 );
00680 
00681     for( int i = 0; i < 3; ++i )
00682         for( int j = 0; j < 3; ++j )
00683         {
00684             vertex_ids[i + 3 * j] = start + i + 5 * j;
00685         }
00686 
00687     // create vertices
00688     int vl_pos = 0;
00689     int startI = 2 * ( rank / 2 );  // so it will be 0, 0, 2, 2, ...)
00690     for( int i = 0; i <= 2; ++i )
00691     {
00692         double coords[9] = {
00693             static_cast< double >( i + startI ), 2. * bottom,     0,
00694             static_cast< double >( i + startI ), 1 + 2. * bottom, 0,
00695             static_cast< double >( i + startI ), 2 + 2. * bottom, 0,
00696         };
00697         iBase_EntityHandle* ptr = vertices[i];
00698         const int n             = 3;
00699         int junk1 = n, junk2 = n;
00700         iMesh_createVtxArr( imesh, n, iBase_INTERLEAVED, coords, 3 * n, &ptr, &junk1, &junk2, &ierr );CHKERR;
00701         assert( ptr == vertices[i] );
00702         assert( junk1 == n );
00703         assert( junk2 == n );
00704         for( int j = 0; j < n; ++j )
00705             vertex_list[vl_pos++] = vertices[i][j];
00706     }
00707 
00708     // create elements
00709     for( int i = 0; i < num_cols; ++i )
00710     {
00711         iBase_EntityHandle conn[8];
00712         for( int j = 0; j < 2; ++j )
00713         {
00714             conn[4 * j]     = vertices[i][j];
00715             conn[4 * j + 1] = vertices[i][j + 1];
00716             conn[4 * j + 2] = vertices[i + 1][j + 1];
00717             conn[4 * j + 3] = vertices[i + 1][j];
00718         }
00719         iBase_EntityHandle* ptr = elements[i];
00720         const int n             = 2;
00721         int junk1 = n, junk2 = n, junk3 = n, junk4 = n;
00722         int stat[4];
00723         int* ptr2 = stat;
00724         iMesh_createEntArr( imesh, iMesh_QUADRILATERAL, conn, 4 * n, &ptr, &junk1, &junk2, &ptr2, &junk3, &junk4,
00725                             &ierr );CHKERR;
00726         assert( ptr == elements[i] );
00727         assert( junk1 == n );
00728         assert( junk2 == n );
00729         assert( ptr2 == stat );
00730         assert( junk3 == n );
00731         assert( junk4 == n );
00732     }
00733 
00734     // create partition
00735     iMeshP_createPartitionAll( imesh, MPI_COMM_WORLD, &partition, &ierr );CHKERR;
00736 
00737     iMeshP_PartHandle part;
00738     iMeshP_createPart( imesh, partition, &part, &ierr );CHKERR;
00739     iBase_EntityHandle quads[] = { elements[0][0], elements[0][1], elements[1][0], elements[1][1] };
00740     iMesh_addEntArrToSet( imesh, quads, 4, part, &ierr );CHKERR;
00741 
00742     // assign global ids to vertices
00743     iBase_TagHandle id_tag = 0;
00744     iMesh_getTagHandle( imesh, tagname, &id_tag, &ierr, strlen( tagname ) );
00745     if( iBase_SUCCESS == ierr )
00746     {
00747         int tag_size, tag_type;
00748         iMesh_getTagSizeValues( imesh, id_tag, &tag_size, &ierr );CHKERR;
00749         if( tag_size != 1 ) return iBase_TAG_ALREADY_EXISTS;
00750         iMesh_getTagType( imesh, id_tag, &tag_type, &ierr );CHKERR;
00751         if( tag_type != iBase_INTEGER ) return iBase_TAG_ALREADY_EXISTS;
00752     }
00753     else
00754     {
00755         iMesh_createTag( imesh, tagname, 1, iBase_INTEGER, &id_tag, &ierr, strlen( tagname ) );CHKERR;
00756     }
00757     iMesh_setIntArrData( imesh, &vertex_list[0], num_vtx, id_tag, &vertex_ids[0], num_vtx, &ierr );CHKERR;
00758 
00759     // some mesh sync
00760     iMeshP_syncPartitionAll( imesh, partition, &ierr );CHKERR;
00761     iMeshP_syncMeshAll( imesh, partition, &ierr );CHKERR;
00762 
00763     ierr = map.build_map( imesh, partition, num_parts );CHKERR;
00764     return 0;
00765 }
00766 // generate unique for each vertex from coordinates.
00767 // Assume integer coordinate values with x in [0,inf] and y in [0,4]
00768 // as generated by create_mean(..).
00769 int vertex_tag( iMesh_Instance imesh, iBase_EntityHandle vertex, int& tag )
00770 {
00771     int ierr;
00772     double x, y, z;
00773     iMesh_getVtxCoord( imesh, vertex, &x, &y, &z, &ierr );CHKERR;
00774 
00775     int xc = (int)round( x );
00776     int yc = (int)round( y );
00777     tag    = 5 * xc + yc + 1;
00778     return ierr;
00779 }
00780 
00781 /**************************************************************************
00782                            Test  Implementations
00783  **************************************************************************/
00784 
00785 int test_load( iMesh_Instance imesh, iMeshP_PartitionHandle prtn, PartMap& map, int proc_size )
00786 {
00787     int ierr;
00788 
00789     iBase_EntitySetHandle root_set;
00790     iMesh_getRootSet( imesh, &root_set, &ierr );
00791     const char* opt = "moab:PARTITION=PARALLEL_PARTITION";
00792     iMeshP_loadAll( imesh, prtn, root_set, FILENAME, opt, &ierr, strlen( FILENAME ), strlen( opt ) );
00793     PCHECK;
00794 
00795     ierr = map.build_map( imesh, prtn, proc_size );CHKERR;
00796     return iBase_SUCCESS;
00797 }
00798 
00799 /**\brief Test partition query methods
00800  *
00801  * Test:
00802  * - iMeshP_getPartitionComm
00803  * - iMeshP_getNumPartitions
00804  * - iMeshP_getPartitions
00805  */
00806 int test_get_partitions( iMesh_Instance imesh, iMeshP_PartitionHandle prtn, const PartMap& )
00807 {
00808     int ierr;
00809 
00810     // test iMeshP_getPartitionCom
00811     MPI_Comm comm = MPI_COMM_SELF;
00812     iMeshP_getPartitionComm( imesh, prtn, &comm, &ierr );
00813     PCHECK;
00814     ASSERT( comm == MPI_COMM_WORLD );
00815 
00816     // test iMeshP_getPartitions
00817     iMeshP_PartitionHandle* array = 0;
00818     int alloc = 0, size = -1;
00819     iMeshP_getPartitions( imesh, &array, &alloc, &size, &ierr );
00820     PCHECK;
00821     ASSERT( array != 0 );
00822     ASSERT( alloc == size );
00823     ASSERT( size > 0 );
00824     int idx = std::find( array, array + size, prtn ) - array;
00825     free( array );
00826     ASSERT( idx < size );
00827 
00828     // test iMesP_getNumPartitions
00829     int size2 = -1;
00830     iMeshP_getNumPartitions( imesh, &size2, &ierr );
00831     PCHECK;
00832     ASSERT( size2 == size );
00833     return 0;
00834 }
00835 
00836 /**\brief Test part quyery methods
00837  *
00838  * Test:
00839  * - iMeshP_getNumGlobalParts
00840  * - iMeshP_getNumLocalParts
00841  * - iMeshP_getLocalParts
00842  */
00843 int test_get_parts( iMesh_Instance imesh, iMeshP_PartitionHandle prtn, const PartMap& map )
00844 {
00845     int size, rank, ierr;
00846     MPI_Comm_rank( MPI_COMM_WORLD, &rank );
00847     MPI_Comm_size( MPI_COMM_WORLD, &size );
00848 
00849     int num_part_g;
00850     iMeshP_getNumGlobalParts( imesh, prtn, &num_part_g, &ierr );
00851     PCHECK;
00852     ASSERT( num_part_g == map.num_parts() );
00853 
00854     int num_part_l;
00855     iMeshP_getNumLocalParts( imesh, prtn, &num_part_l, &ierr );
00856     PCHECK;
00857     ASSERT( num_part_l == map.count_from_rank( rank ) );
00858 
00859     std::vector< iMeshP_PartHandle > parts( num_part_l );
00860     iMeshP_PartHandle* ptr = &parts[0];
00861     int junk1 = num_part_l, count = -1;
00862     iMeshP_getLocalParts( imesh, prtn, &ptr, &junk1, &count, &ierr );
00863     PCHECK;
00864     assert( ptr == &parts[0] );
00865     assert( junk1 == num_part_l );
00866     ASSERT( count == num_part_l );
00867 
00868     return iBase_SUCCESS;
00869 }
00870 
00871 static int test_get_by_type_topo_all( iMesh_Instance imesh, iMeshP_PartitionHandle prtn, bool test_type, int num_parts )
00872 {
00873     // calculate number of quads and vertices in entire mesh
00874     // from number of parts (see create_mesh(..) function.)
00875     const int expected_global_quad_count = 4 * num_parts;
00876     const int num_col                    = 2 * ( num_parts / 2 + num_parts % 2 );
00877     const int expected_global_vtx_count  = num_parts == 1 ? 9 : num_parts % 2 ? 1 + 5 * num_col : 5 + 5 * num_col;
00878 
00879     // test getNumOf*All for root set
00880     int ierr, count;
00881     iBase_EntitySetHandle root;
00882     iMesh_getRootSet( imesh, &root, &ierr );
00883     if( test_type )
00884         iMeshP_getNumOfTypeAll( imesh, prtn, root, iBase_VERTEX, &count, &ierr );
00885     else
00886         iMeshP_getNumOfTopoAll( imesh, prtn, root, iMesh_POINT, &count, &ierr );
00887     PCHECK;
00888     ASSERT( count == expected_global_vtx_count );
00889     if( test_type )
00890         iMeshP_getNumOfTypeAll( imesh, prtn, root, iBase_FACE, &count, &ierr );
00891     else
00892         iMeshP_getNumOfTopoAll( imesh, prtn, root, iMesh_QUADRILATERAL, &count, &ierr );
00893     PCHECK;
00894     ASSERT( count == expected_global_quad_count );
00895 
00896     // create an entity set containing half of the quads
00897     std::vector< iBase_EntityHandle > all_quads, half_quads;
00898     ierr = get_entities( imesh, root, iBase_FACE, iMesh_QUADRILATERAL, all_quads );
00899     assert( 0 == all_quads.size() % 2 );
00900     half_quads.resize( all_quads.size() / 2 );
00901     for( size_t i = 0; i < all_quads.size() / 2; ++i )
00902         half_quads[i] = all_quads[2 * i];
00903     iBase_EntitySetHandle set;
00904     iMesh_createEntSet( imesh, 1, &set, &ierr );CHKERR;
00905     iMesh_addEntArrToSet( imesh, &half_quads[0], half_quads.size(), set, &ierr );CHKERR;
00906 
00907     // test getNumOf*All with defined set
00908     if( test_type )
00909         iMeshP_getNumOfTypeAll( imesh, prtn, set, iBase_VERTEX, &count, &ierr );
00910     else
00911         iMeshP_getNumOfTopoAll( imesh, prtn, set, iMesh_POINT, &count, &ierr );
00912     PCHECK;
00913     ASSERT( count == 0 );
00914     if( test_type )
00915         iMeshP_getNumOfTypeAll( imesh, prtn, set, iBase_FACE, &count, &ierr );
00916     else
00917         iMeshP_getNumOfTopoAll( imesh, prtn, set, iMesh_QUADRILATERAL, &count, &ierr );
00918     PCHECK;
00919     ASSERT( count == expected_global_quad_count / 2 );
00920 
00921     return 0;
00922 }
00923 
00924 static int test_get_by_type_topo_local( iMesh_Instance imesh, iMeshP_PartitionHandle prtn, bool test_type )
00925 {
00926     int ierr;
00927     iBase_EntitySetHandle root;
00928     iMesh_getRootSet( imesh, &root, &ierr );
00929 
00930     // select a single part
00931     std::vector< iMeshP_PartHandle > parts;
00932     ierr = get_local_parts( imesh, prtn, parts );CHKERR;
00933     iMeshP_PartHandle part = parts.front();
00934 
00935     // get the entities contained in the part
00936     std::vector< iBase_EntityHandle > part_quads, part_all;
00937     ierr = get_entities( imesh, part, iBase_FACE, iMesh_QUADRILATERAL, part_quads );CHKERR;
00938     ierr = get_entities( imesh, part, iBase_ALL_TYPES, iMesh_ALL_TOPOLOGIES, part_all );CHKERR;
00939 
00940     // compare local counts (using root set)
00941 
00942     int count;
00943     if( test_type )
00944         iMeshP_getNumOfType( imesh, prtn, part, root, iBase_FACE, &count, &ierr );
00945     else
00946         iMeshP_getNumOfTopo( imesh, prtn, part, root, iMesh_QUADRILATERAL, &count, &ierr );CHKERR;
00947     ASSERT( count == (int)part_quads.size() );
00948 
00949     if( test_type )
00950         iMeshP_getNumOfType( imesh, prtn, part, root, iBase_ALL_TYPES, &count, &ierr );
00951     else
00952         iMeshP_getNumOfTopo( imesh, prtn, part, root, iMesh_ALL_TOPOLOGIES, &count, &ierr );CHKERR;
00953     ASSERT( count == (int)part_all.size() );
00954 
00955     // compare local contents (using root set)
00956 
00957     iBase_EntityHandle* ptr = 0;
00958     int num_ent, junk1 = 0;
00959     iMeshP_getEntities( imesh, prtn, part, root, test_type ? iBase_FACE : iBase_ALL_TYPES,
00960                         test_type ? iMesh_ALL_TOPOLOGIES : iMesh_QUADRILATERAL, &ptr, &junk1, &num_ent, &ierr );CHKERR;
00961     std::vector< iBase_EntityHandle > act_quads( ptr, ptr + num_ent );
00962     free( ptr );
00963     junk1 = num_ent = 0;
00964     ptr             = 0;
00965     iMeshP_getEntities( imesh, prtn, part, root, iBase_ALL_TYPES, iMesh_ALL_TOPOLOGIES, &ptr, &junk1, &num_ent, &ierr );CHKERR;
00966     std::vector< iBase_EntityHandle > act_all( ptr, ptr + num_ent );
00967     free( ptr );
00968     std::sort( part_quads.begin(), part_quads.end() );
00969     std::sort( part_all.begin(), part_all.end() );
00970     std::sort( act_quads.begin(), act_quads.end() );
00971     std::sort( act_all.begin(), act_all.end() );
00972     ASSERT( part_quads == act_quads );
00973     ASSERT( part_all == act_all );
00974 
00975     // create an entity set containing half of the quads from the part
00976     std::vector< iBase_EntityHandle > half_quads( part_quads.size() / 2 );
00977     for( size_t i = 0; i < half_quads.size(); ++i )
00978         half_quads[i] = part_quads[2 * i];
00979     iBase_EntitySetHandle set;
00980     iMesh_createEntSet( imesh, 1, &set, &ierr );CHKERR;
00981     iMesh_addEntArrToSet( imesh, &half_quads[0], half_quads.size(), set, &ierr );CHKERR;
00982 
00983     // check if there exists any quads not in the part that we
00984     // can add to the set
00985     std::vector< iBase_EntityHandle > all_quads, other_quads;
00986     ierr = get_entities( imesh, root, iBase_FACE, iMesh_QUADRILATERAL, all_quads );CHKERR;
00987     std::sort( all_quads.begin(), all_quads.end() );
00988     std::sort( part_quads.begin(), part_quads.end() );
00989     std::set_difference( all_quads.begin(), all_quads.end(), part_quads.begin(), part_quads.end(),
00990                          std::back_inserter( other_quads ) );
00991     iMesh_addEntArrToSet( imesh, &other_quads[0], other_quads.size(), set, &ierr );CHKERR;
00992 
00993     // compare local counts (using non-root set)
00994 
00995     if( test_type )
00996         iMeshP_getNumOfType( imesh, prtn, part, set, iBase_FACE, &count, &ierr );
00997     else
00998         iMeshP_getNumOfTopo( imesh, prtn, part, set, iMesh_QUADRILATERAL, &count, &ierr );CHKERR;
00999     ASSERT( count == (int)half_quads.size() );
01000 
01001     if( test_type )
01002         iMeshP_getNumOfType( imesh, prtn, part, set, iBase_VERTEX, &count, &ierr );
01003     else
01004         iMeshP_getNumOfTopo( imesh, prtn, part, set, iMesh_POINT, &count, &ierr );CHKERR;
01005     ASSERT( count == 0 );
01006 
01007     // compare local contents (using non-root set)
01008 
01009     junk1   = 0;
01010     num_ent = 0;
01011     ptr     = 0;
01012     iMeshP_getEntities( imesh, prtn, part, set, test_type ? iBase_FACE : iBase_ALL_TYPES,
01013                         test_type ? iMesh_ALL_TOPOLOGIES : iMesh_QUADRILATERAL, &ptr, &junk1, &num_ent, &ierr );CHKERR;
01014     act_quads.resize( num_ent );
01015     std::copy( ptr, ptr + num_ent, act_quads.begin() );
01016     free( ptr );
01017     std::sort( half_quads.begin(), half_quads.end() );
01018     std::sort( act_quads.begin(), act_quads.end() );
01019     ASSERT( act_quads == half_quads );
01020 
01021     return iBase_SUCCESS;
01022 }
01023 
01024 /**\brief Test query by entity type
01025  *
01026  * Test:
01027  * - iMeshP_getNumOfTypeAll
01028  * - iMeshP_getNumOfType
01029  * - iMeshP_getEntities
01030  * -
01031  */
01032 int test_get_by_type( iMesh_Instance imesh, iMeshP_PartitionHandle prtn, const PartMap& map )
01033 {
01034     int ierr;
01035     ierr = test_get_by_type_topo_all( imesh, prtn, true, map.num_parts() );
01036     PCHECK;
01037     ierr = test_get_by_type_topo_local( imesh, prtn, true );
01038     PCHECK;
01039     return 0;
01040 }
01041 
01042 /**\brief Test query by entity topology
01043  *
01044  * Test:
01045  * - iMeshP_getNumOfTopoAll
01046  * - iMeshP_getNumOfTopo
01047  * - iMeshP_getEntities
01048  * -
01049  */
01050 int test_get_by_topo( iMesh_Instance imesh, iMeshP_PartitionHandle prtn, const PartMap& map )
01051 {
01052     int ierr;
01053     ierr = test_get_by_type_topo_all( imesh, prtn, false, map.num_parts() );
01054     PCHECK;
01055     ierr = test_get_by_type_topo_local( imesh, prtn, false );
01056     PCHECK;
01057     return 0;
01058 }
01059 
01060 /**\brief Test mapping from part id to part handle
01061  *
01062  * Test:
01063  * - iMeshP_getPartIdFromPartHandle
01064  * - iMeshP_getPartIdsFromPartHandlesArr
01065  * - iMeshP_getPartHandleFromPartId
01066  * - iMeshP_getPartHandlesFromPartsIdsArr
01067  */
01068 int test_part_id_handle( iMesh_Instance imesh, iMeshP_PartitionHandle prtn, const PartMap& map )
01069 {
01070     // get local part ids
01071     int rank, ierr;
01072     MPI_Comm_rank( MPI_COMM_WORLD, &rank );
01073     std::vector< iMeshP_Part > ids;
01074     map.part_id_from_rank( rank, ids );
01075 
01076     // check single-part functions and build list of part handles
01077     std::vector< iMeshP_PartHandle > handles( ids.size() );
01078     size_t i;
01079     for( i = 0; i < ids.size(); ++i )
01080     {
01081         iMeshP_getPartHandleFromPartId( imesh, prtn, ids[i], &handles[i], &ierr );CHKERR;
01082         iMeshP_Part id;
01083         iMeshP_getPartIdFromPartHandle( imesh, prtn, handles[i], &id, &ierr );CHKERR;
01084         if( id != ids[i] ) break;
01085     }
01086     ASSERT( i == ids.size() );
01087 
01088     // test iMeshP_getPartIdsFromPartHandlesArr
01089     std::vector< iMeshP_Part > ids2( ids.size() );
01090     int junk1 = ids.size(), junk2 = 0;
01091     iMeshP_Part* ptr = &ids2[0];
01092     iMeshP_getPartIdsFromPartHandlesArr( imesh, prtn, &handles[0], handles.size(), &ptr, &junk1, &junk2, &ierr );
01093     PCHECK;
01094     ASSERT( ptr == &ids2[0] );
01095     ASSERT( junk2 == (int)ids2.size() );
01096     ASSERT( ids == ids2 );
01097 
01098     // test iMeshP_getPartHandlesFromPartsIdsArr
01099     std::vector< iMeshP_PartHandle > handles2( handles.size() );
01100     junk1                   = handles.size();
01101     junk2                   = 0;
01102     iMeshP_PartHandle* ptr2 = &handles2[0];
01103     iMeshP_getPartHandlesFromPartsIdsArr( imesh, prtn, &ids[0], ids.size(), &ptr2, &junk1, &junk2, &ierr );
01104     PCHECK;
01105     ASSERT( ptr2 == &handles2[0] );
01106     ASSERT( junk2 == (int)handles2.size() );
01107     ASSERT( handles == handles2 );
01108 
01109     return 0;
01110 }
01111 
01112 /**\brief Test get part rank
01113  *
01114  * Tests:
01115  * - iMeshP_getRankOfPart
01116  * - iMeshP_getRankOfPartArr
01117  */
01118 int test_part_rank( iMesh_Instance imesh, iMeshP_PartitionHandle prtn, const PartMap& map )
01119 {
01120     int ierr = 0, rank;
01121     std::vector< iMeshP_Part > invalid, failed;
01122     MPI_Comm_rank( MPI_COMM_WORLD, &rank );
01123 
01124     // test iMeshP_getRankOfPart
01125     for( size_t i = 0; i < map.get_parts().size(); ++i )
01126     {
01127         int pr;
01128         iMeshP_getRankOfPart( imesh, prtn, map.get_parts()[i], &pr, &ierr );
01129         if( iBase_SUCCESS != ierr )
01130             failed.push_back( map.get_parts()[i] );
01131         else if( pr != map.get_ranks()[i] )
01132             invalid.push_back( map.get_parts()[i] );
01133     }
01134     if( !failed.empty() )
01135     {
01136         std::cerr << "Processor " << rank << ": iMeshP_getRankOfPart failed for " << failed.size() << " parts."
01137                   << std::endl;
01138         ierr = iBase_FAILURE;
01139     }
01140     if( !invalid.empty() )
01141     {
01142         std::cerr << "Processor " << rank << ": iMeshP_getRankOfPart was incorrect for " << invalid.size() << " parts."
01143                   << std::endl;
01144         ierr = iBase_FAILURE;
01145     }
01146     PCHECK;
01147 
01148     // test iMeshP_getRankOfPartArr
01149     std::vector< int > ranks( map.get_parts().size() );
01150     int junk1 = ranks.size(), junk2, *ptr = &ranks[0];
01151     iMeshP_getRankOfPartArr( imesh, prtn, &map.get_parts()[0], map.get_parts().size(), &ptr, &junk1, &junk2, &ierr );
01152     PCHECK;
01153     assert( ptr == &ranks[0] );
01154     assert( junk1 == (int)ranks.size() );
01155     ASSERT( junk2 == (int)ranks.size() );
01156     for( size_t i = 0; i < map.get_parts().size(); ++i )
01157     {
01158         if( ranks[i] != map.get_ranks()[i] ) invalid.push_back( map.get_parts()[i] );
01159     }
01160     if( !invalid.empty() )
01161     {
01162         std::cerr << "Processor " << rank << ": iMeshP_getRankOfPartArr was incorrect for " << invalid.size()
01163                   << " parts." << std::endl;
01164         ierr = iBase_FAILURE;
01165     }
01166     PCHECK;
01167 
01168     return 0;
01169 }
01170 
01171 // see create_mesh(..)
01172 static void get_part_neighbors( int logical_part_id, int num_parts, int neighbors[5], int& num_neighbors )
01173 {
01174     num_neighbors = 0;
01175     if( logical_part_id + 1 < num_parts ) neighbors[num_neighbors++] = logical_part_id + 1;
01176     if( logical_part_id + 2 < num_parts ) neighbors[num_neighbors++] = logical_part_id + 2;
01177     if( logical_part_id % 2 )
01178     {
01179         neighbors[num_neighbors++] = logical_part_id - 1;
01180         if( logical_part_id > 2 )
01181         {
01182             neighbors[num_neighbors++] = logical_part_id - 3;
01183             neighbors[num_neighbors++] = logical_part_id - 2;
01184         }
01185     }
01186     else
01187     {
01188         if( logical_part_id + 3 < num_parts ) neighbors[num_neighbors++] = logical_part_id + 3;
01189         if( logical_part_id > 1 )
01190         {
01191             neighbors[num_neighbors++] = logical_part_id - 1;
01192             neighbors[num_neighbors++] = logical_part_id - 2;
01193         }
01194     }
01195 }
01196 
01197 /**\brief Test querying of part neighbors
01198  *
01199  * Test:
01200  * - iMeshP_getNumPartNbors
01201  * - iMeshP_getNumPartNborsArr
01202  * - iMeshP_getPartNbors
01203  * - iMeshP_getPartNborsArr
01204  */
01205 int test_get_neighbors( iMesh_Instance imesh, iMeshP_PartitionHandle prtn, const PartMap& map )
01206 {
01207     int ierr, rank;
01208     MPI_Comm_rank( MPI_COMM_WORLD, &rank );
01209 
01210     std::vector< iMeshP_Part > local_parts;
01211     map.part_id_from_rank( rank, local_parts );
01212 
01213     // get handles for local parts
01214     std::vector< iMeshP_PartHandle > handles( local_parts.size() );
01215     iMeshP_PartHandle* ptr = &handles[0];
01216     int junk1 = handles.size(), junk2 = 0;
01217     iMeshP_getPartHandlesFromPartsIdsArr( imesh, prtn, &local_parts[0], local_parts.size(), &ptr, &junk1, &junk2,
01218                                           &ierr );
01219     PCHECK;
01220     assert( ptr == &handles[0] );
01221     assert( junk2 == (int)handles.size() );
01222 
01223     // get logical ids for local parts
01224     std::vector< int > logical_ids;
01225     map.local_id_from_rank( rank, logical_ids );
01226 
01227     // get neighbors for each local part
01228     std::vector< std::vector< iMeshP_Part > > neighbors( logical_ids.size() );
01229     for( size_t i = 0; i < logical_ids.size(); ++i )
01230     {
01231         int logical_neighbors[5], num_neighbors;
01232         get_part_neighbors( logical_ids[i], map.num_parts(), logical_neighbors, num_neighbors );
01233         neighbors[i].resize( num_neighbors );
01234         for( int j = 0; j < num_neighbors; ++j )
01235             neighbors[i][j] = map.part_id_from_local_id( logical_neighbors[j] );
01236         std::sort( neighbors[i].begin(), neighbors[i].end() );
01237     }
01238 
01239     // test iMeshP_getNumPartNbors
01240     std::vector< iMeshP_Part > invalid, failed;
01241     for( size_t i = 0; i < local_parts.size(); ++i )
01242     {
01243         int count;
01244         iMeshP_getNumPartNbors( imesh, prtn, handles[i], iBase_VERTEX, &count, &ierr );
01245         if( ierr )
01246             failed.push_back( local_parts[i] );
01247         else if( count != (int)neighbors[i].size() )
01248             invalid.push_back( local_parts[i] );
01249     }
01250     if( !failed.empty() )
01251     {
01252         std::cerr << "Processor " << rank << ": iMeshP_getNumPartNbors failed for " << failed.size() << " parts."
01253                   << std::endl;
01254         ierr = iBase_FAILURE;
01255         PCHECK;
01256     }
01257     if( !invalid.empty() )
01258     {
01259         std::cerr << "Processor " << rank << ": iMeshP_getNumPartNbors was incorrect for " << invalid.size()
01260                   << " parts." << std::endl;
01261         ierr = iBase_FAILURE;
01262         PCHECK;
01263     }
01264 
01265     // test iMeshP_getPartNbors
01266     ierr = 0;
01267     for( size_t i = 0; i < local_parts.size(); ++i )
01268     {
01269         int count, junk = 0, another_count;
01270         iMeshP_Part* list = 0;
01271         iMeshP_getPartNbors( imesh, prtn, handles[i], iBase_VERTEX, &another_count, &list, &junk, &count, &ierr );
01272         assert( count == another_count );
01273         if( ierr )
01274             failed.push_back( local_parts[i] );
01275         else
01276         {
01277             std::sort( list, list + count );
01278             std::vector< iMeshP_Part > cpy( list, list + count );
01279             if( cpy != neighbors[i] ) invalid.push_back( local_parts[i] );
01280             free( list );
01281         }
01282     }
01283     if( !failed.empty() )
01284     {
01285         std::cerr << "Processor " << rank << ": iMeshP_getPartNbors failed for " << failed.size() << " parts."
01286                   << std::endl;
01287         ierr = iBase_FAILURE;
01288     }
01289     if( !invalid.empty() )
01290     {
01291         std::cerr << "Processor " << rank << ": iMeshP_getPartNbors was incorrect for " << invalid.size() << " parts."
01292                   << std::endl;
01293         ierr = iBase_FAILURE;
01294     }
01295     PCHECK;
01296 
01297     // test iMeshP_getNumPartNborsArr
01298     std::vector< int > count_vect( handles.size() );
01299     int* count_arr = &count_vect[0];
01300     junk1          = handles.size();
01301     iMeshP_getNumPartNborsArr( imesh, prtn, &handles[0], handles.size(), iBase_VERTEX, &count_arr, &junk1, &junk2,
01302                                &ierr );
01303     PCHECK;
01304     assert( count_arr == &count_vect[0] );
01305     assert( junk2 == (int)handles.size() );
01306     for( size_t i = 0; i < local_parts.size(); ++i )
01307     {
01308         if( count_arr[i] != (int)neighbors[i].size() ) invalid.push_back( local_parts[i] );
01309     }
01310     if( !invalid.empty() )
01311     {
01312         std::cerr << "Processor " << rank << ": iMeshP_getNumPartNborsArr was incorrect for " << invalid.size()
01313                   << " parts." << std::endl;
01314         ierr = iBase_FAILURE;
01315     }
01316     PCHECK;
01317 
01318     // test iMeshP_getPartNborsArr
01319     iMeshP_Part* nbor_arr = 0;
01320     junk1 = handles.size(), junk2 = 0;
01321     int junk3 = 0, nbor_size;
01322     iMeshP_getPartNborsArr( imesh, prtn, &handles[0], handles.size(), iBase_VERTEX, &count_arr, &junk1, &junk2,
01323                             &nbor_arr, &junk3, &nbor_size, &ierr );
01324     PCHECK;
01325     assert( count_arr == &count_vect[0] );
01326     assert( junk2 == (int)handles.size() );
01327     std::vector< iMeshP_Part > all_nbors( nbor_arr, nbor_arr + nbor_size );
01328     free( nbor_arr );
01329     std::vector< iMeshP_Part >::iterator j = all_nbors.begin();
01330     bool bad_length                        = false;
01331     for( size_t i = 0; i < local_parts.size(); ++i )
01332     {
01333         if( all_nbors.end() - j > count_arr[i] )
01334         {
01335             bad_length = true;
01336             break;
01337         }
01338         if( count_arr[i] != (int)neighbors[i].size() ) { invalid.push_back( local_parts[i] ); }
01339         else
01340         {
01341             std::vector< iMeshP_Part >::iterator e = j + count_arr[i];
01342             std::sort( j, e );
01343             if( !std::equal( j, e, neighbors[i].begin() ) ) invalid.push_back( local_parts[i] );
01344         }
01345     }
01346     if( bad_length )
01347     {
01348         std::cerr << "Processor " << rank << ": iMeshP_getPartNborsArr had inconsistent result array lengths."
01349                   << std::endl;
01350         ierr = iBase_FAILURE;
01351     }
01352     if( !invalid.empty() )
01353     {
01354         std::cerr << "Processor " << rank << ": iMeshP_getPartNborsArr was incorrect for " << invalid.size()
01355                   << " parts." << std::endl;
01356         ierr = iBase_FAILURE;
01357     }
01358     PCHECK;
01359 
01360     return 0;
01361 }
01362 
01363 // Determine the expected vertices on the interface between two parts.
01364 // Returns no vertices for non-adjacient parts and fails if both parts
01365 // are the same.
01366 // See create_mesh(..) for the assumed mesh.
01367 static int interface_verts( iMesh_Instance imesh, iMeshP_PartitionHandle prtn, iMeshP_PartHandle local_part,
01368                             iMeshP_Part other_part, const PartMap& map, std::vector< iBase_EntityHandle >& vtx_handles )
01369 {
01370     int ierr, rank;
01371     MPI_Comm_rank( MPI_COMM_WORLD, &rank );
01372 
01373     iMeshP_Part local_id;
01374     iMeshP_getPartIdFromPartHandle( imesh, prtn, local_part, &local_id, &ierr );CHKERR;
01375 
01376     const int local_logical = map.local_id_from_part_id( local_id );
01377     const int other_logical = map.local_id_from_part_id( other_part );
01378 
01379     // get grid of local vertices
01380 
01381     iBase_EntityHandle verts[3][3];
01382     const double xbase = ( local_id / 2 ) * 2;
01383     const double ybase = ( local_id % 2 ) * 2;
01384 
01385     // get quads in partition
01386     iBase_EntityHandle quads[4], *ptr = quads;
01387     int junk1 = 4, junk2;
01388     iMesh_getEntities( imesh, local_part, iBase_FACE, iMesh_QUADRILATERAL, &ptr, &junk1, &junk2, &ierr );CHKERR;
01389     assert( ptr == quads );
01390     assert( junk1 == 4 );
01391     assert( junk2 == 4 );
01392 
01393     // get vertices in quads
01394     iBase_EntityHandle conn[16];
01395     int offsets[5], *off_ptr = offsets, junk3 = 5, junk4;
01396     ptr   = conn;
01397     junk1 = 16;
01398     iMesh_getEntArrAdj( imesh, quads, 4, iBase_VERTEX, &ptr, &junk1, &junk2, &off_ptr, &junk3, &junk4, &ierr );CHKERR;
01399     assert( ptr == conn );
01400     assert( junk1 == 16 );
01401     assert( junk2 == 16 );
01402     assert( off_ptr == offsets );
01403     assert( junk3 == 5 );
01404     assert( junk4 == 5 );
01405 
01406     // make unique vertex list
01407     std::sort( conn, conn + 16 );
01408     const int num_vtx = std::unique( conn, conn + 16 ) - conn;
01409     assert( 9 == num_vtx );
01410 
01411     // get vertex coords
01412     std::vector< double > coords( 27 );
01413     ierr = get_coords( imesh, conn, 9, &coords[0] );CHKERR;
01414 
01415     // use vertex coords to determine logical position
01416     for( int i = 0; i < num_vtx; ++i )
01417     {
01418         int x = (int)round( coords[3 * i] - xbase );
01419         int y = (int)round( coords[3 * i + 1] - ybase );
01420         if( x < 0 || x > 2 || y < 0 || y > 2 )
01421         {
01422             std::cerr << "Processor " << rank << ": Error at " __FILE__ ":" << __LINE__ << std::endl
01423                       << "  Invalid vertex coordinate: (" << coords[3 * i] << ", " << coords[3 * i + 1] << ", "
01424                       << coords[3 * i + 2] << ")" << std::endl
01425                       << "  For logical partition " << local_id << std::endl;
01426             return iBase_FAILURE;
01427         }
01428         verts[x][y] = conn[i];
01429     }
01430 
01431     if( local_logical % 2 )
01432     {
01433         switch( other_logical - local_logical )
01434         {
01435             case 0:
01436                 return iBase_FAILURE;
01437             case 1:  // upper right
01438                 vtx_handles.resize( 1 );
01439                 vtx_handles[0] = verts[2][0];
01440                 break;
01441             case 2:  // right
01442                 vtx_handles.resize( 3 );
01443                 std::copy( verts[2], verts[2] + 3, vtx_handles.begin() );
01444                 break;
01445             case -1:  // above
01446                 vtx_handles.resize( 3 );
01447                 vtx_handles[0] = verts[0][0];
01448                 vtx_handles[1] = verts[1][0];
01449                 vtx_handles[2] = verts[2][0];
01450                 break;
01451             case -2:  // left
01452                 vtx_handles.resize( 3 );
01453                 std::copy( verts[0], verts[0] + 3, vtx_handles.begin() );
01454                 break;
01455             case -3:  // upper left
01456                 vtx_handles.resize( 1 );
01457                 vtx_handles[0] = verts[0][0];
01458                 break;
01459             default:
01460                 vtx_handles.clear();
01461                 break;
01462         }
01463     }
01464     else
01465     {
01466         switch( other_logical - local_logical )
01467         {
01468             case 0:
01469                 return iBase_FAILURE;
01470             case 1:  // below
01471                 vtx_handles.resize( 3 );
01472                 vtx_handles[0] = verts[0][2];
01473                 vtx_handles[1] = verts[1][2];
01474                 vtx_handles[2] = verts[2][2];
01475                 break;
01476             case 2:  // right
01477                 vtx_handles.resize( 3 );
01478                 std::copy( verts[2], verts[2] + 3, vtx_handles.begin() );
01479                 break;
01480             case 3:  // lower right
01481                 vtx_handles.resize( 1 );
01482                 vtx_handles[0] = verts[2][2];
01483                 break;
01484             case -1:  // lower left
01485                 vtx_handles.resize( 1 );
01486                 vtx_handles[0] = verts[0][2];
01487                 break;
01488             case -2:  // left
01489                 vtx_handles.resize( 3 );
01490                 std::copy( verts[0], verts[0] + 3, vtx_handles.begin() );
01491                 break;
01492             default:
01493                 vtx_handles.clear();
01494                 break;
01495         }
01496     }
01497 
01498     return iBase_SUCCESS;
01499 }
01500 
01501 /**\brief Test querying of part boundary entities
01502  *
01503  * Test:
01504  * - iMeshP_getNumPartBdryEnts
01505  * - iMeshP_getPartBdryEnts
01506  */
01507 int test_get_part_boundary( iMesh_Instance imesh, iMeshP_PartitionHandle prtn, const PartMap& map )
01508 {
01509     int ierr, rank;
01510     MPI_Comm_rank( MPI_COMM_WORLD, &rank );
01511 
01512     // get local part handles and part ids, and global part id list
01513     std::vector< iMeshP_PartHandle > local_handles;
01514     std::vector< iMeshP_Part > local_ids;
01515     std::vector< iMeshP_Part > all_parts = map.get_parts();
01516     std::map< iMeshP_PartHandle, std::vector< iBase_EntityHandle > > part_bdry;
01517     ierr = get_local_parts( imesh, prtn, local_handles, &local_ids );CHKERR;
01518 
01519     // for each combination of local part with any other part,
01520     // check for valid function values.
01521     std::vector< std::pair< iMeshP_Part, iMeshP_Part > > num_failed, num_error, list_failed, list_error, error;
01522     for( size_t i = 0; i < local_handles.size(); ++i )
01523     {
01524         iMeshP_PartHandle local_handle = local_handles[i];
01525         iMeshP_Part local_id           = local_ids[i];
01526         for( std::vector< iMeshP_Part >::iterator j = all_parts.begin(); j != all_parts.end(); ++j )
01527         {
01528             iMeshP_Part other_id = *j;
01529             if( other_id == local_id ) continue;
01530 
01531             std::pair< iMeshP_Part, iMeshP_Part > part_pair;
01532             part_pair.first  = local_id;
01533             part_pair.second = other_id;
01534 
01535             // get expected values
01536             std::vector< iBase_EntityHandle > shared_verts;
01537             ierr = interface_verts( imesh, prtn, local_handle, other_id, map, shared_verts );
01538             if( ierr != iBase_SUCCESS )
01539             {
01540                 error.push_back( part_pair );
01541                 continue;
01542             }
01543             std::sort( shared_verts.begin(), shared_verts.end() );
01544 
01545             // test iMeshP_getNumPartBdryEnts
01546             int count;
01547             iMeshP_getNumPartBdryEnts( imesh, prtn, local_handle, iBase_VERTEX, iMesh_POINT, other_id, &count, &ierr );
01548             if( iBase_SUCCESS != ierr )
01549                 num_error.push_back( part_pair );
01550             else if( count != (int)shared_verts.size() )
01551                 num_failed.push_back( part_pair );
01552 
01553             // test iMeshP_getPartBdryEnts
01554             iBase_EntityHandle* ptr = 0;
01555             int junk                = 0;
01556             iMeshP_getPartBdryEnts( imesh, prtn, local_handle, iBase_VERTEX, iMesh_POINT, other_id, &ptr, &junk, &count,
01557                                     &ierr );
01558             if( iBase_SUCCESS != ierr )
01559                 list_error.push_back( part_pair );
01560             else
01561             {
01562                 std::copy( ptr, ptr + count, std::back_inserter( part_bdry[local_handles[i]] ) );
01563                 std::sort( ptr, ptr + count );
01564                 if( (int)shared_verts.size() != count || !std::equal( shared_verts.begin(), shared_verts.end(), ptr ) )
01565                     list_failed.push_back( part_pair );
01566                 free( ptr );
01567             }
01568         }
01569     }
01570 
01571     if( !error.empty() )
01572     {
01573         std::cerr << "Processor " << rank << ": Error at " __FILE__ ":" << __LINE__ << std::endl
01574                   << "  Internal error for " << error.size() << " part pairs." << std::endl;
01575         ierr = iBase_FAILURE;
01576     }
01577     if( !num_error.empty() )
01578     {
01579         std::cerr << "Processor " << rank << ": Error at " __FILE__ ":" << __LINE__ << std::endl
01580                   << "  iMeshP_getNumPartBdryEnts return error for " << num_error.size() << " part pairs." << std::endl;
01581         ierr = iBase_FAILURE;
01582     }
01583     if( !list_error.empty() )
01584     {
01585         std::cerr << "Processor " << rank << ": Error at " __FILE__ ":" << __LINE__ << std::endl
01586                   << "  iMeshP_getPartBdryEnts return error for " << list_error.size() << " part pairs." << std::endl;
01587         ierr = iBase_FAILURE;
01588     }
01589     if( !num_failed.empty() )
01590     {
01591         std::cerr << "Processor " << rank << ": Error at " __FILE__ ":" << __LINE__ << std::endl
01592                   << "  iMeshP_getNumPartBdryEnts return incorrect results for " << num_failed.size() << " part pairs."
01593                   << std::endl;
01594         ierr = iBase_FAILURE;
01595     }
01596     if( !list_failed.empty() )
01597     {
01598         std::cerr << "Processor " << rank << ": Error at " __FILE__ ":" << __LINE__ << std::endl
01599                   << "  iMeshP_getPartBdryEnts return incorrect results for " << list_failed.size() << " part pairs."
01600                   << std::endl;
01601         ierr = iBase_FAILURE;
01602     }
01603 
01604     if( iBase_SUCCESS != ierr ) return ierr;
01605 
01606     // test with iMeshP_ALL_PARTS
01607     for( size_t i = 0; i < local_handles.size(); ++i )
01608     {
01609         std::vector< iBase_EntityHandle >& exp_bdry = part_bdry[local_handles[i]];
01610         std::sort( exp_bdry.begin(), exp_bdry.end() );
01611         exp_bdry.erase( std::unique( exp_bdry.begin(), exp_bdry.end() ), exp_bdry.end() );
01612         std::pair< iMeshP_Part, iMeshP_Part > part_pair;
01613         part_pair.first  = local_ids[i];
01614         part_pair.second = iMeshP_ALL_PARTS;
01615 
01616         int num = 0;
01617         iMeshP_getNumPartBdryEnts( imesh, prtn, local_handles[i], iBase_VERTEX, iMesh_POINT, iMeshP_ALL_PARTS, &num,
01618                                    &ierr );
01619         if( ierr )
01620             num_error.push_back( part_pair );
01621         else if( num != (int)exp_bdry.size() )
01622             num_failed.push_back( part_pair );
01623 
01624         iBase_EntityHandle* bdry = 0;
01625         int junk = num = 0;
01626         iMeshP_getPartBdryEnts( imesh, prtn, local_handles[i], iBase_VERTEX, iMesh_POINT, iMeshP_ALL_PARTS, &bdry,
01627                                 &junk, &num, &ierr );
01628         if( ierr )
01629             list_error.push_back( part_pair );
01630         else
01631         {
01632             std::sort( bdry, bdry + num );
01633             if( num != (int)exp_bdry.size() || !std::equal( bdry, bdry + num, exp_bdry.begin() ) )
01634                 list_failed.push_back( part_pair );
01635             free( bdry );
01636         }
01637     }
01638     if( !num_error.empty() )
01639     {
01640         std::cerr << "Processor " << rank << ": Error at " __FILE__ ":" << __LINE__ << std::endl
01641                   << "  iMeshP_getNumPartBdryEnts return error for " << num_error.size() << " part pairs." << std::endl;
01642         ierr = iBase_FAILURE;
01643     }
01644     if( !list_error.empty() )
01645     {
01646         std::cerr << "Processor " << rank << ": Error at " __FILE__ ":" << __LINE__ << std::endl
01647                   << "  iMeshP_getPartBdryEnts return error for " << list_error.size() << " part pairs." << std::endl;
01648         ierr = iBase_FAILURE;
01649     }
01650     if( !num_failed.empty() )
01651     {
01652         std::cerr << "Processor " << rank << ": Error at " __FILE__ ":" << __LINE__ << std::endl
01653                   << "  iMeshP_getNumPartBdryEnts return incorrect results for " << num_failed.size() << " part pairs."
01654                   << std::endl;
01655         ierr = iBase_FAILURE;
01656     }
01657     if( !list_failed.empty() )
01658     {
01659         std::cerr << "Processor " << rank << ": Error at " __FILE__ ":" << __LINE__ << std::endl
01660                   << "  iMeshP_getPartBdryEnts return incorrect results for " << list_failed.size() << " part pairs."
01661                   << std::endl;
01662         ierr = iBase_FAILURE;
01663     }
01664 
01665     return ierr;
01666 }
01667 
01668 /**\brief Test querying of part boundary entities
01669  *
01670  * Test:
01671  * - iMeshP_initPartBdryEntIter
01672  * - iMeshP_initPartBdryEntArrIter
01673  */
01674 int test_part_boundary_iter( iMesh_Instance imesh, iMeshP_PartitionHandle prtn, const PartMap& map )
01675 {
01676     int ierr, rank, has_data;
01677     MPI_Comm_rank( MPI_COMM_WORLD, &rank );
01678 
01679     // get local part handles and part ids, and global part id list
01680     std::vector< iMeshP_PartHandle > local_handles;
01681     std::vector< iMeshP_Part > local_ids;
01682     std::vector< iMeshP_Part > all_parts = map.get_parts();
01683     ierr                                 = get_local_parts( imesh, prtn, local_handles, &local_ids );CHKERR;
01684 
01685     std::vector< std::pair< iMeshP_Part, iMeshP_Part > > single_failed, single_error, single_step_error, array_failed,
01686         array_error, array_step_error;
01687     for( size_t i = 0; i < local_handles.size(); ++i )
01688     {
01689         iMeshP_PartHandle local_handle = local_handles[i];
01690         iMeshP_Part local_id           = local_ids[i];
01691         for( std::vector< iMeshP_Part >::iterator j = all_parts.begin(); j != all_parts.end(); ++j )
01692         {
01693             iMeshP_Part other_id = *j;
01694             if( other_id == local_id ) continue;
01695 
01696             std::pair< iMeshP_Part, iMeshP_Part > part_pair;
01697             part_pair.first  = local_id;
01698             part_pair.second = other_id;
01699 
01700             // get expected values
01701             std::vector< iBase_EntityHandle > shared_verts;
01702             ierr = interface_verts( imesh, prtn, local_handle, other_id, map, shared_verts );
01703             if( ierr != iBase_SUCCESS || 0 == shared_verts.size() ) continue;
01704             std::sort( shared_verts.begin(), shared_verts.end() );
01705 
01706             // test single entity iterator
01707             iBase_EntityIterator siter;
01708             iMeshP_initPartBdryEntIter( imesh, prtn, local_handle, iBase_VERTEX, iMesh_POINT, other_id, &siter, &ierr );
01709             if( ierr != iBase_SUCCESS ) { single_error.push_back( part_pair ); }
01710             else
01711             {
01712                 std::vector< iBase_EntityHandle > results;
01713                 for( ;; )
01714                 {
01715                     iBase_EntityHandle handle;
01716                     iMesh_getNextEntIter( imesh, siter, &handle, &has_data, &ierr );
01717                     if( ierr != iBase_SUCCESS )
01718                     {
01719                         single_step_error.push_back( part_pair );
01720                         break;
01721                     }
01722                     if( !has_data ) break;
01723                     results.push_back( handle );
01724                 }
01725 
01726                 std::sort( results.begin(), results.end() );
01727                 if( results.size() != shared_verts.size() ||
01728                     !std::equal( results.begin(), results.end(), shared_verts.begin() ) )
01729                     single_failed.push_back( part_pair );
01730             }
01731             iMesh_endEntIter( imesh, siter, &ierr );
01732 
01733             // test array iterator
01734             iBase_EntityArrIterator aiter;
01735             iMeshP_initPartBdryEntArrIter( imesh, prtn, local_handle, iBase_VERTEX, iMesh_POINT, shared_verts.size(),
01736                                            other_id, &aiter, &ierr );
01737             if( ierr != iBase_SUCCESS )
01738             {
01739                 array_error.push_back( part_pair );
01740                 iMesh_endEntArrIter( imesh, aiter, &ierr );
01741                 continue;
01742             }
01743             iBase_EntityHandle results[5], *ptr = results;
01744             int junk = 5, count;
01745             iMesh_getNextEntArrIter( imesh, aiter, &ptr, &junk, &count, &has_data, &ierr );
01746             if( ierr != iBase_SUCCESS || !has_data )
01747             {
01748                 array_step_error.push_back( part_pair );
01749                 iMesh_endEntArrIter( imesh, aiter, &ierr );
01750                 continue;
01751             }
01752             iMesh_endEntArrIter( imesh, aiter, &ierr );
01753             assert( count <= 5 );
01754             assert( ptr == results );
01755             std::sort( ptr, ptr + count );
01756             if( count != (int)shared_verts.size() || !std::equal( shared_verts.begin(), shared_verts.end(), results ) )
01757                 array_failed.push_back( part_pair );
01758         }
01759     }
01760 
01761     if( !single_error.empty() )
01762     {
01763         std::cerr << "Processor " << rank << ": Error at " __FILE__ ":" << __LINE__ << std::endl
01764                   << "  iMeshP_initPartBdryEntIter return error for " << single_error.size() << " part pairs."
01765                   << std::endl;
01766         ierr = iBase_FAILURE;
01767     }
01768     if( !single_step_error.empty() )
01769     {
01770         std::cerr << "Processor " << rank << ": Error at " __FILE__ ":" << __LINE__ << std::endl
01771                   << "  iMesh_getNextEntIter return error for " << single_step_error.size() << " part pairs."
01772                   << std::endl;
01773         ierr = iBase_FAILURE;
01774     }
01775     if( !single_failed.empty() )
01776     {
01777         std::cerr << "Processor " << rank << ": Error at " __FILE__ ":" << __LINE__ << std::endl
01778                   << "  iMeshP_initPartBdryEntIter iterator iterated over invalid entities for " << single_failed.size()
01779                   << " part pairs." << std::endl;
01780         ierr = iBase_FAILURE;
01781     }
01782 
01783     if( !array_error.empty() )
01784     {
01785         std::cerr << "Processor " << rank << ": Error at " __FILE__ ":" << __LINE__ << std::endl
01786                   << "  iMeshP_initPartBdryEntArrIter return error for " << array_error.size() << " part pairs."
01787                   << std::endl;
01788         ierr = iBase_FAILURE;
01789     }
01790     if( !array_step_error.empty() )
01791     {
01792         std::cerr << "Processor " << rank << ": Error at " __FILE__ ":" << __LINE__ << std::endl
01793                   << "  iMesh_getNextEntArrIter return error for " << array_step_error.size() << " part pairs."
01794                   << std::endl;
01795         ierr = iBase_FAILURE;
01796     }
01797     if( !array_failed.empty() )
01798     {
01799         std::cerr << "Processor " << rank << ": Error at " __FILE__ ":" << __LINE__ << std::endl
01800                   << "  iMeshP_initPartBdryEntArrIter iterator iterated over invalid entities for "
01801                   << array_failed.size() << " part pairs." << std::endl;
01802         ierr = iBase_FAILURE;
01803     }
01804 
01805     return ierr;
01806 }
01807 
01808 /**\brief Test adjacent entity query
01809  *
01810  * Test:
01811  * - iMeshP_getAdjEntities
01812  */
01813 int test_get_adjacencies( iMesh_Instance /* imesh */, iMeshP_PartitionHandle /* prtn */, const PartMap& )
01814 {
01815     return iBase_SUCCESS;
01816 }
01817 
01818 /**\brief Test entity iterators
01819  *
01820  * Test:
01821  * - iMeshP_initEntIter
01822  * - iMeshP_initEntArrIter
01823  */
01824 int test_entity_iterator( iMesh_Instance /*imesh */, iMeshP_PartitionHandle /*prtn*/, const PartMap& )
01825 {
01826     return iBase_SUCCESS;
01827 }
01828 
01829 /**\brief Test entity owner queries
01830  *
01831  * Test:
01832  * - iMeshP_getEntOwnerPart
01833  * - iMeshP_getEntOwnerPartArr
01834  * - iMeshP_isEntOwner
01835  * - iMeshP_isEntOwnerArr
01836  */
01837 int test_entity_owner( iMesh_Instance imesh, iMeshP_PartitionHandle prtn, const PartMap& /* map */ )
01838 {
01839     int ierr, rank, size;
01840     MPI_Comm_rank( MPI_COMM_WORLD, &rank );
01841     MPI_Comm_size( MPI_COMM_WORLD, &size );
01842 
01843     // get local part handles and part ids
01844     std::vector< iMeshP_PartHandle > local_handles;
01845     std::vector< iMeshP_Part > local_ids;
01846     ierr = get_local_parts( imesh, prtn, local_handles, &local_ids );
01847     PCHECK;
01848 
01849     // test iMeshP_getEntOwnerPart for quads in each part
01850     std::vector< iBase_EntityHandle > all_quads;
01851     std::vector< iMeshP_Part > quad_owners;
01852     int invalid_count = 0;
01853     for( size_t i = 0; i < local_handles.size(); ++i )
01854     {
01855         std::vector< iBase_EntityHandle > quads;
01856         ierr = get_entities( imesh, local_handles[0], iBase_FACE, iMesh_QUADRILATERAL, quads );
01857         if( ierr ) break;
01858 
01859         for( size_t j = 0; j < quads.size(); ++j )
01860         {
01861             all_quads.push_back( quads[j] );
01862             quad_owners.push_back( local_ids[i] );
01863             iMeshP_Part owner;
01864             iMeshP_getEntOwnerPart( imesh, prtn, quads[j], &owner, &ierr );
01865             if( iBase_SUCCESS != ierr ) break;
01866 
01867             if( owner != local_ids[i] ) ++invalid_count;
01868         }
01869         if( iBase_SUCCESS != ierr ) break;
01870     }
01871     PCHECK;
01872     ASSERT( 0 == invalid_count );
01873 
01874     // test iMeshP_getEntOwnerPartArr for quads in each part
01875     invalid_count = 0;
01876     for( size_t i = 0; i < local_handles.size(); ++i )
01877     {
01878         std::vector< iBase_EntityHandle > quads;
01879         ierr = get_entities( imesh, local_handles[0], iBase_FACE, iMesh_QUADRILATERAL, quads );
01880         if( ierr ) break;
01881 
01882         std::vector< iMeshP_Part > owners( quads.size() ), expected( quads.size(), local_ids[i] );
01883         int junk         = owners.size(), count;
01884         iMeshP_Part* ptr = &owners[0];
01885         iMeshP_getEntOwnerPartArr( imesh, prtn, &quads[0], quads.size(), &ptr, &junk, &count, &ierr );
01886         if( ierr ) break;
01887         assert( ptr == &owners[0] );
01888         assert( junk == (int)owners.size() );
01889         assert( count == (int)quads.size() );
01890         if( owners != expected ) ++invalid_count;
01891     }
01892     PCHECK;
01893     ASSERT( 0 == invalid_count );
01894 
01895     // get all vertices
01896     iBase_EntityHandle* vtx_arr = 0;
01897     int junk1                   = 0, num_vtx;
01898     int *junk2 = 0, junk3 = 0, junk4;
01899     iMesh_getEntArrAdj( imesh, &all_quads[0], all_quads.size(), iBase_VERTEX, &vtx_arr, &junk1, &num_vtx, &junk2,
01900                         &junk3, &junk4, &ierr );
01901     PCHECK;
01902     free( junk2 );
01903     std::sort( vtx_arr, vtx_arr + num_vtx );
01904     num_vtx = std::unique( vtx_arr, vtx_arr + num_vtx ) - vtx_arr;
01905     std::vector< iBase_EntityHandle > all_verts( vtx_arr, vtx_arr + num_vtx );
01906     free( vtx_arr );
01907 
01908     // check consistency between iMeshP_getEntOwnerPart and iMeshP_getEntOwnerPartArr
01909     // for all vertices
01910     std::vector< iMeshP_Part > vert_owners( all_verts.size() );
01911     junk1              = vert_owners.size();
01912     iMeshP_Part* junk5 = &vert_owners[0];
01913     iMeshP_getEntOwnerPartArr( imesh, prtn, &all_verts[0], all_verts.size(), &junk5, &junk1, &junk3, &ierr );
01914     PCHECK;
01915     assert( junk5 == &vert_owners[0] );
01916     assert( junk1 == (int)vert_owners.size() );
01917     assert( junk3 == (int)all_verts.size() );
01918 
01919     invalid_count = 0;
01920     for( size_t i = 0; i < all_verts.size(); ++i )
01921     {
01922         iMeshP_Part owner;
01923         iMeshP_getEntOwnerPart( imesh, prtn, all_verts[i], &owner, &ierr );
01924         if( iBase_SUCCESS != ierr || owner != vert_owners[i] ) ++invalid_count;
01925     }
01926     ASSERT( 0 == invalid_count );
01927 
01928     // get lists for all entities
01929     std::vector< iBase_EntityHandle > all_entities( all_verts );
01930     std::copy( all_quads.begin(), all_quads.end(), std::back_inserter( all_entities ) );
01931     std::vector< iMeshP_Part > all_owners( vert_owners );
01932     std::copy( quad_owners.begin(), quad_owners.end(), std::back_inserter( all_owners ) );
01933 
01934     // check consistency of iMeshP_isEntOwner for all entities
01935     invalid_count = 0;
01936     ierr          = iBase_SUCCESS;
01937     for( size_t i = 0; i < local_handles.size(); ++i )
01938     {
01939         for( size_t j = 0; ierr == iBase_SUCCESS && j < all_entities.size(); ++j )
01940         {
01941             int is_owner;
01942             iMeshP_isEntOwner( imesh, prtn, local_handles[i], all_entities[j], &is_owner, &ierr );
01943             if( ierr != iBase_SUCCESS ) break;
01944             if( !is_owner == ( local_ids[i] == all_owners[j] ) ) ++invalid_count;
01945         }
01946     }
01947     PCHECK;
01948     ASSERT( 0 == invalid_count );
01949 
01950     // check consistency of iMeshP_isEntOwnerArr for all entities
01951     for( size_t i = 0; i < local_handles.size(); ++i )
01952     {
01953         std::vector< int > is_owner_list( all_entities.size() );
01954         junk1      = is_owner_list.size();
01955         int* junk6 = &is_owner_list[0];
01956         iMeshP_isEntOwnerArr( imesh, prtn, local_handles[i], &all_entities[0], all_entities.size(), &junk6, &junk1,
01957                               &junk3, &ierr );
01958         if( iBase_SUCCESS != ierr ) break;
01959         assert( junk6 == &is_owner_list[0] );
01960         assert( junk1 == (int)is_owner_list.size() );
01961         assert( junk3 == (int)all_entities.size() );
01962         invalid_count = 0;
01963         for( size_t j = 0; j < all_entities.size(); ++j )
01964         {
01965             if( !( is_owner_list[j] ) == ( local_ids[0] == all_owners[j] ) ) ++invalid_count;
01966         }
01967     }
01968     PCHECK;
01969     ASSERT( 0 == invalid_count );
01970 
01971     // check globally consistent owners for all vertices
01972 
01973     // first communicate total number of vertex entries to be sent to root proc
01974     int local_count = all_verts.size(), global_count = 0;
01975     ierr = MPI_Reduce( &local_count, &global_count, 1, MPI_INT, MPI_SUM, 0, MPI_COMM_WORLD );CHKERR;
01976 
01977     // for each vertex, store { (x << 2) | y, owning part id }
01978     std::vector< int > vtxdata( 2 * all_verts.size() );
01979     std::vector< double > coords( 3 * all_verts.size() );
01980     ierr = get_coords( imesh, &all_verts[0], all_verts.size(), &coords[0] );CHKERR;
01981     for( size_t i = 0; i < all_verts.size(); ++i )
01982     {
01983         int x              = (int)round( coords[3 * i] );
01984         int y              = (int)round( coords[3 * i + 1] );
01985         vtxdata[2 * i]     = ( x << 3 ) | y;
01986         vtxdata[2 * i + 1] = vert_owners[i];
01987     }
01988 
01989     // collect all data on root procesor
01990     std::vector< int > all_data( 2 * global_count );
01991     std::vector< int > displ( size ), counts( size );
01992     for( int i = 0; i < size; i++ )
01993     {
01994         counts[i] = vtxdata.size();
01995         displ[i]  = i * vtxdata.size();
01996     }
01997 
01998     // we could have used a simple gather, because all sequences are the same
01999     ierr = MPI_Gatherv( &vtxdata[0], vtxdata.size(), MPI_INT, &all_data[0], &counts[0], &displ[0], MPI_INT, 0,
02000                         MPI_COMM_WORLD );CHKERR;
02001 
02002     if( rank == 0 )
02003     {
02004         // map from vertex tag to indices into data
02005         std::multimap< int, int > data_map;
02006         for( int i = 0; i < global_count; ++i )
02007         {
02008             std::pair< int, int > p;
02009             p.first  = all_data[2 * i];
02010             p.second = i;
02011             data_map.insert( p );
02012         }
02013 
02014         // check consistent data for each vtx
02015         std::multimap< int, int >::const_iterator a, b;
02016         for( a = data_map.begin(); a != data_map.end(); a = b )
02017         {
02018             for( b = a; b != data_map.end() && a->first == b->first; ++b )
02019             {
02020                 int idx1 = a->second;
02021                 int idx2 = b->second;
02022                 if( all_data[2 * idx1 + 1] == all_data[2 * idx2 + 1] ) continue;
02023 
02024                 ierr = iBase_FAILURE;
02025 
02026                 int proc1 = std::lower_bound( displ.begin(), displ.end(), 2 * idx1 ) - displ.begin();
02027                 if( displ[proc1] != 2 * idx1 ) ++proc1;
02028                 int proc2 = std::lower_bound( displ.begin(), displ.end(), 2 * idx2 ) - displ.begin();
02029                 if( displ[proc2] != 2 * idx2 ) ++proc2;
02030 
02031                 std::cerr << "Error at " __FILE__ ":" << __LINE__ << " : " << std::endl
02032                           << "  For vertex at (" << ( a->first >> 2 ) << ", " << ( a->first & 3 ) << ") :" << std::endl
02033                           << "  Processor " << proc1 << " has " << all_data[2 * idx1 + 1] << " as the owning part"
02034                           << std::endl
02035                           << "  Processor " << proc2 << " has " << all_data[2 * idx2 + 1] << " as the owning part"
02036                           << std::endl;
02037             }
02038         }
02039     }
02040 
02041     return ierr;
02042 }
02043 
02044 static int get_part_boundary_verts( iMesh_Instance imesh, iMeshP_PartitionHandle prtn, const PartMap& map,
02045                                     iMeshP_PartHandle part, std::vector< iBase_EntityHandle >& boundary )
02046 {
02047     int ierr, logical_id;
02048     ierr = map.part_from_coords( imesh, part, logical_id );CHKERR;
02049 
02050     int neighbors[5], num_neighbors;
02051     get_part_neighbors( logical_id, map.get_parts().size(), neighbors, num_neighbors );
02052 
02053     for( int j = 0; j < num_neighbors; ++j )
02054     {
02055         std::vector< iBase_EntityHandle > iface;
02056         ierr = interface_verts( imesh, prtn, part, neighbors[j], map, iface );CHKERR;
02057         std::copy( iface.begin(), iface.end(), std::back_inserter( boundary ) );
02058     }
02059 
02060     std::sort( boundary.begin(), boundary.end() );
02061     boundary.erase( std::unique( boundary.begin(), boundary.end() ), boundary.end() );
02062     return iBase_SUCCESS;
02063 }
02064 
02065 /**\brief Test entity status
02066  *
02067  * Test:
02068  * - iMeshP_getEntStatus
02069  * - iMeshP_getEntStatusArr
02070  */
02071 int test_entity_status( iMesh_Instance imesh, iMeshP_PartitionHandle prtn, const PartMap& map )
02072 {
02073     int ierr, rank, size;
02074     MPI_Comm_rank( MPI_COMM_WORLD, &rank );
02075     MPI_Comm_size( MPI_COMM_WORLD, &size );
02076 
02077     // get local part handles
02078     std::vector< iMeshP_PartHandle > parts;
02079     ierr = get_local_parts( imesh, prtn, parts );
02080     PCHECK;
02081 
02082     // for each part
02083     int num_quad_ent_incorrect = 0, num_quad_ent_error = 0;
02084     int num_quad_arr_incorrect = 0, num_quad_arr_error = 0;
02085     int num_vert_ent_incorrect = 0, num_vert_ent_error = 0;
02086     int num_vert_arr_incorrect = 0, num_vert_arr_error = 0;
02087     for( size_t i = 0; i < parts.size(); ++i )
02088     {
02089         const iMeshP_PartHandle part = parts[i];
02090 
02091         // get quads and vertices
02092         std::vector< iBase_EntityHandle > quads, verts;
02093         ierr = get_part_quads_and_verts( imesh, part, quads, verts );
02094         if( ierr ) break;
02095 
02096         // check quad status (no ghosting yet)
02097         for( size_t j = 0; j < quads.size(); ++j )
02098         {
02099             int status;
02100             iMeshP_getEntStatus( imesh, prtn, part, quads[j], &status, &ierr );
02101             if( ierr != iBase_SUCCESS )
02102             {
02103                 ++num_quad_ent_error;
02104                 ierr = iBase_SUCCESS;
02105                 continue;
02106             }
02107 
02108             if( status != iMeshP_INTERNAL ) ++num_quad_ent_incorrect;
02109         }
02110 
02111         // check quad status using iMeshP_getEntStatusArr
02112         std::vector< int > stat_list( quads.size() );
02113         int* junk1 = &stat_list[0];
02114         int junk2  = stat_list.size(), count;
02115         iMeshP_getEntStatusArr( imesh, prtn, part, &quads[0], quads.size(), &junk1, &junk2, &count, &ierr );
02116         if( ierr != iBase_SUCCESS )
02117         {
02118             ++num_quad_arr_error;
02119             ierr = iBase_SUCCESS;
02120             continue;
02121         }
02122         assert( junk1 == &stat_list[0] );
02123         assert( junk2 == (int)stat_list.size() );
02124         assert( count == (int)quads.size() );
02125         for( size_t j = 0; j < quads.size(); ++j )
02126             if( stat_list[j] != iMeshP_INTERNAL ) ++num_quad_arr_incorrect;
02127 
02128         // figure out which vertices are on the boundary
02129         std::vector< iBase_EntityHandle > boundary;
02130         ierr = get_part_boundary_verts( imesh, prtn, map, part, boundary );
02131         if( ierr ) break;
02132         std::sort( boundary.begin(), boundary.end() );
02133 
02134         // check vertex status (no ghosting yet)
02135         for( size_t j = 0; j < verts.size(); ++j )
02136         {
02137             int status;
02138             iMeshP_getEntStatus( imesh, prtn, part, verts[j], &status, &ierr );
02139             if( ierr != iBase_SUCCESS )
02140             {
02141                 ++num_vert_ent_error;
02142                 ierr = iBase_SUCCESS;
02143                 continue;
02144             }
02145             bool on_boundary = std::binary_search( boundary.begin(), boundary.end(), verts[j] );
02146             if( status != ( on_boundary ? iMeshP_BOUNDARY : iMeshP_INTERNAL ) ) ++num_vert_ent_incorrect;
02147         }
02148 
02149         // check vert status using iMeshP_getEntStatusArr
02150         stat_list.resize( verts.size() );
02151         junk1 = &stat_list[0];
02152         junk2 = stat_list.size();
02153         iMeshP_getEntStatusArr( imesh, prtn, part, &verts[0], verts.size(), &junk1, &junk2, &count, &ierr );
02154         if( ierr != iBase_SUCCESS )
02155         {
02156             ++num_vert_arr_error;
02157             ierr = iBase_SUCCESS;
02158             continue;
02159         }
02160         assert( junk1 == &stat_list[0] );
02161         assert( junk2 == (int)stat_list.size() );
02162         assert( count == (int)verts.size() );
02163         for( size_t j = 0; j < verts.size(); ++j )
02164         {
02165             bool on_boundary = std::binary_search( boundary.begin(), boundary.end(), verts[j] );
02166             if( stat_list[j] != ( on_boundary ? iMeshP_BOUNDARY : iMeshP_INTERNAL ) ) ++num_vert_arr_incorrect;
02167         }
02168     }
02169     PCHECK;  // check if loop interrupted by any internal errors
02170 
02171     ASSERT( 0 == num_quad_ent_error );
02172     ASSERT( 0 == num_quad_arr_error );
02173     ASSERT( 0 == num_vert_ent_error );
02174     ASSERT( 0 == num_vert_arr_error );
02175     ASSERT( 0 == num_quad_ent_incorrect );
02176     ASSERT( 0 == num_quad_arr_incorrect );
02177     ASSERT( 0 == num_vert_ent_incorrect );
02178     ASSERT( 0 == num_vert_arr_incorrect );
02179 
02180     return iBase_SUCCESS;
02181 }
02182 
02183 /**\brief Test information about entity copies for interface entities
02184  *
02185  * Test:
02186  * - iMeshP_getNumCopies
02187  * - iMeshP_getCopyParts
02188  */
02189 int test_entity_copy_parts( iMesh_Instance imesh, iMeshP_PartitionHandle prtn, const PartMap& map )
02190 {
02191     int ierr, rank, size;
02192     MPI_Comm_rank( MPI_COMM_WORLD, &rank );
02193     MPI_Comm_size( MPI_COMM_WORLD, &size );
02194 
02195     // get local part handles
02196     std::vector< iMeshP_PartHandle > parts;
02197     ierr = get_local_parts( imesh, prtn, parts );
02198     PCHECK;
02199     ASSERT( !parts.empty() );
02200 
02201     // select a singe part to test
02202     const iMeshP_PartHandle part = parts[0];
02203     int logical_id;
02204     ierr = map.part_from_coords( imesh, part, logical_id );CHKERR;
02205     const iMeshP_Part part_id = map.part_id_from_local_id( logical_id );
02206 
02207     // get vertices in part
02208     std::vector< iBase_EntityHandle > quads, verts;
02209     ierr = get_part_quads_and_verts( imesh, part, quads, verts );
02210     PCHECK;
02211 
02212     // get neighbors
02213     int neighbors[5], num_neighbors;
02214     get_part_neighbors( logical_id, map.get_parts().size(), neighbors, num_neighbors );
02215 
02216     // build map of sharing data for each vertex
02217     std::map< iBase_EntityHandle, std::vector< iMeshP_Part > > vert_sharing;
02218     for( int j = 0; j < num_neighbors; ++j )
02219     {
02220         std::vector< iBase_EntityHandle > iface;
02221         ierr = interface_verts( imesh, prtn, part, neighbors[j], map, iface );CHKERR;
02222         for( size_t k = 0; k < iface.size(); ++k )
02223             vert_sharing[iface[k]].push_back( map.part_id_from_local_id( neighbors[j] ) );
02224     }
02225 
02226     // test getNumCopies for each vertex
02227     std::map< iBase_EntityHandle, std::vector< iMeshP_Part > >::iterator i;
02228     int num_failed = 0, num_incorrect = 0;
02229     for( i = vert_sharing.begin(); i != vert_sharing.end(); ++i )
02230     {
02231         int count;
02232         iBase_EntityHandle vtx = i->first;
02233         iMeshP_getNumCopies( imesh, prtn, vtx, &count, &ierr );
02234         if( ierr )
02235             ++num_failed;
02236         else if( (unsigned)count != i->second.size() + 1 )  // add one for the part we queried from
02237             ++num_incorrect;
02238     }
02239     ASSERT( 0 == num_failed );
02240     ASSERT( 0 == num_incorrect );
02241 
02242     // get getCopyParts for each vertex
02243     num_failed = num_incorrect = 0;
02244     for( i = vert_sharing.begin(); i != vert_sharing.end(); ++i )
02245     {
02246         iMeshP_Part* list = 0;
02247         int junk          = 0, count;
02248         iMeshP_getCopyParts( imesh, prtn, i->first, &list, &junk, &count, &ierr );
02249         if( iBase_SUCCESS != ierr )
02250         {
02251             ++num_failed;
02252             continue;
02253         }
02254         if( (unsigned)count != i->second.size() + 1 )
02255         {  // add one for the part we queried from
02256             ++num_incorrect;
02257             free( list );
02258             continue;
02259         }
02260 
02261         std::vector< iMeshP_Part > expected( i->second );
02262         expected.push_back( part_id );
02263         std::sort( list, list + count );
02264         std::sort( expected.begin(), expected.end() );
02265         bool eq = std::equal( list, list + count, expected.begin() );
02266         free( list );
02267         if( !eq ) ++num_incorrect;
02268     }
02269     ASSERT( 0 == num_failed );
02270     ASSERT( 0 == num_incorrect );
02271 
02272     return iBase_SUCCESS;
02273 }
02274 
02275 // store remote handle data for a vertex
02276 struct VtxCopyData
02277 {
02278     std::vector< iMeshP_Part > parts;
02279     std::vector< iBase_EntityHandle > handles;
02280 };
02281 
02282 /**\brief Test information about entity copies for interface entities
02283  *
02284  * Test:
02285  * - iMeshP_getCopies
02286  * - iMeshP_getCopyOnPart
02287  * - iMeshP_getOwnerCopy
02288  */
02289 int test_entity_copies( iMesh_Instance imesh, iMeshP_PartitionHandle prtn, const PartMap& /* map */ )
02290 {
02291     int ierr, rank, size;
02292     MPI_Comm_rank( MPI_COMM_WORLD, &rank );
02293     MPI_Comm_size( MPI_COMM_WORLD, &size );
02294 
02295     // generate a unique ID for each vertex using the coordinates.
02296     // see create_mesh(..): each vertex has integer coordinates (x,y,0)
02297     //                      with x in [0,inf] and y in [0,4]
02298     // then to an Allgatherv to exchange handles for each processor
02299 
02300     // cast everything to iBase_EntityHandle so we can pack it all in one communication
02301     MPI_Datatype tmp_type;
02302     if( sizeof( iBase_EntityHandle ) == sizeof( unsigned ) )
02303         tmp_type = MPI_UNSIGNED;
02304     else if( sizeof( iBase_EntityHandle ) == sizeof( unsigned long ) )
02305         tmp_type = MPI_UNSIGNED_LONG;
02306     else if( sizeof( iBase_EntityHandle ) == sizeof( unsigned long long ) )
02307         tmp_type = MPI_UNSIGNED_LONG_LONG;
02308     else
02309         return iBase_FAILURE;
02310     const MPI_Datatype type = tmp_type;  // make it const
02311 
02312     // get local part handles
02313     std::vector< iMeshP_PartHandle > parts;
02314     ierr = get_local_parts( imesh, prtn, parts );
02315     PCHECK;
02316     std::vector< iMeshP_Part > part_ids( parts.size() );
02317     iMeshP_Part* junk1 = &part_ids[0];
02318     int junk2          = part_ids.size(), junk3;
02319     iMeshP_getPartIdsFromPartHandlesArr( imesh, prtn, &parts[0], parts.size(), &junk1, &junk2, &junk3, &ierr );
02320     PCHECK;
02321     assert( junk1 == &part_ids[0] );
02322     assert( junk2 == (int)part_ids.size() );
02323     assert( junk3 == (int)parts.size() );
02324 
02325     // build list of {vtx_id, part_id, handle} tuples to send
02326     // also build list of local vertex handles
02327     std::vector< iBase_EntityHandle > local_data, local_vertices;
02328     for( size_t i = 0; i < parts.size(); ++i )
02329     {
02330         // get vertices
02331         std::vector< iBase_EntityHandle > quads, verts;
02332         ierr = get_part_quads_and_verts( imesh, parts[i], quads, verts );
02333         if( ierr ) break;
02334 
02335         // add all vertices to local_data
02336         for( size_t j = 0; j < verts.size(); ++j )
02337         {
02338             int tag = 0;
02339             ierr    = vertex_tag( imesh, verts[j], tag );
02340             if( ierr ) break;
02341             long tmp_h = tag;
02342             local_data.push_back( (iBase_EntityHandle)tmp_h );
02343             tmp_h = part_ids[i];
02344             local_data.push_back( (iBase_EntityHandle)tmp_h );
02345             local_data.push_back( verts[j] );
02346         }
02347         if( ierr ) break;
02348 
02349         std::copy( verts.begin(), verts.end(), std::back_inserter( local_vertices ) );
02350     }
02351 
02352     // build list of local vertices
02353     std::sort( local_vertices.begin(), local_vertices.end() );
02354     local_vertices.erase( std::unique( local_vertices.begin(), local_vertices.end() ), local_vertices.end() );
02355     std::vector< int > local_vtx_tags( local_vertices.size() );CHKERR;
02356     for( size_t i = 0; i < local_vertices.size(); ++i )
02357     {
02358         ierr = vertex_tag( imesh, local_vertices[i], local_vtx_tags[i] );
02359         if( ierr ) break;
02360     }
02361     CHKERR;
02362 
02363     // communicate data
02364     std::vector< int > gcounts( size ), gdisp( size );
02365     int local_data_size = local_data.size();
02366     ierr                = MPI_Allgather( &local_data_size, 1, MPI_INT, &gcounts[0], 1, MPI_INT, MPI_COMM_WORLD );CHKERR;
02367     gdisp[0] = 0;
02368     for( int i = 1; i < size; ++i )
02369         gdisp[i] = gdisp[i - 1] + gcounts[i - 1];
02370     std::vector< iBase_EntityHandle > global_data( gdisp[size - 1] + gcounts[size - 1] );
02371     ierr = MPI_Allgatherv( &local_data[0], local_data_size, type, &global_data[0], &gcounts[0], &gdisp[0], type,
02372                            MPI_COMM_WORLD );CHKERR;
02373 
02374     // arrange global data in a more useful way
02375     std::map< int, VtxCopyData > vtx_sharing;
02376     assert( global_data.size() % 3 == 0 );
02377     for( size_t i = 0; i < global_data.size(); i += 3 )
02378     {
02379         int tag                   = (int)(size_t)global_data[i];
02380         iMeshP_Part part          = ( iMeshP_Part )(size_t)global_data[i + 1];
02381         iBase_EntityHandle handle = global_data[i + 2];
02382         vtx_sharing[tag].parts.push_back( part );
02383         vtx_sharing[tag].handles.push_back( handle );
02384     }
02385 
02386     // test iMeshP_getCopies for each local vertex
02387     int num_error = 0, num_incorrect = 0, junk4;
02388     for( size_t i = 0; i < local_vertices.size(); ++i )
02389     {
02390         int num_copies = -1;
02391         // iMeshP_Part* part_ids = 0;
02392         iMeshP_Part* ptr_part_ids  = 0;  // Use ptr_part_ids to avoid shadowing std::vector<iMeshP_Part> part_ids
02393         iBase_EntityHandle* copies = 0;
02394         junk2 = junk3 = junk4 = 0;
02395         iMeshP_getCopies( imesh, prtn, local_vertices[i], &ptr_part_ids, &junk2, &num_copies, &copies, &junk3, &junk4,
02396                           &ierr );
02397         if( iBase_SUCCESS != ierr )
02398         {
02399             ++num_error;
02400             continue;
02401         }
02402         assert( junk4 == num_copies );
02403 
02404         VtxCopyData& expected = vtx_sharing[local_vtx_tags[i]];
02405         if( num_copies != (int)expected.parts.size() )
02406             ++num_incorrect;
02407         else
02408             for( size_t j = 0; j < expected.parts.size(); ++j )
02409             {
02410                 int idx = std::find( ptr_part_ids, ptr_part_ids + num_copies, expected.parts[j] ) - ptr_part_ids;
02411                 if( idx == num_copies || copies[idx] != expected.handles[j] )
02412                 {
02413                     ++num_incorrect;
02414                     break;
02415                 }
02416             }
02417         free( ptr_part_ids );
02418         free( copies );
02419     }
02420     ASSERT( 0 == num_error );
02421     ASSERT( 0 == num_incorrect );
02422 
02423     // test iMeshP_getCopyOnPart for each local vertex
02424     num_error = num_incorrect = 0;
02425     for( size_t i = 0; i < local_vertices.size(); ++i )
02426     {
02427         VtxCopyData& expected = vtx_sharing[local_vtx_tags[i]];
02428         for( size_t j = 0; j < expected.parts.size(); ++j )
02429         {
02430             iBase_EntityHandle copy;
02431             iMeshP_getCopyOnPart( imesh, prtn, local_vertices[i], expected.parts[j], &copy, &ierr );
02432             if( iBase_SUCCESS != ierr )
02433                 ++num_error;
02434             else if( expected.handles[j] != copy )
02435                 ++num_incorrect;
02436         }
02437     }
02438     ASSERT( 0 == num_error );
02439     ASSERT( 0 == num_incorrect );
02440 
02441     // test iMeshP_getOwnerCopy for each local vertex
02442     num_error = num_incorrect = 0;
02443     for( size_t i = 0; i < local_vertices.size(); ++i )
02444     {
02445         VtxCopyData& expected = vtx_sharing[local_vtx_tags[i]];
02446         iMeshP_Part owner_id  = 0;
02447         iMeshP_getEntOwnerPart( imesh, prtn, local_vertices[i], &owner_id, &ierr );
02448         if( iBase_SUCCESS != ierr ) continue;  // not testing getEntOwnerPart here
02449 
02450         size_t idx = std::find( expected.parts.begin(), expected.parts.end(), owner_id ) - expected.parts.begin();
02451         if( idx == expected.parts.size() ) continue;  // not testing getEntOwnerPart here
02452 
02453         iMeshP_Part owner_id_2  = 0;
02454         iBase_EntityHandle copy = 0;
02455         iMeshP_getOwnerCopy( imesh, prtn, local_vertices[i], &owner_id_2, &copy, &ierr );
02456         if( iBase_SUCCESS != ierr )
02457             ++num_error;
02458         else if( owner_id_2 != owner_id && copy != expected.handles[idx] )
02459             ++num_incorrect;
02460     }
02461     ASSERT( 0 == num_error );
02462     ASSERT( 0 == num_incorrect );
02463 
02464     return iBase_SUCCESS;
02465 }
02466 
02467 int get_num_adj_quads( iMesh_Instance imesh, iBase_EntityHandle vtx, int& num )
02468 {
02469     iBase_EntityHandle* list = 0;
02470     int ierr, junk = 0;
02471     iMesh_getEntAdj( imesh, vtx, iBase_FACE, &list, &junk, &num, &ierr );
02472     if( iBase_SUCCESS == ierr ) free( list );
02473     return ierr;
02474 }
02475 
02476 int get_adj( iMesh_Instance imesh, iBase_EntityHandle ent, int type, std::vector< iBase_EntityHandle >& adj )
02477 {
02478     iBase_EntityHandle* list = 0;
02479     int ierr, num, junk = 0;
02480     iMesh_getEntAdj( imesh, ent, type, &list, &junk, &num, &ierr );
02481     if( iBase_SUCCESS == ierr )
02482     {
02483         std::copy( list, list + num, std::back_inserter( adj ) );
02484         free( list );
02485     }
02486     return ierr;
02487 }
02488 
02489 // assume regular quad mesh
02490 int get_boundary_vertices( iMesh_Instance imesh, std::vector< iBase_EntityHandle >& bdry )
02491 {
02492     int ierr, n;
02493     iBase_EntitySetHandle root;
02494     iMesh_getRootSet( imesh, &root, &ierr );CHKERR;
02495     std::vector< iBase_EntityHandle > all_verts;
02496     ierr = get_entities( imesh, root, iBase_VERTEX, iMesh_POINT, all_verts );CHKERR;
02497     bdry.clear();
02498     for( size_t i = 0; i < all_verts.size(); ++i )
02499     {
02500         ierr = get_num_adj_quads( imesh, all_verts[i], n );CHKERR;
02501         if( n != 4 ) bdry.push_back( all_verts[i] );
02502     }
02503     return iBase_SUCCESS;
02504 }
02505 
02506 int check_one_layer( iMesh_Instance imesh, iBase_EntityHandle vtx,
02507                      const std::vector< iBase_EntityHandle >& sorted_vertices )
02508 {
02509     int ierr;
02510     if( std::binary_search( sorted_vertices.begin(), sorted_vertices.end(), vtx ) ) return iBase_SUCCESS;
02511     std::vector< iBase_EntityHandle > quads, verts;
02512     ierr = get_adj( imesh, vtx, iBase_FACE, quads );CHKERR;
02513     for( size_t i = 0; i < quads.size(); ++i )
02514     {
02515         verts.clear();
02516         ierr = get_adj( imesh, quads[i], iBase_VERTEX, verts );CHKERR;
02517         for( size_t j = 0; j < verts.size(); ++j )
02518         {
02519             if( std::binary_search( sorted_vertices.begin(), sorted_vertices.end(), verts[j] ) ) return iBase_SUCCESS;
02520         }
02521     }
02522 
02523     return iBase_FAILURE;
02524 }
02525 
02526 // get number of adjacent quads to each vertex, both on the local
02527 // processor and in the entire mesh
02528 int get_num_adj_all( iMesh_Instance imesh, const std::vector< iBase_EntityHandle >& verts,
02529                      std::vector< int >& num_local_adj, std::vector< int >& num_all_adj )
02530 {
02531     int ierr, size;
02532     MPI_Comm_size( MPI_COMM_WORLD, &size );
02533 
02534     std::vector< int > vtx_tags( verts.size() );
02535     num_local_adj.resize( verts.size() );
02536     for( size_t i = 0; i < verts.size(); ++i )
02537     {
02538         ierr = get_num_adj_quads( imesh, verts[i], num_local_adj[i] );CHKERR;
02539         ierr = vertex_tag( imesh, verts[i], vtx_tags[i] );CHKERR;
02540     }
02541 
02542     std::vector< int > counts( size ), displ( size );
02543     int num_vtx = verts.size();
02544     ierr        = MPI_Allgather( &num_vtx, 1, MPI_INT, &counts[0], 1, MPI_INT, MPI_COMM_WORLD );CHKERR;
02545     displ[0] = 0;
02546     for( int i = 1; i < size; ++i )
02547         displ[i] = displ[i - 1] + counts[i - 1];
02548     int total = displ[size - 1] + counts[size - 1];
02549     std::vector< int > all_tags( total ), all_adj_counts( total );
02550     ierr = MPI_Allgatherv( &vtx_tags[0], vtx_tags.size(), MPI_INT, &all_tags[0], &counts[0], &displ[0], MPI_INT,
02551                            MPI_COMM_WORLD );CHKERR;
02552     ierr = MPI_Allgatherv( &num_local_adj[0], num_local_adj.size(), MPI_INT, &all_adj_counts[0], &counts[0], &displ[0],
02553                            MPI_INT, MPI_COMM_WORLD );CHKERR;
02554 
02555     num_all_adj.clear();
02556     num_all_adj.resize( total, 0 );
02557     for( int i = 0; i < total; ++i )
02558     {
02559         std::vector< int >::iterator it = std::find( vtx_tags.begin(), vtx_tags.end(), all_tags[i] );
02560         if( it == vtx_tags.end() ) continue;
02561         int idx = it - vtx_tags.begin();
02562         num_all_adj[idx] += all_adj_counts[i];
02563     }
02564 
02565     return iBase_SUCCESS;
02566 }
02567 
02568 /**\brief Test creation of ghost entities
02569  *
02570  * Test:
02571  * - iMeshP_createGhostEntsAll
02572  */
02573 int test_create_ghost_ents( iMesh_Instance imesh, iMeshP_PartitionHandle prtn, const PartMap& /* map */ )
02574 {
02575     int ierr;
02576 
02577     // get boundary vertices
02578     std::vector< iBase_EntityHandle > bdry;
02579     ierr = get_boundary_vertices( imesh, bdry );
02580     PCHECK;
02581     // get counts of adjacent entities
02582     std::vector< int > num_local_adj, num_global_adj;
02583     ierr = get_num_adj_all( imesh, bdry, num_local_adj, num_global_adj );
02584     PCHECK;
02585     // create one layer of ghost entities
02586     iMeshP_createGhostEntsAll( imesh, prtn, iBase_FACE, iBase_VERTEX, 1, 0, &ierr );
02587     PCHECK;
02588     // check that each vertex has the correct number of adjacent entities
02589     int num_incorrect = 0;
02590     for( size_t i = 0; i < bdry.size(); ++i )
02591     {
02592         int n;
02593         ierr = get_num_adj_quads( imesh, bdry[i], n );
02594         if( iBase_SUCCESS != ierr || num_global_adj[i] != n ) ++num_incorrect;
02595     }
02596     ASSERT( 0 == num_incorrect );
02597     // get new the new boundary
02598     std::vector< iBase_EntityHandle > new_bdry;
02599     ierr = get_boundary_vertices( imesh, new_bdry );
02600     PCHECK;
02601     // check that each vertex on the new boundary is separated by
02602     // at most one layer from the old boundary
02603     std::sort( bdry.begin(), bdry.end() );
02604     num_incorrect = 0;
02605     for( size_t i = 0; i < new_bdry.size(); ++i )
02606     {
02607         ierr = check_one_layer( imesh, new_bdry[i], bdry );
02608         if( ierr ) ++num_incorrect;
02609     }
02610     ASSERT( 0 == num_incorrect );
02611     // make another layer of ghost entiites
02612     bdry.swap( new_bdry );
02613     new_bdry.clear();
02614     ierr = get_num_adj_all( imesh, bdry, num_local_adj, num_global_adj );
02615     PCHECK;
02616     iMeshP_createGhostEntsAll( imesh, prtn, iBase_FACE, iBase_VERTEX, 2, 0, &ierr );
02617     PCHECK;
02618     // check that each vertex has the correct number of adjacent entities
02619     num_incorrect = 0;
02620     for( size_t i = 0; i < bdry.size(); ++i )
02621     {
02622         int n;
02623         ierr = get_num_adj_quads( imesh, bdry[i], n );
02624         if( iBase_SUCCESS != ierr || num_global_adj[i] != n ) ++num_incorrect;
02625     }
02626     // check that each vertex on the new boundary is separated by
02627     // at most one layer from the old boundary
02628     std::sort( bdry.begin(), bdry.end() );
02629     num_incorrect = 0;
02630     for( size_t i = 0; i < new_bdry.size(); ++i )
02631     {
02632         ierr = check_one_layer( imesh, new_bdry[i], bdry );
02633         if( ierr ) ++num_incorrect;
02634     }
02635     ASSERT( 0 == num_incorrect );
02636 
02637     return iBase_SUCCESS;
02638 }
02639 
02640 /**\brief Test exchange entities
02641  *
02642  * Test:
02643  * - iMeshP_exchEntArrToPartsAll
02644  */
02645 int test_exchange_ents( iMesh_Instance imesh, iMeshP_PartitionHandle prtn, const PartMap& map )
02646 {
02647     int ierr, rank, size;
02648     int num_err = 0;
02649     iMeshP_RequestHandle request;
02650     MPI_Comm_rank( MPI_COMM_WORLD, &rank );
02651     MPI_Comm_size( MPI_COMM_WORLD, &size );
02652 
02653     std::vector< iBase_EntityHandle > all_elems;
02654     std::vector< iMeshP_Part > all_ids;
02655     std::vector< iBase_EntityHandle > quads;
02656 
02657     // get local part handles and part ids
02658     std::vector< iMeshP_PartHandle > local_handles;
02659     std::vector< iMeshP_Part > local_ids;
02660     ierr = get_local_parts( imesh, prtn, local_handles, &local_ids );
02661     PCHECK;
02662 
02663     // get loacal quads before exchange
02664     quads.clear();
02665     ierr = get_entities( imesh, local_handles[0], iBase_FACE, iMesh_QUADRILATERAL, quads );CHKERR;
02666     int n_quads = quads.size();
02667 
02668     // send all elements in local processor to all other processors
02669     for( size_t i = 0; i < map.get_parts().size(); ++i )
02670     {
02671         if( map.get_parts()[i] == (unsigned int)rank ) continue;  // skip own rank
02672 
02673         for( int j = 0; j < n_quads; j++ )
02674         {
02675             all_elems.push_back( quads[j] );
02676             all_ids.push_back( map.get_parts()[i] );
02677         }
02678     }
02679 
02680     // exchange entities
02681     iMeshP_exchEntArrToPartsAll( imesh, prtn, &all_elems[0], all_elems.size(), &all_ids[0], 0, 0, &request, &ierr );
02682     if( iBase_SUCCESS != ierr ) ++num_err;
02683 
02684     // get local quads after exchange
02685     quads.clear();
02686     ierr = get_entities( imesh, local_handles[0], iBase_FACE, iMesh_QUADRILATERAL, quads );CHKERR;
02687 
02688     // # of elements should be # of quads * # of processors
02689     ASSERT( quads.size() == (unsigned int)n_quads * size );
02690 
02691     ASSERT( 0 == num_err );
02692 
02693     return iBase_SUCCESS;
02694 }
02695 
02696 /**\brief Test commuinication of tag data
02697  *
02698  * Test:
02699  * - iMeshP_pushTags
02700  * - iMeshP_pushTagsEnt
02701  */
02702 int test_push_tag_data_common( iMesh_Instance imesh, iMeshP_PartitionHandle prtn, int num_ghost_layers )
02703 {
02704     const char* src_name = "test_src";
02705     const char* dst_name = "test_dst";
02706     int ierr, rank;
02707     MPI_Comm_rank( MPI_COMM_WORLD, &rank );
02708 
02709     if( num_ghost_layers )
02710     {
02711         iMeshP_createGhostEntsAll( imesh, prtn, iBase_FACE, iBase_VERTEX, num_ghost_layers, 0, &ierr );
02712         PCHECK;
02713     }
02714 
02715     iBase_TagHandle src_tag, dst_tag;
02716     iMesh_createTag( imesh, src_name, 1, iBase_INTEGER, &src_tag, &ierr, strlen( src_name ) );CHKERR;
02717     iMesh_createTag( imesh, dst_name, 1, iBase_INTEGER, &dst_tag, &ierr, strlen( dst_name ) );CHKERR;
02718 
02719     iBase_EntitySetHandle root;
02720     iMesh_getRootSet( imesh, &root, &ierr );CHKERR;
02721 
02722     std::vector< iBase_EntityHandle > verts;
02723     ierr = get_entities( imesh, root, iBase_VERTEX, iMesh_POINT, verts );CHKERR;
02724 
02725     // test iMeshP_pushTags
02726     // each processor writes its rank on all vertices
02727     // after push, each vertex should be tagged with the rank of its owner
02728 
02729     std::vector< int > tag_vals( verts.size(), rank );
02730     iMesh_setIntArrData( imesh, &verts[0], verts.size(), src_tag, &tag_vals[0], tag_vals.size(), &ierr );CHKERR;
02731 
02732     iMeshP_pushTags( imesh, prtn, src_tag, dst_tag, iBase_VERTEX, iMesh_POINT, &ierr );
02733     PCHECK;
02734 
02735     tag_vals.clear();
02736     tag_vals.resize( verts.size(), -1 );
02737     iBase_TagHandle id_tag;
02738     iMesh_getTagHandle( imesh, "GLOBAL_ID", &id_tag, &ierr, strlen( "GLOBAL_ID" ) );
02739     std::vector< int > ids( verts.size() );
02740     int *junk1 = &ids[0], junk2 = ids.size(), junk3;
02741     iMesh_getIntArrData( imesh, &verts[0], verts.size(), id_tag, &junk1, &junk2, &junk3, &ierr );
02742     PCHECK;
02743     int errcount = 0;
02744     for( size_t i = 0; i < verts.size(); ++i )
02745     {
02746         iMesh_getIntData( imesh, verts[i], dst_tag, &tag_vals[i], &ierr );
02747         if( ierr != iBase_SUCCESS )
02748         {
02749             std::cerr << "Rank " << rank << " : getIntData failed for vertex " << ids[i] << std::endl;
02750             std::cerr.flush();
02751             ++errcount;
02752         }
02753     }
02754     ASSERT( 0 == errcount );
02755 
02756     //  int *junk1 = &tag_vals[0], junk2 = tag_vals.size(), junk3;
02757     //  iMesh_getIntArrData( imesh, &verts[0], verts.size(), dst_tag, &junk1, &junk2, &junk3, &ierr
02758     //  ); PCHECK; assert( junk1 == &tag_vals[0] ); assert( junk2 == (int)tag_vals.size() ); assert(
02759     //  junk3 == (int)verts.size() );
02760 
02761     std::vector< int > expected( verts.size() );
02762     std::vector< iMeshP_Part > parts( verts.size() );
02763     iMeshP_Part* junk4 = &parts[0];
02764     junk2              = parts.size();
02765     iMeshP_getEntOwnerPartArr( imesh, prtn, &verts[0], verts.size(), &junk4, &junk2, &junk3, &ierr );
02766     PCHECK;
02767     assert( junk4 == &parts[0] );
02768     assert( junk2 == (int)parts.size() );
02769     assert( junk3 == (int)verts.size() );
02770     junk1 = &expected[0];
02771     junk2 = expected.size();
02772     iMeshP_getRankOfPartArr( imesh, prtn, &parts[0], parts.size(), &junk1, &junk2, &junk3, &ierr );
02773     PCHECK;
02774     assert( junk1 == &expected[0] );
02775     assert( junk2 == (int)expected.size() );
02776     assert( junk3 == (int)parts.size() );
02777 
02778     ASSERT( tag_vals == expected );
02779 
02780     // test iMeshP_pushTagsEnt
02781     // write -1 on all vertices
02782     // For each vertex owned by this processor and shared with more than
02783     // two others, write the rank of the owning processor.
02784 
02785     tag_vals.clear();
02786     tag_vals.resize( verts.size(), -1 );
02787     iMesh_setIntArrData( imesh, &verts[0], verts.size(), src_tag, &tag_vals[0], tag_vals.size(), &ierr );
02788     PCHECK;
02789     tag_vals.resize( verts.size(), -1 );
02790     iMesh_setIntArrData( imesh, &verts[0], verts.size(), dst_tag, &tag_vals[0], tag_vals.size(), &ierr );
02791     PCHECK;
02792 
02793     std::vector< iBase_EntityHandle > some;
02794     for( size_t i = 0; i < verts.size(); ++i )
02795     {
02796         int num;
02797         iMeshP_getNumCopies( imesh, prtn, verts[i], &num, &ierr );
02798         if( iBase_SUCCESS != ierr ) break;
02799         if( num > 2 )
02800             some.push_back( verts[i] );
02801         else
02802             expected[i] = -1;
02803     }
02804 
02805     tag_vals.clear();
02806     tag_vals.resize( some.size(), rank );
02807     iMesh_setIntArrData( imesh, &some[0], some.size(), src_tag, &tag_vals[0], tag_vals.size(), &ierr );
02808     PCHECK;
02809 
02810     iMeshP_pushTagsEnt( imesh, prtn, src_tag, dst_tag, &some[0], some.size(), &ierr );
02811     PCHECK;
02812 
02813     tag_vals.clear();
02814     tag_vals.resize( verts.size(), -1 );
02815     junk1 = &tag_vals[0];
02816     junk2 = tag_vals.size();
02817     iMesh_getIntArrData( imesh, &verts[0], verts.size(), dst_tag, &junk1, &junk2, &junk3, &ierr );CHKERR;
02818     assert( junk1 == &tag_vals[0] );
02819     assert( junk2 == (int)tag_vals.size() );
02820     assert( junk3 == (int)verts.size() );
02821 
02822     ASSERT( tag_vals == expected );
02823     return iBase_SUCCESS;
02824 }
02825 
02826 /**\brief Test commuinication of tag data
02827  *
02828  * Test:
02829  * - iMeshP_pushTags
02830  * - iMeshP_pushTagsEnt
02831  */
02832 int test_push_tag_data_iface( iMesh_Instance imesh, iMeshP_PartitionHandle prtn, const PartMap& )
02833 {
02834     return test_push_tag_data_common( imesh, prtn, 0 );
02835 }
02836 
02837 /**\brief Test commuinication of tag data
02838  *
02839  * Test:
02840  * - iMeshP_pushTags
02841  * - iMeshP_pushTagsEnt
02842  */
02843 int test_push_tag_data_ghost( iMesh_Instance imesh, iMeshP_PartitionHandle prtn, const PartMap& )
02844 {
02845     return test_push_tag_data_common( imesh, prtn, 1 );
02846 }
02847 
02848 /**************************************************************************
02849                           PartMap class
02850  **************************************************************************/
02851 
02852 int PartMap::build_map( iMesh_Instance imesh, iMeshP_PartitionHandle prtn, int num_expected_parts )
02853 {
02854     int ierr, rank, size;
02855     MPI_Comm_rank( MPI_COMM_WORLD, &rank );
02856     MPI_Comm_size( MPI_COMM_WORLD, &size );
02857 
02858     // get local parts
02859     std::vector< iMeshP_PartHandle > local_parts;
02860     std::vector< iMeshP_Part > imesh_ids;
02861     ierr = get_local_parts( imesh, prtn, local_parts, &imesh_ids );CHKERR;
02862 
02863     // get logical ids for local parts
02864     std::vector< int > local_ids( local_parts.size() );
02865     for( size_t i = 0; i < local_parts.size(); ++i )
02866     {
02867         ierr = part_from_coords( imesh, local_parts[i], local_ids[i] );CHKERR;
02868     }
02869 
02870     // get total number of parts
02871     int num_global = 0, num_local = local_parts.size();
02872     ierr = MPI_Allreduce( &num_local, &num_global, 1, MPI_INT, MPI_SUM, MPI_COMM_WORLD );CHKERR;
02873     if( num_global != num_expected_parts )
02874     {
02875         std::cerr << "Invalid/unexpected global part count at " __FILE__ ":" << __LINE__ << " (proc " << rank
02876                   << "): " << std::endl
02877                   << "  Expected: " << num_expected_parts << std::endl
02878                   << "  Actual:   " << num_global << std::endl;
02879         return 1;
02880     }
02881 
02882     // get counts and displacements for Allgatherv calls
02883     std::vector< int > dspls( size ), counts( size );
02884     ierr = MPI_Allgather( &num_local, 1, MPI_INT, &counts[0], 1, MPI_INT, MPI_COMM_WORLD );CHKERR;
02885     dspls[0] = 0;
02886     for( int i = 1; i < size; ++i )
02887         dspls[i] = dspls[i - 1] + counts[i - 1];
02888 
02889     // gather iMeshP_Part list from each processor
02890     std::vector< unsigned > global_part_ids( num_expected_parts );
02891     assert( sizeof( iMeshP_Part ) == sizeof( int ) );
02892     ierr = MPI_Allgatherv( &imesh_ids[0], num_local, MPI_UNSIGNED, &global_part_ids[0], &counts[0], &dspls[0],
02893                            MPI_UNSIGNED, MPI_COMM_WORLD );CHKERR;
02894 
02895     // gather local ids from each processor
02896     std::vector< int > global_id_list( num_expected_parts );
02897     ierr = MPI_Allgatherv( &local_ids[0], num_local, MPI_INT, &global_id_list[0], &counts[0], &dspls[0], MPI_INT,
02898                            MPI_COMM_WORLD );CHKERR;
02899 
02900     // build owner list
02901     std::vector< int > global_owners( num_expected_parts );
02902     for( int i = 0; i < size; ++i )
02903         for( int j = 0; j < counts[i]; ++j )
02904             global_owners[dspls[i] + j] = i;
02905 
02906     // populate member lists
02907     sortedPartList = global_part_ids;
02908     std::sort( sortedPartList.begin(), sortedPartList.end() );
02909     partLocalIds.resize( num_expected_parts );
02910     partRanks.resize( num_expected_parts );
02911     for( int i = 0; i < num_expected_parts; ++i )
02912     {
02913         int idx = std::lower_bound( sortedPartList.begin(), sortedPartList.end(), global_part_ids[i] ) -
02914                   sortedPartList.begin();
02915         partLocalIds[idx] = global_id_list[i];
02916         partRanks[idx]    = global_owners[i];
02917     }
02918 
02919     // do some consistency checking
02920     if( std::unique( sortedPartList.begin(), sortedPartList.end() ) != sortedPartList.end() )
02921     {
02922         if( rank == 0 )
02923         { std::cerr << "ERROR: Duplicate iMeshP_Part values detected at " __FILE__ ":" << __LINE__ << std::endl; }
02924         return 1;
02925     }
02926 
02927     // build revesre local id map and check for duplicates
02928     localIdReverseMap.clear();
02929     localIdReverseMap.resize( num_expected_parts, -1 );
02930     for( int i = 0; i < num_expected_parts; ++i )
02931     {
02932         int idx = partLocalIds[i];
02933         if( localIdReverseMap[idx] != -1 )
02934         {
02935             if( rank == 0 )
02936             {
02937                 std::cerr << "ERROR: Part mesh has been duplicated in multiple parts." << std::endl
02938                           << "  Detected at " __FILE__ ":" << __LINE__ << std::endl
02939                           << "  See PartMap::part_from_coords" << std::endl;
02940             }
02941             return 1;
02942         }
02943         if( idx >= num_expected_parts )
02944         {
02945             if( rank == 0 )
02946             {
02947                 std::cerr << "ERROR: Part mesh invalid/incorrect mesh." << std::endl
02948                           << "  Detected at " __FILE__ ":" << __LINE__ << std::endl
02949                           << "  See PartMap::part_from_coords" << std::endl;
02950             }
02951             return 1;
02952         }
02953 
02954         localIdReverseMap[idx] = i;
02955     }
02956 
02957     return 0;
02958 }
02959 
02960 void PartMap::part_id_from_rank( int rank, std::vector< iMeshP_Part >& parts ) const
02961 {
02962     for( size_t i = 0; i < sortedPartList.size(); ++i )
02963         if( partRanks[i] == rank ) parts.push_back( sortedPartList[i] );
02964 }
02965 
02966 void PartMap::local_id_from_rank( int rank, std::vector< int >& ids ) const
02967 {
02968     for( size_t i = 0; i < sortedPartList.size(); ++i )
02969         if( partRanks[i] == rank ) ids.push_back( partLocalIds[i] );
02970 }
02971 
02972 int PartMap::part_from_coords( iMesh_Instance imesh, iMeshP_PartHandle part, int& id )
02973 {
02974     int ierr, rank;
02975     MPI_Comm_rank( MPI_COMM_WORLD, &rank );
02976 
02977     // get elements
02978     const int num_elem = 4;
02979     iBase_EntityHandle array[num_elem];
02980     iBase_EntityHandle* ptr = array;
02981     int junk1 = num_elem, n = -1;
02982     iMesh_getEntities( imesh, part, iBase_FACE, iMesh_QUADRILATERAL, &ptr, &junk1, &n, &ierr );CHKERR;
02983     assert( ptr == array );
02984     assert( junk1 == num_elem );
02985     if( n != num_elem )
02986     {
02987         std::cerr << "Internal error at " __FILE__ ":" << __LINE__ << " (proc " << rank
02988                   << "): Expected all parts to have " << num_elem << " elements.  Found one with " << n << std::endl;
02989         return 1;
02990     }
02991 
02992     // get vertices
02993     iBase_EntityHandle adj_array[4 * num_elem];
02994     int junk2, junk3, offset_array[5];
02995     ptr       = adj_array;
02996     junk1     = sizeof( adj_array ) / sizeof( adj_array[0] );
02997     junk2     = sizeof( offset_array ) / sizeof( offset_array[0] );
02998     int* ptr2 = offset_array;
02999     iMesh_getEntArrAdj( imesh, array, num_elem, iBase_VERTEX, &ptr, &junk1, &n, &ptr2, &junk2, &junk3, &ierr );CHKERR;
03000     assert( ptr == adj_array );
03001     assert( ptr2 == offset_array );
03002     assert( junk1 == sizeof( adj_array ) / sizeof( adj_array[0] ) );
03003     assert( junk2 == sizeof( offset_array ) / sizeof( offset_array[0] ) );
03004     assert( n == 4 * num_elem );
03005     assert( offset_array[0] == 0 );
03006     for( int i = 1; i < junk3; ++i )
03007         assert( offset_array[i] - offset_array[i - 1] == 4 );
03008 
03009     // find center vertex
03010     iBase_EntityHandle vtx;
03011     bool all_match;
03012     for( int i = 0; i < 4; ++i )
03013     {
03014         vtx       = adj_array[i];
03015         all_match = true;
03016         for( int j = 1; j < 4; ++j )
03017         {
03018             iBase_EntityHandle* mvtx = adj_array + 4 * j;
03019             int k;
03020             for( k = 0; k < 4; ++k )
03021                 if( mvtx[k] == vtx ) break;
03022             if( k == 4 ) all_match = false;
03023         }
03024         if( all_match ) break;
03025     }
03026     assert( all_match );
03027 
03028     // get center vertex coordinates
03029     double x, y, z;
03030     iMesh_getVtxCoord( imesh, vtx, &x, &y, &z, &ierr );CHKERR;
03031     assert( 0.0 == z );
03032     const int xi = ( (int)round( x ) - 1 ) / 2;
03033     const int yi = ( (int)round( y ) - 1 ) / 2;
03034     assert( xi >= 0 );
03035     assert( yi >= 0 );
03036     assert( fabs( x - 2 * xi - 1 ) < 1e-12 );
03037     assert( fabs( y - 2 * yi - 1 ) < 1e-12 );
03038 
03039     id = 2 * xi + yi;
03040     return 0;
03041 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines