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