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