MOAB: Mesh Oriented datABase  (version 5.4.1)
ReadUtil.cpp
Go to the documentation of this file.
00001 /**
00002  * MOAB, a Mesh-Oriented datABase, is a software component for creating,
00003  * storing and accessing finite element mesh data.
00004  *
00005  * Copyright 2004 Sandia Corporation.  Under the terms of Contract
00006  * DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government
00007  * retains certain rights in this software.
00008  *
00009  * This library is free software; you can redistribute it and/or
00010  * modify it under the terms of the GNU Lesser General Public
00011  * License as published by the Free Software Foundation; either
00012  * version 2.1 of the License, or (at your option) any later version.
00013  *
00014  */
00015 
00016 #ifdef WIN32
00017 #ifdef _DEBUG
00018 // turn off warnings that say they debugging identifier has been truncated
00019 // this warning comes up when using some STL containers
00020 #pragma warning( disable : 4786 )
00021 #endif
00022 #endif
00023 
00024 #include "ReadUtil.hpp"
00025 #include "moab/Core.hpp"
00026 #include "AEntityFactory.hpp"
00027 #include "moab/Error.hpp"
00028 #include "SequenceManager.hpp"
00029 #include "VertexSequence.hpp"
00030 #include "ElementSequence.hpp"
00031 
00032 namespace moab
00033 {
00034 
00035 #define RR \
00036     if( MB_SUCCESS != result ) return result
00037 
00038 ReadUtil::ReadUtil( Core* mdb, Error* /*error_handler*/ ) : ReadUtilIface(), mMB( mdb ) {}
00039 
00040 ErrorCode ReadUtil::get_node_coords( const int /*num_arrays*/,
00041                                      const int num_nodes,
00042                                      const int preferred_start_id,
00043                                      EntityHandle& actual_start_handle,
00044                                      std::vector< double* >& arrays,
00045                                      int sequence_size )
00046 {
00047     ErrorCode error;
00048     EntitySequence* seq = 0;
00049 
00050     if( num_nodes < 1 )
00051     {
00052         actual_start_handle = 0;
00053         arrays.clear();
00054         return MB_INDEX_OUT_OF_RANGE;
00055     }
00056 
00057     // Create an entity sequence for these nodes
00058     error = mMB->sequence_manager()->create_entity_sequence( MBVERTEX, num_nodes, 0, preferred_start_id,
00059                                                              actual_start_handle, seq, sequence_size );
00060 
00061     if( error != MB_SUCCESS ) return error;
00062 
00063     if( seq->start_handle() > actual_start_handle || seq->end_handle() < actual_start_handle ||
00064         seq->end_handle() - actual_start_handle + 1 < (unsigned)num_nodes )
00065         return MB_FAILURE;
00066 
00067     arrays.resize( 3 );
00068 
00069     error = static_cast< VertexSequence* >( seq )->get_coordinate_arrays( arrays[0], arrays[1], arrays[2] );
00070     for( unsigned i = 0; i < arrays.size(); ++i )
00071         if( arrays[i] ) arrays[i] += ( actual_start_handle - seq->start_handle() );
00072 
00073     return error;
00074 }
00075 
00076 ErrorCode ReadUtil::get_element_connect( const int num_elements,
00077                                          const int verts_per_element,
00078                                          const EntityType mdb_type,
00079                                          const int preferred_start_id,
00080                                          EntityHandle& actual_start_handle,
00081                                          EntityHandle*& array,
00082                                          int sequence_size )
00083 {
00084     ErrorCode error;
00085     EntitySequence* seq;
00086 
00087     if( num_elements < 1 )
00088     {
00089         actual_start_handle = 0;
00090         array               = 0;
00091         return MB_INDEX_OUT_OF_RANGE;
00092     }
00093 
00094     // if (mdb_type <= MBVERTEX || mdb_type >= MBPOLYHEDRON || mdb_type == MBPOLYGON)
00095     // return MB_TYPE_OUT_OF_RANGE;
00096 
00097     // Make an entity sequence to hold these elements
00098     error =
00099         mMB->sequence_manager()->create_entity_sequence( mdb_type, num_elements, verts_per_element, preferred_start_id,
00100                                                          actual_start_handle, seq, sequence_size );
00101     if( MB_SUCCESS != error ) return error;
00102 
00103     if( seq->start_handle() > actual_start_handle || seq->end_handle() < actual_start_handle ||
00104         seq->end_handle() - actual_start_handle + 1 < (unsigned)num_elements )
00105         return MB_FAILURE;
00106 
00107     // Get an array for the connectivity
00108     array = static_cast< ElementSequence* >( seq )->get_connectivity_array();
00109     if( !array ) return MB_FAILURE;
00110     array +=
00111         ( actual_start_handle - seq->start_handle() ) * static_cast< ElementSequence* >( seq )->nodes_per_element();
00112 
00113     return error;
00114 }
00115 
00116 ErrorCode ReadUtil::create_entity_sets( EntityID num_sets,
00117                                         const unsigned* flags,
00118                                         EntityID start_id,
00119                                         EntityHandle& start_handle )
00120 {
00121     if( num_sets < 1 )
00122     {
00123         start_handle = 0;
00124         return MB_INDEX_OUT_OF_RANGE;
00125     }
00126 
00127     ErrorCode error;
00128     EntitySequence* seq;
00129     error = mMB->sequence_manager()->create_meshset_sequence( num_sets, start_id, flags, start_handle, seq );
00130     if( MB_SUCCESS != error ) return error;
00131 
00132     if( seq->start_handle() > start_handle || seq->end_handle() < start_handle ||
00133         seq->end_handle() - start_handle + 1 < (EntityHandle)num_sets )
00134         return MB_FAILURE;
00135 
00136     return MB_SUCCESS;
00137 }
00138 
00139 ErrorCode ReadUtil::update_adjacencies( const EntityHandle start_handle,
00140                                         const int number_elements,
00141                                         const int number_vertices_per_element,
00142                                         const EntityHandle* conn_array )
00143 {
00144     EntityHandle tmp_hndl    = start_handle;
00145     AEntityFactory* adj_fact = mMB->a_entity_factory();
00146 
00147     // Iterate over the elements and update adjacency information
00148     if( adj_fact != NULL && adj_fact->vert_elem_adjacencies() )
00149     {
00150         int j = 0;
00151         for( int i = 0; i < number_elements; i++ )
00152         {
00153             adj_fact->notify_create_entity( tmp_hndl, ( conn_array + j ), number_vertices_per_element );
00154             tmp_hndl++;
00155             j += number_vertices_per_element;
00156         }
00157     }
00158 
00159     return MB_SUCCESS;
00160 }
00161 
00162 ErrorCode ReadUtil::gather_related_ents( Range& partition, Range& related_ents, EntityHandle* file_set )
00163 {
00164     // Loop over any sets, getting contained ents
00165     std::pair< Range::const_iterator, Range::const_iterator > pair_it = partition.equal_range( MBENTITYSET );
00166 
00167     ErrorCode result = MB_SUCCESS;
00168     for( Range::const_iterator rit = pair_it.first; rit != pair_it.second; ++rit )
00169     {
00170         ErrorCode tmp_result = mMB->get_entities_by_handle( *rit, related_ents, Interface::UNION );
00171         if( MB_SUCCESS != tmp_result ) result = tmp_result;
00172     }
00173     RR;
00174 
00175     // Gather adjacent ents of other dimensions
00176     Range tmp_ents;
00177     for( int dim = 3; dim >= 0; dim-- )
00178     {
00179         tmp_ents.clear();
00180         ErrorCode tmp_result = mMB->get_adjacencies( related_ents, dim, false, tmp_ents, Interface::UNION );
00181         if( MB_SUCCESS != tmp_result )
00182             result = tmp_result;
00183         else
00184             related_ents.merge( tmp_ents );
00185     }
00186     RR;
00187 
00188     // Related ents includes the partition itself
00189     related_ents.merge( partition );
00190 
00191     // Get contains-related sets
00192     Range tmp_ents3, last_related;
00193     if( file_set )
00194         result = mMB->get_entities_by_type( *file_set, MBENTITYSET, tmp_ents3 );
00195     else
00196         result = mMB->get_entities_by_type( 0, MBENTITYSET, tmp_ents3 );RR;
00197 
00198     while( related_ents.size() != last_related.size() )
00199     {
00200         last_related = related_ents;
00201         for( Range::iterator rit = tmp_ents3.begin(); rit != tmp_ents3.end(); ++rit )
00202         {
00203             if( related_ents.find( *rit ) != related_ents.end() ) continue;
00204 
00205             tmp_ents.clear();
00206             result = mMB->get_entities_by_handle( *rit, tmp_ents, true );RR;
00207             Range tmp_ents2 = intersect( tmp_ents, related_ents );
00208 
00209             // If the intersection is not empty, set is related
00210             if( !tmp_ents2.empty() ) related_ents.insert( *rit );
00211         }
00212     }
00213 
00214     // Get child-related sets
00215     last_related.clear();
00216     while( related_ents.size() != last_related.size() )
00217     {
00218         last_related                                                      = related_ents;
00219         std::pair< Range::const_iterator, Range::const_iterator > it_pair = last_related.equal_range( MBENTITYSET );
00220 
00221         for( Range::const_iterator rit = it_pair.first; rit != it_pair.second; ++rit )
00222         {
00223             // Get all children and add to related ents
00224             tmp_ents.clear();
00225             result = mMB->get_child_meshsets( *rit, tmp_ents, 0 );RR;
00226             related_ents.merge( tmp_ents );
00227         }
00228     }
00229 
00230     // Get parent-related sets
00231     last_related.clear();
00232     while( related_ents.size() != last_related.size() )
00233     {
00234         last_related                                                      = related_ents;
00235         std::pair< Range::const_iterator, Range::const_iterator > it_pair = last_related.equal_range( MBENTITYSET );
00236 
00237         for( Range::const_iterator rit = it_pair.first; rit != it_pair.second; ++rit )
00238         {
00239             // Get all parents and add to related ents
00240             tmp_ents.clear();
00241             result = mMB->get_parent_meshsets( *rit, tmp_ents, 0 );RR;
00242             related_ents.merge( tmp_ents );
00243         }
00244     }
00245 
00246     return MB_SUCCESS;
00247 }
00248 
00249 ErrorCode ReadUtil::get_ordered_vertices( EntityHandle* bound_ents,
00250                                           int* sense,
00251                                           int bound_size,
00252                                           int dim,
00253                                           EntityHandle* bound_verts,
00254                                           EntityType& etype )
00255 {
00256     // Get dimension of bounding entities
00257     int bound_dim = CN::Dimension( TYPE_FROM_HANDLE( bound_ents[0] ) );
00258     int indices[MAX_SUB_ENTITY_VERTICES];
00259     const EntityHandle* connect = NULL;
00260     std::vector< EntityHandle > tmp_connect;
00261 
00262     // Find the right entity type based on # bounding ents
00263     int numv = 0, num_connect = 0;
00264     ErrorCode result;
00265     for( EntityType t = MBEDGE; t < MBENTITYSET; t++ )
00266     {
00267         int nindex = CN::NumSubEntities( t, bound_dim );
00268         if( CN::Dimension( t ) != dim || nindex != bound_size ) continue;
00269 
00270         // Fill in vertices from bounding entity vertices
00271         int nverts = CN::VerticesPerEntity( t );
00272         std::fill( bound_verts, bound_verts + nverts, 0 );
00273         for( int index = 0; index < nindex; index++ )
00274         {
00275             result = mMB->get_connectivity( bound_ents[index], connect, num_connect, false, &tmp_connect );
00276             if( MB_SUCCESS != result ) return result;
00277 
00278             CN::SubEntityVertexIndices( t, bound_dim, index, indices );
00279 
00280             for( int c = 0; c < num_connect; c++ )
00281             {
00282                 if( !bound_verts[indices[c]] )
00283                 {
00284                     bound_verts[indices[c]] = ( sense[index] > 0 ) ? connect[c] : connect[num_connect - c - 1];
00285                     numv++;
00286                 }
00287             }
00288             if( numv == nverts )
00289             {
00290                 etype = t;
00291                 return MB_SUCCESS;
00292             }
00293         }
00294     }
00295 
00296     // If we get here, we didn't get full connectivity
00297     etype = MBMAXTYPE;
00298     return MB_FAILURE;
00299 }
00300 
00301 static ErrorCode check_int_tag( Interface* mb, Tag tag )
00302 {
00303     int size;
00304     DataType type;
00305     ErrorCode rval = mb->tag_get_bytes( tag, size );
00306     if( MB_SUCCESS != rval ) return rval;
00307     if( size != sizeof( int ) ) return MB_TYPE_OUT_OF_RANGE;
00308     rval = mb->tag_get_data_type( tag, type );
00309     if( type != MB_TYPE_OPAQUE && type != MB_TYPE_INTEGER ) return MB_TYPE_OUT_OF_RANGE;
00310 
00311     return MB_SUCCESS;
00312 }
00313 
00314 ErrorCode ReadUtil::assign_ids( Tag id_tag, const Range& ents, int start )
00315 {
00316     ErrorCode rval = check_int_tag( mMB, id_tag );
00317     if( MB_SUCCESS != rval ) return rval;
00318 
00319     Range tmp_range;
00320     std::vector< int > data;
00321     for( Range::const_pair_iterator i = ents.pair_begin(); i != ents.pair_end(); ++i )
00322     {
00323         data.resize( i->second + 1 - i->first );
00324         for( std::vector< int >::iterator j = data.begin(); j != data.end(); ++j )
00325             *j = start++;
00326         tmp_range.clear();
00327         tmp_range.insert( i->first, i->second );
00328         rval = mMB->tag_set_data( id_tag, tmp_range, &data[0] );
00329         if( MB_SUCCESS != rval ) return rval;
00330     }
00331 
00332     return MB_SUCCESS;
00333 }
00334 
00335 ErrorCode ReadUtil::assign_ids( Tag id_tag, const EntityHandle* ents, size_t num_ents, int start )
00336 {
00337     ErrorCode rval = check_int_tag( mMB, id_tag );
00338     if( MB_SUCCESS != rval ) return rval;
00339 
00340     std::vector< int > data;
00341     const EntityHandle* const end = ents + num_ents;
00342     const EntityHandle* i         = ents;
00343     while( i != end )
00344     {
00345         const EntityHandle* next = std::find( i, end, 0u );
00346         size_t size              = next - i;
00347         if( !size )
00348         {
00349             ++i;
00350             continue;
00351         }
00352 
00353         int id = start + ( i - ents );
00354         data.resize( size );
00355         for( std::vector< int >::iterator j = data.begin(); j != data.end(); ++j )
00356             *j = id++;
00357 
00358         rval = mMB->tag_set_data( id_tag, i, size, &data[0] );
00359         if( MB_SUCCESS != rval ) return rval;
00360     }
00361 
00362     return MB_SUCCESS;
00363 }
00364 
00365 ErrorCode ReadUtil::create_gather_set( EntityHandle& gather_set )
00366 {
00367     ErrorCode rval = mMB->create_meshset( MESHSET_SET, gather_set );
00368     if( MB_SUCCESS != rval ) return rval;
00369 
00370     Tag gather_set_tag;
00371     rval = mMB->tag_get_handle( "GATHER_SET", 1, MB_TYPE_INTEGER, gather_set_tag, MB_TAG_CREAT | MB_TAG_SPARSE );
00372     if( MB_SUCCESS != rval ) return rval;
00373 
00374     int gather_val = 1;
00375     rval           = mMB->tag_set_data( gather_set_tag, &gather_set, 1, &gather_val );
00376     if( MB_SUCCESS != rval ) return rval;
00377 
00378     return MB_SUCCESS;
00379 }
00380 
00381 ErrorCode ReadUtil::get_gather_set( EntityHandle& gather_set )
00382 {
00383     Tag gather_set_tag;
00384     ErrorCode rval = mMB->tag_get_handle( "GATHER_SET", 1, MB_TYPE_INTEGER, gather_set_tag, MB_TAG_SPARSE );
00385     if( MB_SUCCESS != rval ) return rval;
00386 
00387     int gather_val = 1;
00388     void* vals[]   = { &gather_val };
00389     Range gather_sets;
00390     rval = mMB->get_entities_by_type_and_tag( 0, MBENTITYSET, &gather_set_tag, vals, 1, gather_sets );
00391     if( MB_SUCCESS != rval ) return rval;
00392 
00393     if( gather_sets.empty() ) return MB_ENTITY_NOT_FOUND;
00394 
00395     gather_set = gather_sets[0];
00396 
00397     return MB_SUCCESS;
00398 }
00399 
00400 }  // namespace moab
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines