BulkDataParallelVerify.cpp

00001 /*------------------------------------------------------------------------*/
00002 /*                 Copyright 2010 Sandia Corporation.                     */
00003 /*  Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive   */
00004 /*  license for use of this work by or on behalf of the U.S. Government.  */
00005 /*  Export of this program may require a license from the                 */
00006 /*  United States Government.                                             */
00007 /*------------------------------------------------------------------------*/
00008 
00009 
00014 #include <set>
00015 #include <stdexcept>
00016 #include <iostream>
00017 #include <sstream>
00018 #include <algorithm>
00019 #include <assert.h>
00020 
00021 #include <stk_util/parallel/ParallelComm.hpp>
00022 #include <stk_util/parallel/ParallelReduce.hpp>
00023 
00024 #include <stk_mesh/base/Ghosting.hpp>
00025 #include <stk_mesh/base/BulkData.hpp>
00026 #include <stk_mesh/base/MetaData.hpp>
00027 #include <stk_mesh/base/FieldData.hpp>
00028 #include <stk_mesh/base/EntityComm.hpp>
00029 #include <stk_mesh/base/Comm.hpp>
00030 
00031 namespace stk {
00032 namespace mesh {
00033 
00034 namespace {
00035 
00036 bool verify_parallel_attributes( BulkData & M , std::ostream & error_log );
00037 
00038 void pack_owned_verify( CommAll & all , const BulkData & mesh );
00039 
00040 bool unpack_not_owned_verify( CommAll & comm_all ,
00041                               const BulkData & mesh ,
00042                               std::ostream & error_log );
00043 
00044 }
00045 
00046 bool comm_mesh_verify_parallel_consistency(
00047   BulkData & M ,
00048   std::ostream & error_log )
00049 {
00050   int result = 1 ;
00051 
00052   // Verify consistency of parallel attributes
00053 
00054   result = verify_parallel_attributes( M , error_log );
00055 
00056   all_reduce( M.parallel() , ReduceMin<1>( & result ) );
00057 
00058   // Verify entities against owner.
00059 
00060   if ( result ) {
00061     CommAll all( M.parallel() );
00062 
00063     pack_owned_verify( all , M );
00064 
00065     all.allocate_buffers( all.parallel_size() / 4 );
00066 
00067     pack_owned_verify( all , M );
00068 
00069     all.communicate();
00070 
00071     result = unpack_not_owned_verify( all , M , error_log );
00072 
00073     all_reduce( M.parallel() , ReduceMin<1>( & result ) );
00074   }
00075 
00076 #if 0
00077   // Verify sharing against parallel index
00078   if ( result ) {
00079     std::vector< stk::parallel::DistributedIndex::KeyProc > index ;
00080 
00081     m_entities_index.query( index );
00082 
00083     std::vector< stk::parallel::DistributedIndex::KeyProc >::iterator
00084       j_index = index.begin();
00085     std::vector< Entity * >::iterator i_entity = m_entity_comm.begin();
00086 
00087     while ( j_index != index.end() ) {
00088       std::vector< stk::parallel::DistributedIndex::KeyProc >::iterator
00089         i_index = j_index ;
00090 
00091       while ( j_index != index.end() && i_index->first == j_index->first ) {
00092         ++j_index ;
00093       }
00094 
00095       const EntityKey key( & i_index->first );
00096 
00097       for ( ; i_entity != m_entity_comm.end() &&
00098               (**i_entity).key() < key ; ++i_entity );
00099 
00100       bool this_result = i_entity != m_entity_comm.end() &&
00101                          (**i_entity).key() == key ;
00102 
00103       if ( this_result ) {
00104         Entity & entity = **i_entity ;
00105         PairIterEntityComm ec = entity.sharing();
00106         for ( ; ! ec.empty() && i_index != j_index &&
00107                 ec->proc == i_index->second ; ++ec , ++i_index );
00108         this_result = ! ec.empty() || i_index != j_index ;
00109       }
00110 
00111       if ( ! this_result ) {
00112         error_log << "P" << m_parallel_rank ;
00113         error_log << ": Parallel index " ;
00114         print_entity_key( error_log , mesh_meta_data(), key );
00115         error_log << " != " ;
00116         if ( i_entity == m_entity_comm.end() ) {
00117           error_log << " NULL" ;
00118         }
00119         else {
00120           print_entity_key( error_log , mesh_meta_data(), (**i_entity).key() );
00121         }
00122         result = 0 ;
00123       }
00124     }
00125     all_reduce( M.parallel() , ReduceMin<1>( & result ) );
00126   }
00127 #endif
00128 
00129   return result == 1 ;
00130 }
00131 
00132 namespace {
00133 
00134 bool in_owned_closure( const Entity & entity , unsigned p_local )
00135 {
00136   if ( entity.owner_rank() == p_local ) { return true ; }
00137 
00138   const unsigned erank = entity.entity_rank();
00139 
00140   for ( PairIterRelation rel = entity.relations(); ! rel.empty() ; ++rel ) {
00141     if ( erank < rel->entity_rank() &&
00142          p_local == rel->entity()->owner_rank() ) {
00143       return true ;
00144     }
00145   }
00146   return false ;
00147 }
00148 
00149 bool ordered_comm( const Entity & entity )
00150 {
00151   const PairIterEntityComm ec = entity.comm();
00152   const size_t n = ec.size();
00153   for ( size_t i = 1 ; i < n ; ++i ) {
00154     if ( ! ( ec[i-1] < ec[i] ) ) {
00155       return false ;
00156     }
00157   }
00158   return true ;
00159 }
00160 
00161 bool verify_parallel_attributes( BulkData & M , std::ostream & error_log )
00162 {
00163   bool result = true ;
00164 
00165   const MetaData & S = M.mesh_meta_data();
00166   Part & uses_part = S.locally_used_part();
00167   Part & owns_part = S.locally_owned_part();
00168 
00169   const unsigned p_rank = M.parallel_rank();
00170 
00171   const size_t EntityTypeEnd = M.mesh_meta_data().entity_type_count();
00172 
00173   size_t comm_count = 0 ;
00174 
00175   for ( size_t itype = 0 ; itype < EntityTypeEnd ; ++itype ) {
00176     const std::vector< Bucket * > & all_buckets = M.buckets( itype );
00177 
00178     const std::vector<Bucket*>::const_iterator i_end = all_buckets.end();
00179           std::vector<Bucket*>::const_iterator i     = all_buckets.begin();
00180 
00181     while ( i != i_end ) {
00182       Bucket & bucket = **i ; ++i ;
00183 
00184       const bool uses = has_superset( bucket , uses_part );
00185       const bool owns = has_superset( bucket , owns_part );
00186 
00187       const Bucket::iterator j_end = bucket.end();
00188             Bucket::iterator j     = bucket.begin();
00189 
00190       while ( j != j_end ) {
00191         Entity & entity = *j ; ++j ;
00192 
00193         bool this_result = true ;
00194 
00195         const unsigned p_owner    = entity.owner_rank();
00196         const bool     ordered    = ordered_comm( entity );
00197         const bool     shares     = in_shared( entity );
00198         const bool     recv_ghost = in_receive_ghost( entity );
00199         const bool     send_ghost = in_send_ghost( entity );
00200         const bool     owned_closure = in_owned_closure( entity , p_rank );
00201 
00202         if ( ! ordered ) { this_result = false ; }
00203 
00204         // Owner consistency:
00205 
00206         if (   owns && p_owner != p_rank ) { this_result = false ; }
00207         if ( ! owns && p_owner == p_rank ) { this_result = false ; }
00208 
00209         // Owns is defined to be a subset of uses
00210 
00211         if ( owns && ! uses ) { this_result = false ; }
00212 
00215         // if ( uses != owned_closure ) { this_result = false ; }
00216 
00217         // Shared is defined to be a subset of uses.
00218 
00219         if ( shares && ! uses ) { this_result = false ; }
00220 
00221         // Must be either uses or recv_ghost but not both.
00222 
00223         if (   uses &&   recv_ghost ) { this_result = false ; }
00224         if ( ! uses && ! recv_ghost ) { this_result = false ; }
00225 
00226         // If sending as a ghost then I must own it
00227 
00228         if ( ! owns && send_ghost ) { this_result = false ; }
00229  
00230         // If shared then I am owner or owner is in the shared list
00231 
00232         if ( shares && p_owner != p_rank ) {
00233           PairIterEntityComm ip = entity.sharing();
00234           for ( ; ! ip.empty() && p_owner != ip->proc ; ++ip );
00235           if ( ip.empty() ) { this_result = false ; }
00236         }
00237 
00238         if ( shares || recv_ghost || send_ghost ) { ++comm_count ; }
00239 
00240         if ( ! this_result ) {
00241           result = false ;
00242           error_log << "P" << M.parallel_rank() << ": " ;
00243           print_entity_key( error_log , M.mesh_meta_data(), entity.key() );
00244           error_log << " ERROR: owner(" << p_owner
00245                     << ") owns(" << owns
00246                     << ") uses(" << uses
00247                     << ") owned_closure(" << owned_closure
00248                     << ") recv_ghost(" << recv_ghost
00249                     << ") send_ghost(" << send_ghost
00250                     << ") comm(" ;
00251           PairIterEntityComm ip = entity.comm();
00252           for ( ; ! ip.empty() ; ++ip ) {
00253             error_log << " " << ip->ghost_id << ":" << ip->proc ;
00254           }
00255           error_log << " )" << std::endl ;
00256         }
00257       }
00258     }
00259   }
00260 
00261   for ( std::vector<Entity*>::const_iterator
00262         i =  M.entity_comm().begin() ;
00263         i != M.entity_comm().end() ; ++i ) {
00264 
00265     const PairIterEntityComm ec = (*i)->comm();
00266 
00267     if ( ec.empty() ) {
00268       print_entity_key( error_log , M.mesh_meta_data(), (*i)->key() );
00269       error_log << " ERROR: in entity_comm but has no comm info" << std::endl ;
00270       result = false ;
00271     }
00272   }
00273 
00274   if ( M.entity_comm().size() != comm_count ) {
00275     error_log << " ERROR: entity_comm.size() = " << M.entity_comm().size();
00276     error_log << " != " << comm_count << " = entities with comm info" ;
00277     error_log << std::endl ;
00278     result = false ;
00279   }
00280 
00281   return result ;
00282 }
00283 
00284 //----------------------------------------------------------------------------
00285 // Packing my owned entities.
00286 
00287 void insert( std::vector<unsigned> & vec , unsigned val )
00288 {
00289   std::vector<unsigned>::iterator j =
00290     std::lower_bound( vec.begin() , vec.end() , val );
00291   if ( j == vec.end() || *j != val ) {
00292     vec.insert( j , val );
00293   }
00294 }
00295 
00296 void pack_owned_verify( CommAll & all , const BulkData & mesh )
00297 {
00298   const std::vector<Entity*> & entity_comm = mesh.entity_comm();
00299   const unsigned p_rank = all.parallel_rank();
00300 
00301   for ( std::vector<Entity*>::const_iterator
00302         i = entity_comm.begin() ; i != entity_comm.end() ; ++i ) {
00303 
00304     Entity & entity = **i ;
00305 
00306     if ( entity.owner_rank() == p_rank ) {
00307 
00308       std::vector<unsigned> share_proc ;
00309       std::vector<unsigned> ghost_proc ;
00310 
00311       const PairIterEntityComm comm = entity.comm();
00312 
00313       for ( size_t j = 0 ; j < comm.size() ; ++j ) {
00314         if ( comm[j].ghost_id == 0 ) {
00315           // Will be ordered by proc
00316           share_proc.push_back( comm[j].proc );
00317         }
00318         else {
00319           // No guarantee of ordering by proc
00320           insert( ghost_proc , comm[j].proc );
00321         }
00322       }
00323 
00324       const unsigned share_count = share_proc.size();
00325 
00326       for ( size_t j = 0 ; j < share_proc.size() ; ++j ) {
00327 
00328         // Sharing process, send sharing process list
00329 
00330         const unsigned p = share_proc[j] ;
00331 
00332         CommBuffer & buf = all.send_buffer( p );
00333 
00334         pack_entity_info( buf , entity );
00335 
00336         buf.pack<unsigned>( share_count );
00337 
00338         // Pack what the receiver should have:
00339         // My list, remove receiver, add myself
00340         size_t k = 0 ;
00341         for ( ; k < share_count && share_proc[k] < p_rank ; ++k ) {
00342           if ( k != j ) { buf.pack<unsigned>( share_proc[k] ); }
00343         }
00344         buf.pack<unsigned>( p_rank );
00345         for ( ; k < share_count ; ++k ) {
00346           if ( k != j ) { buf.pack<unsigned>( share_proc[k] ); }
00347         }
00348       }
00349 
00350       for ( size_t j = 0 ; j < ghost_proc.size() ; ++j ) {
00351         const unsigned p = ghost_proc[j] ;
00352 
00353         CommBuffer & buf = all.send_buffer( p );
00354 
00355         pack_entity_info( buf , entity );
00356 
00357         // What ghost subsets go to this process?
00358         unsigned count = 0 ;
00359         for ( size_t k = 0 ; k < comm.size() ; ++k ) {
00360           if ( comm[k].ghost_id != 0 && comm[k].proc == p ) {
00361             ++count ;
00362           }
00363         }
00364         buf.pack<unsigned>( count );
00365         for ( size_t k = 0 ; k < comm.size() ; ++k ) {
00366           if ( comm[k].ghost_id != 0 && comm[k].proc == p ) {
00367             buf.pack<unsigned>( comm[k].ghost_id );
00368           }
00369         }
00370       }
00371     }
00372   }
00373 }
00374 
00375 //----------------------------------------------------------------------------
00376 // Unpacking all of my not-owned entities.
00377 
00378 bool unpack_not_owned_verify( CommAll & comm_all ,
00379                               const BulkData & mesh ,
00380                               std::ostream & error_log )
00381 {
00382   const MetaData & meta = mesh.mesh_meta_data();
00383   Part * const       owns_part  = & meta.locally_owned_part();
00384   Part * const       uses_part  = & meta.locally_used_part();
00385   const PartVector & mesh_parts = meta.get_parts();
00386   const unsigned     p_rank = mesh.parallel_rank();
00387   const std::vector<Entity*> & entity_comm = mesh.entity_comm();
00388 
00389   bool result = true ;
00390 
00391   EntityKey             recv_entity_key ;
00392   unsigned              recv_owner_rank = 0 ;
00393   unsigned              recv_comm_count = 0 ;
00394   std::vector<Part*>    recv_parts ;
00395   std::vector<Relation> recv_relations ;
00396   std::vector<unsigned> recv_comm ;
00397 
00398   for ( std::vector<Entity*>::const_iterator
00399         i = entity_comm.begin() ; i != entity_comm.end() ; ++i ) {
00400 
00401     Entity & entity = **i ;
00402 
00403     if ( entity.owner_rank() != p_rank ) {
00404 
00405       const Bucket & bucket = entity.bucket();
00406 
00407       std::pair<const unsigned *,const unsigned *>
00408         part_ordinals = bucket.superset_part_ordinals();
00409 
00410       CommBuffer & buf = comm_all.recv_buffer( entity.owner_rank() );
00411 
00412       unpack_entity_info( buf , mesh ,
00413                           recv_entity_key , recv_owner_rank ,
00414                           recv_parts , recv_relations );
00415 
00416       recv_comm_count = 0 ;
00417       buf.unpack<unsigned>( recv_comm_count );
00418       recv_comm.resize( recv_comm_count );
00419       buf.unpack<unsigned>( & recv_comm[0] , recv_comm_count );
00420 
00421       // Match key and owner
00422 
00423       const bool bad_key = entity.key()        != recv_entity_key ;
00424       const bool bad_own = entity.owner_rank() != recv_owner_rank ;
00425       bool bad_part = false ;
00426       bool bad_rel  = false ;
00427       bool bad_comm = false ;
00428 
00429       // Compare communication information:
00430 
00431       if ( ! bad_key && ! bad_own ) {
00432         const PairIterEntityComm ec = entity.comm();
00433         const unsigned ec_size = ec.size();
00434         bad_comm = ec_size != recv_comm.size();
00435         if ( ! bad_comm ) {
00436           size_t j = 0 ;
00437           if ( in_shared( entity ) ) {
00438             for ( ; j < ec_size &&
00439                     ec[j].ghost_id == 0 &&
00440                     ec[j].proc   == recv_comm[j] ; ++j );
00441             bad_comm = j != ec_size ;
00442           }
00443           else {
00444             for ( ; j < ec_size &&
00445                     ec[j].ghost_id == recv_comm[j] &&
00446                     ec[j].proc   == entity.owner_rank() ; ++j );
00447             bad_comm = j != ec_size ;
00448           }
00449         }
00450       }
00451 
00452       // Compare everything but the owns part and uses part
00453 
00454       if ( ! bad_key && ! bad_own && ! bad_comm ) {
00455 
00456         const unsigned * k = part_ordinals.first ;
00457 
00458         std::vector<Part*>::iterator ip = recv_parts.begin();
00459 
00460         for ( ; ! bad_part && ip != recv_parts.end() ; ++ip ) {
00461           if ( owns_part != *ip ) {
00462             if ( uses_part != *ip ) {
00463               // All not-owned and not-uses parts must match:
00464               bad_part = k == part_ordinals.second ||
00465                          (*ip)->mesh_meta_data_ordinal() != *k ;
00466               ++k ;
00467             }
00468             else if ( k != part_ordinals.second &&
00469                      *k == uses_part->mesh_meta_data_ordinal() ) {
00470               // Uses-part matches
00471               ++k ;
00472             }
00473           }
00474         }
00475       }
00476 
00477       // Compare the closure relations:
00478       if ( ! bad_key && ! bad_own && ! bad_comm && ! bad_part ) {
00479 
00480         PairIterRelation ir = entity.relations();
00481 
00482         std::vector<Relation>::iterator jr = recv_relations.begin() ;
00483 
00484         for ( ; ! bad_rel && jr != recv_relations.end() &&
00485                 jr->entity_rank() < entity.entity_rank() ; ++jr , ++ir ) {
00486           bad_rel = ir.empty() || *jr != *ir ;
00487         }
00488       }
00489 
00490       if ( bad_key || bad_own || bad_comm || bad_part || bad_rel ) {
00491         error_log << "P" << p_rank << ": " ;
00492         print_entity_key( error_log , meta, entity.key() );
00493         error_log << " owner(" << entity.owner_rank() << ")" ;
00494 
00495         if ( bad_key || bad_own ) {
00496           error_log << " != received " ;
00497           print_entity_key( error_log , meta, recv_entity_key );
00498           error_log << " owner(" << recv_owner_rank
00499                     << ")" << std::endl ;
00500         }
00501         else if ( bad_comm ) {
00502           const PairIterEntityComm ec = entity.comm();
00503           if ( in_shared( entity ) ) {
00504             error_log << " sharing(" ;
00505             for ( size_t j = 0 ; j < ec.size() &&
00506                                  ec[j].ghost_id == 0 ; ++j ) {
00507               error_log << " " << ec[j].proc ;
00508             }
00509             error_log << " ) != received sharing(" ;
00510             for ( size_t j = 0 ; j < recv_comm.size() ; ++j ) {
00511               error_log << " " << recv_comm[j] ;
00512             }
00513             error_log << " )" << std::endl ;
00514           }
00515           else {
00516             error_log << " ghosting(" ;
00517             for ( size_t j = 0 ; j < ec.size() ; ++j ) {
00518               error_log << " (g" << ec[j].ghost_id ;
00519               error_log << ",p" << ec[j].proc ;
00520               error_log << ")" ;
00521             }
00522             error_log << " ) != received ghosting(" ;
00523             for ( size_t j = 0 ; j < recv_comm.size() ; ++j ) {
00524               error_log << " (g" << recv_comm[j] ;
00525               error_log << ",p" << entity.owner_rank();
00526               error_log << ")" ;
00527             }
00528             error_log << " )" << std::endl ;
00529           }
00530         }
00531         else if ( bad_part ) {
00532           error_log << " Parts( " ;
00533 
00534           for ( const unsigned * k = part_ordinals.first ;
00535                                  k < part_ordinals.second ; ++k ) {
00536             error_log << " \"" << mesh_parts[ *k ]->name() << "\"" ;
00537           }
00538           error_log << " ) != received Parts( " ;
00539 
00540           for ( std::vector<Part*>::iterator
00541                 ip =  recv_parts.begin();
00542                 ip != recv_parts.end() ; ++ip ) {
00543             error_log << " \"" << (*ip)->name() << "\"" ;
00544           }
00545           error_log << " )" << std::endl ;
00546         }
00547         else if ( bad_rel ) {
00548           error_log << " Relations(" ;
00549           PairIterRelation ir = entity.relations();
00550           for ( ; ! ir.empty() &&
00551                   ir->entity_rank() < entity.entity_rank() ; ++ir ) {
00552             error_log << " " << *ir ;
00553           }
00554           error_log << " ) != received Relations(" ;
00555           std::vector<Relation>::iterator jr = recv_relations.begin() ;
00556           for ( ; jr != recv_relations.end() &&
00557                   jr->entity_rank() < entity.entity_rank() ; ++jr ) {
00558             error_log << " " << *jr ;
00559           }
00560           error_log << " )" << std::endl ;
00561         }
00562         result = false ;
00563       }
00564     }
00565   }
00566 
00567   return result ;
00568 }
00569 
00570 } // namespace<>
00571 
00572 } // namespace mesh
00573 } // namespace stk
00574 

Generated on Tue Jul 13 09:27:31 2010 for Sierra Toolkit by  doxygen 1.4.7