MOAB: Mesh Oriented datABase
(version 5.2.1)
|
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 #ifndef MOAB_IMPL_GENERAL_HPP 00017 #define MOAB_IMPL_GENERAL_HPP 00018 00019 #include "moab/Interface.hpp" 00020 #include "moab/ReaderIface.hpp" 00021 #include <map> 00022 #include <list> 00023 00024 namespace moab 00025 { 00026 00027 class WriteUtil; 00028 class ReadUtil; 00029 class ScdInterface; 00030 class AEntityFactory; 00031 class SequenceManager; 00032 class Error; 00033 class HomCoord; 00034 class ReaderWriterSet; 00035 class EntitySequence; 00036 class FileOptions; 00037 class SetIterator; 00038 00039 #ifdef MOAB_HAVE_AHF 00040 class HalfFacetRep; 00041 #endif 00042 00043 /**\class Core 00044 * \brief Implementation of MOAB Interface 00045 * Implementation of the MOAB Interface class. You shouldn't call functions directly 00046 * on an object of type Core (use Interface instead), unless you really have to access 00047 * non-API functionality. 00048 */ 00049 class Core : public Interface 00050 { 00051 00052 public: 00053 friend class SetIterator; 00054 00055 //! constructor 00056 Core(); 00057 00058 //! destructor 00059 ~Core(); 00060 00061 //! Get a pointer to an internal MOAB interface 00062 //!\return NULL if not found, iterface pointer otherwise 00063 virtual ErrorCode query_interface_type( const std::type_info& iface_type, void*& iface ); 00064 00065 //! Release reference to MB interface 00066 virtual ErrorCode release_interface_type( const std::type_info& iface_type, void* iface ); 00067 00068 virtual int QueryInterface( const MBuuid& uuid, UnknownInterface** iface ); 00069 00070 //! Returns the major.minor version number of the implementation 00071 /** 00072 \param iface_name If non-NULL, will be filled in with a string, possibly 00073 containing implementation-specific information 00074 */ 00075 virtual float impl_version( std::string* version_string = NULL ); 00076 00077 //! get the type from a handle, returns type 00078 virtual EntityType type_from_handle( const EntityHandle handle ) const; 00079 00080 //! get the id from a handle, returns id 00081 virtual EntityID id_from_handle( const EntityHandle handle ) const; 00082 00083 //! get a handle from an id and type 00084 virtual ErrorCode handle_from_id( const EntityType type, const EntityID id, EntityHandle& handle ) const; 00085 00086 virtual int dimension_from_handle( const EntityHandle ) const; 00087 00088 //! load mesh from data in file 00089 //! NOTE: if there is mesh already present, the new mesh will be added 00090 virtual ErrorCode load_mesh( const char* file_name, const int* active_block_id_list = NULL, 00091 const int num_blocks = 0 ); 00092 00093 /**Load or import a file. */ 00094 virtual ErrorCode load_file( const char* file_name, const EntityHandle* file_set = 0, const char* options = 0, 00095 const char* set_tag_name = 0, const int* set_tag_values = 0, 00096 int num_set_tag_values = 0 ); 00097 00098 /**Load or import a file. */ 00099 ErrorCode serial_load_file( const char* file_name, const EntityHandle* file_set, const FileOptions& opts, 00100 const ReaderIface::SubsetList* subset_list = 0, const Tag* file_id_tag = 0 ); 00101 00102 ErrorCode serial_read_tag( const char* file_name, const char* tag_name, const FileOptions& opts, 00103 std::vector< int >& tag_vals, const ReaderIface::SubsetList* subset_list = 0 ); 00104 00105 virtual ErrorCode write_mesh( const char* file_name, const EntityHandle* output_list = NULL, 00106 const int num_sets = 0 ); 00107 /** Write or export a file. */ 00108 virtual ErrorCode write_file( const char* file_name, const char* file_type = 0, const char* options = 0, 00109 const EntityHandle* output_sets = 0, int num_output_sets = 0, const Tag* tag_list = 0, 00110 int num_tags = 0 ); 00111 00112 /** Write or export a file */ 00113 virtual ErrorCode write_file( const char* file_name, const char* file_type, const char* options, 00114 const Range& output_sets, const Tag* tag_list = 0, int num_tags = 0 ); 00115 00116 //! deletes all mesh entities from this datastore 00117 virtual ErrorCode delete_mesh(); 00118 00119 //! get overall geometric dimension 00120 virtual ErrorCode get_dimension( int& dim ) const; 00121 00122 //! set overall geometric dimension 00123 /** Returns error if setting to 3 dimensions, mesh has been created, and 00124 * there are only 2 dimensions on that mesh 00125 */ 00126 virtual ErrorCode set_dimension( const int dim ); 00127 00128 //! get blocked vertex coordinates for all vertices 00129 /** Blocked = all x, then all y, etc. 00130 */ 00131 virtual ErrorCode get_vertex_coordinates( std::vector< double >& coords ) const; 00132 00133 //! get pointers to coordinate data 00134 virtual ErrorCode coords_iterate( Range::const_iterator iter, Range::const_iterator end, double*& xcoords_ptr, 00135 double*& ycoords_ptr, double*& zcoords_ptr, int& count ); 00136 00137 //! get the coordinate information for this handle if it is of type Vertex 00138 //! otherwise, return an error 00139 virtual ErrorCode get_coords( const Range& entity_handles, double* coords ) const; 00140 00141 virtual ErrorCode get_coords( const EntityHandle* entity_handles, const int num_entities, double* coords ) const; 00142 00143 virtual ErrorCode get_coords( const EntityHandle entity_handle, const double*& x, const double*& y, 00144 const double*& z ) const; 00145 00146 virtual ErrorCode get_coords( const Range& entity_handles, double* x_coords, double* y_coords, 00147 double* z_coords ) const; 00148 00149 //! set the coordinate information for this handle if it is of type Vertex 00150 //! otherwise, return an error 00151 virtual ErrorCode set_coords( const EntityHandle* entity_handles, const int num_entities, const double* coords ); 00152 00153 //! set the coordinate information for this handle if it is of type Vertex 00154 //! otherwise, return an error 00155 virtual ErrorCode set_coords( Range entity_handles, const double* coords ); 00156 00157 //! get global connectivity array for specified entity type 00158 /** Assumes just vertices, no higher order nodes 00159 */ 00160 virtual ErrorCode get_connectivity_by_type( const EntityType type, std::vector< EntityHandle >& connect ) const; 00161 00162 //! get pointer to connectivity data 00163 virtual ErrorCode connect_iterate( Range::const_iterator iter, Range::const_iterator end, EntityHandle*& connect, 00164 int& verts_per_entity, int& count ); 00165 00166 //! Gets the connectivity for an element EntityHandle. 00167 /** For non-element handles (ie, MeshSets), 00168 * returns an error. Connectivity data is copied from the database into the vector 00169 * <em>connectivity</em>. The nodes in <em>connectivity</em> are properly ordered. 00170 * \param entity_handle EntityHandle to get connectivity of. 00171 * \param connectivity Vector in which connectivity of <em>entity_handle</em> is returned. 00172 * Should contain MeshVertices. 00173 * \param corners_only If true, returns only corner vertices, otherwise returns all of them 00174 * (including any higher-order vertices) 00175 * 00176 * Example: \code 00177 * std::vector<EntityHandle> conn; 00178 * get_connectivity( entity_handle, conn ); \endcode 00179 */ 00180 virtual ErrorCode get_connectivity( const EntityHandle* entity_handles, const int num_handles, 00181 std::vector< EntityHandle >& connectivity, bool corners_only = false, 00182 std::vector< int >* offsets = NULL ) const; 00183 00184 //! Gets the connectivity for a vector of elements 00185 /** Same as vector-based version except range is returned (unordered!) 00186 */ 00187 virtual ErrorCode get_connectivity( const EntityHandle* entity_handles, const int num_handles, Range& connectivity, 00188 bool corners_only = false ) const; 00189 00190 //! Gets the connectivity for elements 00191 /** Same as vector-based version except range is returned (unordered!) 00192 */ 00193 virtual ErrorCode get_connectivity( const Range& entity_handles, Range& connectivity, 00194 bool corners_only = false ) const; 00195 00196 //! Gets a pointer to constant connectivity data of <em>entity_handle</em> 00197 /** Sets <em>number_nodes</em> equal to the number of nodes of the <em> 00198 entity_handle </em>. Faster then the other <em>get_connectivity</em> function. 00199 The nodes in 'connectivity' are properly ordered. 00200 \param entity_handle EntityHandle to get connectivity of. 00201 \param connectivity Array in which connectivity of <em>entity_handle</em> is returned. 00202 Should contain MBVERTEX's. 00203 \param num_nodes Number of MeshVertices in array <em>connectivity</em>. 00204 00205 Example: \code 00206 const EntityHandle* conn; 00207 int number_nodes = 0; 00208 get_connectivity( entity_handle, conn, number_nodes ); \endcode 00209 00210 Example2: \code 00211 std::vector<EntityHandle> sm_storage; 00212 const EntityHandle* conn; 00213 int number_nodes; 00214 get_connectivity( handle, conn, number_nodes, false, &sm_storage ); 00215 if (conn == &sm_storage[0]) 00216 std::cout << "Structured mesh element" << std::endl; 00217 \endcode 00218 */ 00219 virtual ErrorCode get_connectivity( const EntityHandle entity_handle, const EntityHandle*& connectivity, 00220 int& num_nodes, bool corners_only = false, 00221 std::vector< EntityHandle >* storage = 0 ) const; 00222 00223 //! Sets the connectivity for an EntityHandle. For non-element handles, return an error. 00224 /** Connectivity is stored exactly as it is ordered in vector <em>connectivity</em>. 00225 \param entity_handle EntityHandle to set connectivity of. 00226 \param connect Vector containing new connectivity of <em>entity_handle</em>. 00227 \param num_connect Number of vertices in <em>connect</em> 00228 00229 Example: \code 00230 std::vector<EntityHandle> conn(3); 00231 conn[0] = node1; 00232 conn[1] = node2; 00233 conn[2] = node3; 00234 set_connectivity( entity_handle, conn, 3 ); \endcode */ 00235 virtual ErrorCode set_connectivity( const EntityHandle entity_handle, EntityHandle* connect, 00236 const int num_connect ); 00237 00238 //! get the adjacencies associated with a set of entities 00239 /** \param from_entities vector of EntityHandle to get adjacencies of. 00240 \param to_dimension Dimension of desired adjacency information. 00241 \param adj_entities Vector in which adjacent EntityHandles are returned. 00242 \param operation_type enum of INTERSECT or UNION. Defines whether to take 00243 the intersection or union of the set of adjacencies recovered for the from_entities. 00244 00245 The adjacent entities in vector <em>adjacencies</em> are not in any particular 00246 order. 00247 00248 Example: \code 00249 // get the set of edges that are adjacent to all entities in the from_entities list 00250 std::vector<EntityHandle> from_entities = {hex1, hex2}; 00251 std::vector<EntityHandle> adjacencies; 00252 get_adjacencies( from_entities, MB_1D_ENTITY, adjacencies ); 00253 \endcode */ 00254 00255 virtual ErrorCode get_adjacencies( const EntityHandle* from_entities, const int num_entities, 00256 const int to_dimension, const bool create_if_missing, 00257 std::vector< EntityHandle >& adj_entities, 00258 const int operation_type = Interface::INTERSECT ); 00259 00260 virtual ErrorCode get_adjacencies( const EntityHandle* from_entities, const int num_entities, 00261 const int to_dimension, const bool create_if_missing, Range& adj_entities, 00262 const int operation_type = Interface::INTERSECT ); 00263 00264 virtual ErrorCode get_adjacencies( const Range& from_entities, const int to_dimension, const bool create_if_missing, 00265 Range& adj_entities, const int operation_type = Interface::INTERSECT ); 00266 00267 /**\brief Get a ptr to adjacency lists 00268 * Get a pointer to adjacency lists. These lists are std::vector<EntityHandle>, which are 00269 * pointed to by adjs[i]. Adjacencies are not guaranteed to be in order of increasing 00270 * dimension. Only a const version of this function is given, because adjacency data is managed 00271 * more carefully in MOAB and should be treated as read-only by applications. If adjacencies 00272 * have not yet been initialized, adjs_ptr will be NULL (i.e. adjs_ptr == NULL). There may also 00273 * be NULL entries for individual entities, i.e. adjs_ptr[i] == NULL. \param iter Iterator to 00274 * beginning of entity range desired \param end End iterator for which adjacencies are requested 00275 * \param adjs_ptr Pointer to pointer to const std::vector<EntityHandle>; each member of that 00276 * array is the vector of adjacencies for this entity \param count Number of entities in the 00277 * contiguous chunk starting from *iter 00278 */ 00279 ErrorCode adjacencies_iterate( Range::const_iterator iter, Range::const_iterator end, 00280 const std::vector< EntityHandle >**& adjs_ptr, int& count ); 00281 00282 /**\brief Get all vertices for input entities 00283 * 00284 * Special case of get_adjacencies where to_dimension == 0 00285 * and operation_type == Interface::UNION. 00286 *\Note This is not a variation of get_connectivity because 00287 * the behavior is different for polyhedra. 00288 */ 00289 virtual ErrorCode get_vertices( const Range& from_entities, Range& vertices ); 00290 00291 //! Adds adjacencies 00292 /** \param from_handle entities 00293 \param both_ways add the adjacency information to both the 00294 to_handle and and the from_from :handle 00295 00296 Example: \code 00297 */ 00298 virtual ErrorCode add_adjacencies( const EntityHandle from_handle, const EntityHandle* to_handles, 00299 const int num_handles, bool both_ways ); 00300 00301 //! Adds adjacencies; same as vector-based, but with range instead 00302 virtual ErrorCode add_adjacencies( const EntityHandle from_handle, Range& adjacencies, bool both_ways ); 00303 00304 //! Removes adjacencies 00305 /** \param handle EntityHandle to get adjacencies of. 00306 00307 Example: \code 00308 */ 00309 virtual ErrorCode remove_adjacencies( const EntityHandle from_handle, const EntityHandle* to_handles, 00310 const int num_handles ); 00311 00312 //! Retrieves all entities in the database of given dimension. 00313 /** \param dimension Dimension of entities desired. 00314 \param entities Range in which entities of dimension <em>dimension</em> are returned. 00315 00316 Example: \code 00317 int dimension = 2; 00318 Range entities; 00319 get_entities_by_dimension( dimension, entities ); //get 2D EntityHandles in the database 00320 \endcode */ 00321 virtual ErrorCode get_entities_by_dimension( const EntityHandle meshset, const int dimension, Range& entities, 00322 const bool recursive = false ) const; 00323 00324 virtual ErrorCode get_entities_by_dimension( const EntityHandle meshset, const int dimension, 00325 std::vector< EntityHandle >& entities, 00326 const bool recursive = false ) const; 00327 00328 //! Retrieves all entities in the data base of given type. 00329 /** \param type EntityType of entities desired (ie, MBHEX, MBEDGE, MBTRI, etc ) 00330 \param entities Range in which entities of EntityType <em>type</em> are returned. 00331 00332 Example: \code 00333 EntityType type = MBTET; 00334 Range entities; 00335 get_entities_by_dimension( type, entities ); //get MBTET type EntityHandles in the database 00336 \endcode */ 00337 virtual ErrorCode get_entities_by_type( const EntityHandle meshset, const EntityType type, Range& entities, 00338 const bool recursive = false ) const; 00339 00340 //! Retrieves all entities in the data base of given type. 00341 /** \param type EntityType of entities desired (ie, MBHEX, MBEDGE, MBTRI, etc ) 00342 \param entities Range in which entities of EntityType <em>type</em> are returned. 00343 00344 Example: \code 00345 EntityType type = MBTET; 00346 Range entities; 00347 get_entities_by_dimension( type, entities ); //get MBTET type EntityHandles in the database 00348 \endcode */ 00349 virtual ErrorCode get_entities_by_type( const EntityHandle meshset, const EntityType type, 00350 std::vector< EntityHandle >& entities, const bool recursive = false ) const; 00351 00352 virtual ErrorCode get_entities_by_type_and_tag( const EntityHandle meshset, const EntityType type, 00353 const Tag* tag_handles, const void* const* values, 00354 const int num_tags, Range& entities, 00355 const int condition = Interface::INTERSECT, 00356 const bool recursive = false ) const; 00357 00358 //! Retrieves all entities in the data base 00359 /** \param entities Range in which entities of EntityType <em>type</em> are returned. 00360 00361 Example: \code 00362 Range entities; 00363 get_entities( entities ); //get MBTET type EntityHandles in the database 00364 \endcode */ 00365 virtual ErrorCode get_entities_by_handle( const EntityHandle meshset, Range& entities, 00366 const bool recursive = false ) const; 00367 00368 //! Retrieves all entities in the data base 00369 /** \param entities Range in which entities of EntityType <em>type</em> are returned. 00370 00371 Example: \code 00372 Range entities; 00373 get_entities( entities ); //get MBTET type EntityHandles in the database 00374 \endcode */ 00375 virtual ErrorCode get_entities_by_handle( const EntityHandle meshset, std::vector< EntityHandle >& entities, 00376 const bool recursive = false ) const; 00377 00378 //! Retrieves all entities in the database of given dimension. 00379 /** \param dimension Dimension of entities desired. 00380 \param entities Range in which entities of dimension <em>dimension</em> are returned. 00381 00382 Example: \code 00383 int dimension = 2; 00384 Range entities; 00385 get_entities_by_dimension( dimension, entities ); //get 2D EntityHandles in the database 00386 \endcode */ 00387 virtual ErrorCode get_number_entities_by_dimension( const EntityHandle meshset, const int dimension, 00388 int& num_entities, const bool recursive = false ) const; 00389 00390 //! Retrieves all entities in the data base of given type. 00391 /** \param type EntityType of entities desired (ie, MBHEX, MBEDGE, MBTRI, etc ) 00392 \param entities Range in which entities of EntityType <em>type</em> are returned. 00393 00394 Example: \code 00395 EntityType type = MBTET; 00396 Range entities; 00397 get_entities_by_dimension( type, entities ); //get MBTET type EntityHandles in the database 00398 \endcode */ 00399 virtual ErrorCode get_number_entities_by_type( const EntityHandle meshset, const EntityType type, int& num_entities, 00400 const bool recursive = false ) const; 00401 00402 virtual ErrorCode get_number_entities_by_type_and_tag( const EntityHandle meshset, const EntityType type, 00403 const Tag* tag_handles, const void* const* values, 00404 const int num_tags, int& num_entities, 00405 const int condition = Interface::INTERSECT, 00406 const bool recursive = false ) const; 00407 00408 //! Retrieves all entities in the data base 00409 /** \param entities Range in which entities of EntityType <em>type</em> are returned. 00410 00411 Example: \code 00412 Range entities; 00413 get_entities( entities ); //get MBTET type EntityHandles in the database 00414 \endcode */ 00415 virtual ErrorCode get_number_entities_by_handle( const EntityHandle meshset, int& num_entities, 00416 const bool recursive = false ) const; 00417 00418 //! Creates an element based on the type and connectivity. 00419 /** If connectivity vector is not correct for EntityType <em>type</em> (ie, a vector with 00420 3 vertices is passed in to make an MBQUAD), the function returns MB_FAILURE. 00421 \param type Type of element to create. (MBTET, MBTRI, MBKNIFE, etc.) 00422 \param connectivity Vector containing connectivity of element to create. 00423 \param handle EntityHandle representing the newly created element in the database. 00424 00425 Example: \code 00426 EntityType type = MBQUAD; 00427 std::vector<EntityHandle> connectivity(4); 00428 quad_conn[0] = vertex0; 00429 quad_conn[1] = vertex1; 00430 quad_conn[2] = vertex2; 00431 quad_conn[3] = vertex3; 00432 EntityHandle element_handle; 00433 create_element( type, connectivity, element_handle ); \endcode */ 00434 virtual ErrorCode create_element( const EntityType type, const EntityHandle* connectivity, const int num_nodes, 00435 EntityHandle& element_handle ); 00436 00437 //! Creates a vertex based on coordinates. 00438 /** 00439 \param coordinates Array that has 3 doubles in it. 00440 \param entity_handle EntityHandle representing the newly created vertex in the database. 00441 00442 Example: \code 00443 double *coordinates = double[3]; 00444 coordinates[0] = 1.034; 00445 coordinates[1] = 23.23; 00446 coordinates[2] = -0.432; 00447 EntityHandle entity_handle = 0; 00448 create_vertex( coordinates, entity_handle ); \endcode */ 00449 virtual ErrorCode create_vertex( const double coordinates[3], EntityHandle& entity_handle ); 00450 00451 //! Create a set of vertices with the specified coordinates 00452 /** 00453 \param coordinates Array that has 3*n doubles in it. 00454 \param nverts Number of vertices to create 00455 \param entity_handles Range passed back with new vertex handles 00456 */ 00457 virtual ErrorCode create_vertices( const double* coordinates, const int nverts, Range& entity_handles ); 00458 00459 //! merges two entities 00460 virtual ErrorCode merge_entities( EntityHandle entity_to_keep, EntityHandle entity_to_remove, bool auto_merge, 00461 bool delete_removed_entity ); 00462 00463 //! Removes entities in a vector from the data base. 00464 /** If any of the entities are contained in any meshsets, it is removed from those meshsets 00465 which were created with MESHSET_TRACK_OWNER option bit set. Tags for <em>entity<\em> are 00466 removed as part of this function. 00467 \param entities 1d vector of entities to delete 00468 \param num_entities Number of entities in 1d vector 00469 */ 00470 virtual ErrorCode delete_entities( const EntityHandle* entities, const int num_entities ); 00471 00472 //! Removes entities in a range from the data base. 00473 /** If any of the entities are contained in any meshsets, it is removed from those meshsets 00474 which were created with MESHSET_TRACK_OWNER option bit set. Tags for <em>entity<\em> are 00475 removed as part of this function. 00476 \param entities Range of entities to delete 00477 */ 00478 virtual ErrorCode delete_entities( const Range& entities ); 00479 00480 virtual ErrorCode list_entities( const Range& entities ) const; 00481 00482 virtual ErrorCode list_entities( const EntityHandle* entities, const int num_entities ) const; 00483 00484 virtual ErrorCode list_entity( const EntityHandle entity ) const; 00485 00486 typedef unsigned long long type_memstorage; 00487 00488 //! function object for recieving events from MB of higher order nodes 00489 //! added to entities 00490 class HONodeAddedRemoved 00491 { 00492 public: 00493 HONodeAddedRemoved() {} 00494 virtual ~HONodeAddedRemoved() {} 00495 //! node_added called when a node was added to an element's connectivity array 00496 //! note: connectivity array of element may be incomplete (corner nodes will exist always) 00497 virtual void node_added( EntityHandle node, EntityHandle element ); 00498 virtual void node_removed( EntityHandle node ); 00499 }; 00500 00501 virtual ErrorCode convert_entities( const EntityHandle meshset, const bool mid_side, const bool mid_face, 00502 const bool mid_volume, Interface::HONodeAddedRemoved* function_object = 0 ); 00503 00504 //! function to get the side number given two elements; returns 00505 //! MB_FAILURE if child not related to parent; does *not* create adjacencies 00506 //! between parent and child 00507 virtual ErrorCode side_number( const EntityHandle parent, const EntityHandle child, int& side_number, int& sense, 00508 int& offset ) const; 00509 00510 //! given an entity and the connectivity and type of one of its subfacets, find the 00511 //! high order node on that subfacet, if any 00512 virtual ErrorCode high_order_node( const EntityHandle parent_handle, const EntityHandle* subfacet_conn, 00513 const EntityType subfacet_type, EntityHandle& high_order_node ) const; 00514 00515 //! given an entity and a target dimension & side number, get that entity 00516 virtual ErrorCode side_element( const EntityHandle source_entity, const int dim, const int side_number, 00517 EntityHandle& target_entity ) const; 00518 00519 //-------------------------Tag Stuff-------------------------------------// 00520 00521 /**\brief Get a tag handle, possibly creating the tag 00522 * 00523 * Get a handle used to associate application-defined values 00524 * with MOAB entities. If the tag does not already exist then 00525 * \c flags should contain exactly one of \c MB_TAG_SPARSE, 00526 * \c MB_TAG_DENSE, \c MB_TAG_MESH unless \c type is MB_TYPE_BIT, 00527 * which implies \c MB_TAG_BIT storage. 00528 * . 00529 *\param name The tag name 00530 *\param size Tag size as number of values of of data type per entity 00531 * (or number of bytes if \c MB_TAG_BYTES is passed in flags). If \c 00532 *MB_TAG_VARLEN is specified, this value is taken to be the size of the default value if one is 00533 *specified and is otherwise ignored. \param type The type of the data (used for IO) 00534 *\param tag_handle Output: the resulting tag handle. 00535 *\param flags Bitwise OR of values from \c TagType 00536 *\param default_value Optional default value for tag. 00537 *\param created Optional returned boolean indicating that the tag 00538 * was created. 00539 *\return - \c MB_ALREADY_ALLOCATED if tag exists and \c MB_TAG_EXCL is specified, or 00540 *default values do not match (and \c MB_TAG_ANY or \c MB_TAG_DFTOK not specified). 00541 * - \c MB_TAG_NOT_FOUND if tag does not exist and \c MB_TAG_CREAT is not 00542 *specified 00543 * - \c MB_INVALID_SIZE if tag value size is not a multiple of the size of the 00544 *data type (and \c MB_TAG_ANY not specified). 00545 * - \c MB_TYPE_OUT_OF_RANGE invalid or inconsistent parameter 00546 * - \c MB_VARIABLE_DATA_LENGTH if \c MB_TAG_VARLEN and \c default_value is non-null and 00547 * \c default_value_size is not specified. 00548 */ 00549 virtual ErrorCode tag_get_handle( const char* name, int size, DataType type, Tag& tag_handle, unsigned flags = 0, 00550 const void* default_value = 0, bool* created = 0 ); 00551 00552 /**\brief same as non-const version, except that TAG_CREAT flag is ignored. */ 00553 virtual ErrorCode tag_get_handle( const char* name, int size, DataType type, Tag& tag_handle, unsigned flags = 0, 00554 const void* default_value = 0 ) const; 00555 00556 //! Gets the tag name string of the tag_handle. 00557 /** \param tag_handle Tag you want the name of. 00558 \param tag_name Name string of <em>tag_handle</em>. 00559 00560 Example: \code 00561 Tag tag_handle = 0; 00562 std::string tag_name = "my_special_tag"; 00563 tag_get_name( tag_handle, tag_name ); //gets the Tag from the tag's name string 00564 \endcode */ 00565 virtual ErrorCode tag_get_name( const Tag tag_handle, std::string& tag_name ) const; 00566 00567 /**\brief Gets the tag handle corresponding to a name 00568 * 00569 * If a tag of that name does not exist, returns MB_TAG_NOT_FOUND 00570 * \param tag_name Name of the desired tag. 00571 * \param tag_handle Tag handle corresponding to <em>tag_name</em> 00572 */ 00573 virtual ErrorCode tag_get_handle( const char* tag_name, Tag& tag_handle ) const; 00574 00575 //! Get handles for all tags defined on this entity 00576 virtual ErrorCode tag_get_tags_on_entity( const EntityHandle entity, std::vector< Tag >& tag_handles ) const; 00577 00578 //! Get the size of the specified tag in bytes 00579 virtual ErrorCode tag_get_bytes( const Tag tag, int& bytes_per_tag ) const; 00580 00581 //! Get the array length of a tag 00582 virtual ErrorCode tag_get_length( const Tag tag, int& length ) const; 00583 00584 //! Get the default value of the specified tag 00585 virtual ErrorCode tag_get_default_value( const Tag tag, void* def_val ) const; 00586 virtual ErrorCode tag_get_default_value( Tag tag, const void*& def_val, int& size ) const; 00587 00588 //! get type of tag (sparse, dense, etc.; 0 = dense, 1 = sparse, 2 = bit, 3 = mesh) 00589 virtual ErrorCode tag_get_type( const Tag, TagType& tag_type ) const; 00590 00591 /** \brief Get data type of tag. 00592 * 00593 * Get the type of the tag data. The tag is data is assumed to 00594 * be a vector of this type. If the tag data vetcor contains 00595 * more than one value, then the tag size must be a multiple of 00596 * the size of this type. 00597 * \param tag The tag 00598 * \param type The type of the specified tag (output). 00599 */ 00600 virtual ErrorCode tag_get_data_type( const Tag tag, DataType& type ) const; 00601 00602 //! get handles for all tags defined 00603 virtual ErrorCode tag_get_tags( std::vector< Tag >& tag_handles ) const; 00604 00605 virtual ErrorCode tag_get_data( const Tag tag_handle, const EntityHandle* entity_handles, const int num_entities, 00606 void* tag_data ) const; 00607 00608 virtual ErrorCode tag_get_data( const Tag tag_handle, const Range& entity_handles, void* tag_data ) const; 00609 00610 //! Sets the data of a given EntityHandle and Tag. 00611 /** If the <em>tag_handle</em> and the entity type of <em>entity_handle</em> are not 00612 compatible, data of <em>entity_handle</em> never existed and MB_FAILURE 00613 is returned. 00614 \param tag_handle Tag indicating what data is to be set. 00615 \param entity_handle EntityHandle on which to set tag's data. 00616 \param tag_data Data to set the <em>entity_handle</em>'s tag data to. 00617 00618 Example: \code 00619 int tag_data = 1004; 00620 tag_set_data( tag_handle, entity_handle, &tag_data ); \endcode */ 00621 virtual ErrorCode tag_set_data( Tag tag_handle, const EntityHandle* entity_handles, int num_entities, 00622 const void* tag_data ); 00623 00624 virtual ErrorCode tag_set_data( Tag tag_handle, const Range& entity_handles, const void* tag_data ); 00625 00626 /**\brief Get pointers to tag data 00627 * 00628 * For a tag, get the values for a list of passed entity handles. 00629 *\note This function may not be used for bit tags. 00630 *\param tag_handle The tag 00631 *\param entity_handles An array of entity handles for which to retreive tag values. 00632 *\param num_entities The length of the 'entity_handles' array. 00633 *\param tag_data An array of 'const void*'. Array must be at least 00634 * 'num_entitities' long. Array is populated (output) 00635 * with pointers to the internal storage for the 00636 * tag value corresponding to each entity handle. 00637 *\param tag_sizes The length of each tag value. Optional for 00638 * fixed-length tags. Required for variable-length tags. 00639 */ 00640 virtual ErrorCode tag_get_by_ptr( const Tag tag_handle, const EntityHandle* entity_handles, int num_entities, 00641 const void** tag_data, int* tag_sizes = 0 ) const; 00642 00643 /**\brief Get pointers to tag data 00644 * 00645 * For a tag, get the values for a list of passed entity handles. 00646 *\note This function may not be used for bit tags. 00647 *\param tag_handle The tag 00648 *\param entity_handles The entity handles for which to retreive tag values. 00649 *\param tag_data An array of 'const void*'. Array is populated (output) 00650 * with pointers to the internal storage for the 00651 * tag value corresponding to each entity handle. 00652 *\param tag_sizes The length of each tag value. Optional for 00653 * fixed-length tags. Required for variable-length tags. 00654 */ 00655 virtual ErrorCode tag_get_by_ptr( const Tag tag_handle, const Range& entity_handles, const void** tag_data, 00656 int* tag_sizes = 0 ) const; 00657 00658 /**\brief Set tag data given an array of pointers to tag values. 00659 * 00660 * For a tag, set the values for a list of passed entity handles. 00661 *\note This function may not be used for bit tags. 00662 *\param tag_handle The tag 00663 *\param entity_handles An array of entity handles for which to set tag values. 00664 *\param num_entities The length of the 'entity_handles' array. 00665 *\param tag_data An array of 'const void*'. Array must be at least 00666 * 'num_entitities' long. Array is expected to 00667 * contain pointers to tag values for the corresponding 00668 * EntityHandle in 'entity_handles'. 00669 *\param tag_sizes The length of each tag value. Optional for 00670 * fixed-length tags. Required for variable-length tags. 00671 */ 00672 virtual ErrorCode tag_set_by_ptr( Tag tag_handle, const EntityHandle* entity_handles, int num_entities, 00673 void const* const* tag_data, const int* tag_sizes = 0 ); 00674 00675 /**\brief Set tag data given an array of pointers to tag values. 00676 * 00677 * For a tag, set the values for a list of passed entity handles. 00678 *\note This function may not be used for bit tags. 00679 *\param tag_handle The tag 00680 *\param entity_handles The entity handles for which to set tag values. 00681 *\param tag_data An array of 'const void*'. Array is expected to 00682 * contain pointers to tag values for the corresponding 00683 * EntityHandle in 'entity_handles'. 00684 *\param tag_sizes The length of each tag value. Optional for 00685 * fixed-length tags. Required for variable-length tags. 00686 */ 00687 virtual ErrorCode tag_set_by_ptr( Tag tag_handle, const Range& entity_handles, void const* const* tag_data, 00688 const int* tag_sizes = 0 ); 00689 00690 /**\brief Set tag data given value. 00691 * 00692 * For a tag, set the values for a list of passed entity handles to 00693 * the same, specified value. 00694 * 00695 *\param tag_handle The tag 00696 *\param entity_handles The entity handles for which to set tag values. 00697 *\param tag_data A pointer to the tag value. 00698 *\param tag_sizes For variable-length tags, the lenght of the 00699 * tag value. This argument will be ignored for 00700 * fixed-length tags. 00701 */ 00702 virtual ErrorCode tag_clear_data( Tag tag_handle, const Range& entity_handles, const void* value, 00703 int value_size = 0 ); 00704 00705 /**\brief Set tag data given value. 00706 * 00707 * For a tag, set the values for a list of passed entity handles to 00708 * the same, specified value. 00709 * 00710 *\param tag_handle The tag 00711 *\param entity_handles The entity handles for which to set tag values. 00712 *\param tag_data A pointer to the tag value. 00713 *\param tag_sizes For variable-length tags, the lenght of the 00714 * tag value. This argument will be ignored for 00715 * fixed-length tags. 00716 */ 00717 virtual ErrorCode tag_clear_data( Tag tag_handle, const EntityHandle* entity_handles, int num_entity_handles, 00718 const void* value, int value_size = 0 ); 00719 00720 //! Delete the data of a vector of entity handles and sparse tag 00721 /** Delete the data of a tag on a vector of entity handles. Only sparse tag data are deleted 00722 with this function; dense tags are deleted by deleting the tag itself using tag_delete. 00723 \param tag_handle Handle of the (sparse) tag being deleted from entity 00724 \param entity_handles 1d vector of entity handles from which the tag is being deleted 00725 \param num_handles Number of entity handles in 1d vector 00726 */ 00727 virtual ErrorCode tag_delete_data( Tag tag_handle, const EntityHandle* entity_handles, int num_handles ); 00728 00729 //! Delete the data of a range of entity handles and sparse tag 00730 /** Delete the data of a tag on a range of entity handles. Only sparse tag data are deleted 00731 with this function; dense tags are deleted by deleting the tag itself using tag_delete. 00732 \param tag_handle Handle of the (sparse) tag being deleted from entity 00733 \param entity_range Range of entities from which the tag is being deleted 00734 */ 00735 virtual ErrorCode tag_delete_data( Tag tag_handle, const Range& entity_range ); 00736 00737 //! Removes the tag from the database and deletes all of its associated data. 00738 virtual ErrorCode tag_delete( Tag tag_handle ); 00739 00740 /**\brief Access tag data via direct pointer into contiguous blocks 00741 * 00742 * Iteratively obtain direct access to contiguous blocks of tag 00743 * storage. This function cannot be used with bit tags because 00744 * of the compressed bit storage. This function cannot be used 00745 * with variable length tags because it does not provide a mechanism 00746 * to determine the length of the value for each entity. This 00747 * function may be used with sparse tags, but if it is used, it 00748 * will return data for a single entity at a time. 00749 * 00750 *\param tag_handle The handle of the tag for which to access data 00751 *\param iter The first entity for which to return data. 00752 *\param end One past the last entity for which data is desired. 00753 *\param count The number of entities for which data was returned 00754 *\param data_ptr Output: pointer to tag storage. 00755 *\param allocate If true, space for this tag will be allocated, if not it wont 00756 * 00757 *\Note If this function is called for entities for which no tag value 00758 * has been set, but for which a default value exists, it will 00759 * force the allocation of explicit storage for each such entity 00760 * even though MOAB would normally not explicitly store tag values 00761 * for such entities. 00762 * 00763 *\Example: 00764 *\code 00765 * Range ents; // range to iterate over 00766 * Tag tag; // tag for which to access data 00767 * int bytes; 00768 * ErrorCode err = mb.tag_get_size( tag, bytes ); 00769 * if (err) { ... } 00770 * 00771 * ... 00772 * Range::iterator iter = ents.begin(); 00773 * while (iter != ents.end()) { 00774 * int count; 00775 * // get contiguous block of tag dat 00776 * void* ptr; 00777 * err = mb.tag_iterate( tag, iter, ents.end(), count, ptr ); 00778 * if (err) { ... } 00779 * // do something with tag data 00780 * process_Data( ptr, count ); 00781 * // advance to next block of data 00782 * iter += count; 00783 * } 00784 *\endcode 00785 */ 00786 virtual ErrorCode tag_iterate( Tag tag_handle, Range::const_iterator begin, Range::const_iterator end, int& count, 00787 void*& data_ptr, bool allocate = true ); 00788 00789 //! creates a mesh set 00790 virtual ErrorCode create_meshset( const unsigned int options, EntityHandle& ms_handle, int start_id = 0 ); 00791 00792 //! Empty a vector of mesh set 00793 /** Empty a mesh set. 00794 \param ms_handles 1d vector of handles of sets being emptied 00795 \param num_meshsets Number of entities in 1d vector 00796 */ 00797 virtual ErrorCode clear_meshset( const EntityHandle* ms_handles, const int num_meshsets ); 00798 00799 //! Empty a range of mesh set 00800 /** Empty a mesh set. 00801 \param ms_handles Range of handles of sets being emptied 00802 */ 00803 virtual ErrorCode clear_meshset( const Range& ms_handles ); 00804 00805 //! get the options of a mesh set 00806 virtual ErrorCode get_meshset_options( const EntityHandle ms_handle, unsigned int& options ) const; 00807 00808 //! set the options of a mesh set 00809 virtual ErrorCode set_meshset_options( const EntityHandle ms_handle, const unsigned int options ); 00810 00811 //! subtracts meshset2 from meshset1 - modifies meshset1 00812 virtual ErrorCode subtract_meshset( EntityHandle meshset1, const EntityHandle meshset2 ); 00813 00814 //! intersects meshset2 with meshset1 - modifies meshset1 00815 virtual ErrorCode intersect_meshset( EntityHandle meshset1, const EntityHandle meshset2 ); 00816 00817 //! unites meshset2 with meshset1 - modifies meshset1 00818 virtual ErrorCode unite_meshset( EntityHandle meshset1, const EntityHandle meshset2 ); 00819 00820 //! add entities to meshset 00821 virtual ErrorCode add_entities( EntityHandle meshset, const Range& entities ); 00822 00823 //! add entities to meshset 00824 virtual ErrorCode add_entities( EntityHandle meshset, const EntityHandle* entities, const int num_entities ); 00825 00826 //! remove entities from meshset 00827 virtual ErrorCode remove_entities( EntityHandle meshset, const Range& entities ); 00828 00829 //! remove entities from meshset 00830 virtual ErrorCode remove_entities( EntityHandle meshset, const EntityHandle* entities, const int num_entities ); 00831 00832 //! return true if all entities are contained in set 00833 virtual bool contains_entities( EntityHandle meshset, const EntityHandle* entities, int num_entities, 00834 const int operation_type = Interface::INTERSECT ); 00835 00836 //! replace entities 00837 virtual ErrorCode replace_entities( EntityHandle meshset, const EntityHandle* old_entities, 00838 const EntityHandle* new_entities, int num_entities ); 00839 00840 //------MeshSet Parent/Child functions------ 00841 00842 //! get parent meshsets 00843 virtual ErrorCode get_parent_meshsets( const EntityHandle meshset, std::vector< EntityHandle >& parents, 00844 const int num_hops = 1 ) const; 00845 00846 //! get parent meshsets 00847 virtual ErrorCode get_parent_meshsets( const EntityHandle meshset, Range& parents, const int num_hops = 1 ) const; 00848 00849 //! get child meshsets 00850 virtual ErrorCode get_child_meshsets( const EntityHandle meshset, std::vector< EntityHandle >& children, 00851 const int num_hops = 1 ) const; 00852 00853 //! get child meshsets 00854 virtual ErrorCode get_child_meshsets( const EntityHandle meshset, Range& children, const int num_hops = 1 ) const; 00855 00856 //! get contained meshsets 00857 virtual ErrorCode get_contained_meshsets( const EntityHandle meshset, std::vector< EntityHandle >& contained, 00858 const int num_hops = 1 ) const; 00859 00860 //! get contained meshsets 00861 virtual ErrorCode get_contained_meshsets( const EntityHandle meshset, Range& contained, 00862 const int num_hops = 1 ) const; 00863 00864 //! gets number of parent meshsets 00865 virtual ErrorCode num_parent_meshsets( const EntityHandle meshset, int* number, const int num_hops = 1 ) const; 00866 00867 //! gets number of child meshsets 00868 virtual ErrorCode num_child_meshsets( const EntityHandle meshset, int* number, const int num_hops = 1 ) const; 00869 00870 //! gets number of contained meshsets 00871 virtual ErrorCode num_contained_meshsets( const EntityHandle meshset, int* number, const int num_hops = 1 ) const; 00872 00873 //! add a parent meshset 00874 virtual ErrorCode add_parent_meshset( EntityHandle meshset, const EntityHandle parent_meshset ); 00875 00876 //! add parent meshsets 00877 virtual ErrorCode add_parent_meshsets( EntityHandle meshset, const EntityHandle* parent_meshsets, 00878 int num_parent_meshsets ); 00879 00880 //! add a child meshset 00881 virtual ErrorCode add_child_meshset( EntityHandle meshset, const EntityHandle child_meshset ); 00882 00883 //! add parent meshsets 00884 virtual ErrorCode add_child_meshsets( EntityHandle meshset, const EntityHandle* child_meshsets, 00885 int num_child_meshsets ); 00886 00887 //! adds 'parent' to child's parent list and adds 'child' to parent's child list 00888 virtual ErrorCode add_parent_child( EntityHandle parent, EntityHandle child ); 00889 00890 //! removes 'parent' to child's parent list and removes 'child' to parent's child list 00891 virtual ErrorCode remove_parent_child( EntityHandle parent, EntityHandle child ); 00892 00893 //! remove parent meshset 00894 virtual ErrorCode remove_parent_meshset( EntityHandle meshset, const EntityHandle parent_meshset ); 00895 00896 //! remove child meshset 00897 virtual ErrorCode remove_child_meshset( EntityHandle meshset, const EntityHandle child_meshset ); 00898 00899 // ************************ tag information *************** 00900 00901 //! return various specific tag handles 00902 Tag material_tag(); 00903 Tag neumannBC_tag(); 00904 Tag dirichletBC_tag(); 00905 Tag globalId_tag(); 00906 Tag geom_dimension_tag(); 00907 00908 //! get/set the number of nodes 00909 // int total_num_nodes() const; 00910 // void total_num_nodes(const int val); 00911 00912 //! get/set the number of elements 00913 // int total_num_elements() const; 00914 // void total_num_elements(const int val); 00915 00916 // ************************ structured sequence information *************** 00917 00918 //! return a reference to the sequence manager 00919 SequenceManager* sequence_manager() 00920 { 00921 return sequenceManager; 00922 } 00923 const SequenceManager* sequence_manager() const 00924 { 00925 return sequenceManager; 00926 } 00927 00928 /// create structured sequence 00929 ErrorCode create_scd_sequence( const HomCoord& coord_min, const HomCoord& coord_max, EntityType type, 00930 EntityID start_id_hint, EntityHandle& first_handle_out, 00931 EntitySequence*& sequence_out ); 00932 00933 ErrorCode add_vsequence( EntitySequence* vert_seq, EntitySequence* elem_seq, const HomCoord& p1, const HomCoord& q1, 00934 const HomCoord& p2, const HomCoord& q2, const HomCoord& p3, const HomCoord& q3, 00935 bool bb_input = false, const HomCoord* bb_min = NULL, const HomCoord* bb_max = NULL ); 00936 00937 //! return the a_entity_factory pointer 00938 AEntityFactory* a_entity_factory() 00939 { 00940 return aEntityFactory; 00941 } 00942 const AEntityFactory* a_entity_factory() const 00943 { 00944 return aEntityFactory; 00945 } 00946 00947 #ifdef MOAB_HAVE_AHF 00948 HalfFacetRep* a_half_facet_rep() 00949 { 00950 return ahfRep; 00951 } 00952 const HalfFacetRep* a_half_facet_rep() const 00953 { 00954 return ahfRep; 00955 } 00956 #endif 00957 00958 //! return set of registered IO tools 00959 ReaderWriterSet* reader_writer_set() 00960 { 00961 return readerWriterSet; 00962 } 00963 00964 //-----------------MeshSet Interface Functions------------------// 00965 00966 void print( const EntityHandle handle, const char* prefix, bool first_call = true ) const; 00967 00968 ErrorCode print_entity_tags( std::string indent_prefix, const EntityHandle handle, TagType tp ) const; 00969 00970 virtual ErrorCode get_last_error( std::string& info ) const; 00971 00972 virtual std::string get_error_string( const ErrorCode code ) const; 00973 00974 //! check all adjacencies for consistency 00975 ErrorCode check_adjacencies(); 00976 00977 //! check some adjacencies for consistency 00978 ErrorCode check_adjacencies( const EntityHandle* ents, int num_ents ); 00979 00980 //! return whether the input handle is valid or not 00981 bool is_valid( const EntityHandle this_ent ) const; 00982 00983 /** \brief Create an iterator over the set 00984 * Create a new iterator that iterates over entities with the specified type or dimension. 00985 * Only one of ent_type or dim can be set; use dim=-1 or ent_type=MBMAXTYPE for the other. 00986 * Iterators for list-type (ordered) sets are stable over set modification, unless entity 00987 * removed or deleted is the one at the current position of the iterator. If the check_valid 00988 * parameter is passed as true, entities are checked for validity before being passed back by 00989 * get_next_entities function (checking entity validity can have a non-negligible cost). 00990 * 00991 * Iterators returned by this function can be deleted using the normal C++ delete function. 00992 * After creating the iterator through this function, further interactions are through methods 00993 * on the SetIterator class. 00994 * \param meshset The entity set associated with this iterator (use 0 for whole instance) 00995 * \param ent_type Entity type associated with this iterator 00996 * \param ent_dim Dimension associated with this iterator 00997 * \param chunk_size Chunk size of the iterator 00998 * \param check_valid If true, entities are checked for validity before being returned 00999 */ 01000 virtual ErrorCode create_set_iterator( EntityHandle meshset, EntityType ent_type, int ent_dim, int chunk_size, 01001 bool check_valid, SetIterator*& set_iter ); 01002 01003 /** \brief Remove the set iterator from the instance's list 01004 * \param set_iter Set iterator being removed 01005 */ 01006 ErrorCode remove_set_iterator( SetIterator* set_iter ); 01007 01008 /** \brief Get all set iterators associated with the set passed in 01009 * \param meshset Meshset for which iterators are requested 01010 * \param set_iters Set iterators for the set 01011 */ 01012 ErrorCode get_set_iterators( EntityHandle meshset, std::vector< SetIterator* >& set_iters ); 01013 01014 //-----------------Memory Functions------------------// 01015 01016 /**\brief Calculate amount of memory used to store MOAB data 01017 * 01018 * This function calculates the amount of memory used to store 01019 * MOAB data. 01020 * 01021 * There are two possible values for each catagory of memory use. 01022 * The exact value and the amortized value. The exact value is the 01023 * amount of memory used to store the data for the specified entities. 01024 * The amortized value includes the exact value and an amortized 01025 * estimate of the memory consumed in overhead for storing the values 01026 * (indexing structures, access structures, etc.) 01027 * 01028 * Note: If ent_array is NULL, the total memory used by MOAB for storing 01029 * data will be returned in the address pointed to by 01030 * total_amortized_storage, if total_amortized_storage is not NULL. 01031 * 01032 *\param ent_array Array of entities for which to estimate the memory use. 01033 * If NULL, estimate is done for all entities. 01034 *\param num_ents The length of ent_array. Not used if ent_rray is NULL. 01035 *\param total_(amortized_)storage The sum of the memory entity, adjacency, and all tag storage. 01036 *\param (amortized_)entity_storage The storage for the entity definitions 01037 * (connectivity arrays for elements, coordinates for vertices, 01038 * list storage within sets, etc.) 01039 *\param (amortized_)adjacency_storage The storage for adjacency data. 01040 *\param tag_array An array of tags for which to calculate the memory use. 01041 *\param num_tags The lenght of tag_array 01042 *\param (amortized_)tag_storage If tag_array is not NULL, then one value 01043 * for each tag specifying the memory used for storing that 01044 * tag. If tag_array is NULL and this value is not, the 01045 * location at which to store the total memory used for 01046 * all tags. 01047 */ 01048 void estimated_memory_use( const EntityHandle* ent_array = 0, unsigned long num_ents = 0, 01049 unsigned long long* total_storage = 0, unsigned long long* total_amortized_storage = 0, 01050 unsigned long long* entity_storage = 0, unsigned long long* amortized_entity_storage = 0, 01051 unsigned long long* adjacency_storage = 0, 01052 unsigned long long* amortized_adjacency_storage = 0, const Tag* tag_array = 0, 01053 unsigned num_tags = 0, unsigned long long* tag_storage = 0, 01054 unsigned long long* amortized_tag_storage = 0 ); 01055 01056 /**\brief Calculate amount of memory used to store MOAB data 01057 * 01058 * This function calculates the amount of memory used to store 01059 * MOAB data. 01060 * 01061 * There are two possible values for each catagory of memory use. 01062 * The exact value and the amortized value. The exact value is the 01063 * amount of memory used to store the data for the specified entities. 01064 * The amortized value includes the exact value and an amortized 01065 * estimate of the memory consumed in overhead for storing the values 01066 * (indexing structures, access structures, etc.) 01067 * 01068 *\param ents Entities for which to estimate the memory use. 01069 *\param total_(amortized_)storage The sum of the memory entity, adjacency, and all tag storage. 01070 *\param (amortized_)entity_storage The storage for the entity definitions 01071 * (connectivity arrays for elements, coordinates for vertices, 01072 * list storage within sets, etc.) 01073 *\param (amortized_)adjacency_storage The storage for adjacency data. 01074 *\param tag_array An array of tags for which to calculate the memory use. 01075 *\param num_tags The lenght of tag_array 01076 *\param (amortized_)tag_storage If tag_array is not NULL, then one value 01077 * for each tag specifying the memory used for storing that 01078 * tag. If tag_array is NULL and this value is not, the 01079 * location at which to store the total memory used for 01080 * all tags. 01081 */ 01082 void estimated_memory_use( const Range& ents, unsigned long long* total_storage = 0, 01083 unsigned long long* total_amortized_storage = 0, unsigned long long* entity_storage = 0, 01084 unsigned long long* amortized_entity_storage = 0, 01085 unsigned long long* adjacency_storage = 0, 01086 unsigned long long* amortized_adjacency_storage = 0, const Tag* tag_array = 0, 01087 unsigned num_tags = 0, unsigned long long* tag_storage = 0, 01088 unsigned long long* amortized_tag_storage = 0 ); 01089 01090 void print_database() const; 01091 01092 /** \name Sequence Option controllers */ 01093 01094 /**@{*/ 01095 01096 /** \brief Interface to control memory allocation for sequences 01097 * Provide a factor that controls the size of the sequence that gets allocated. 01098 * This is typically useful in the parallel setting when a-priori, the number of ghost entities 01099 * and the memory required for them within the same sequence as the owned entities are unknown. 01100 * The default factor is 1.0 but this can be appropriately updated at runtime so that we do not 01101 * have broken sequences. 01102 */ 01103 virtual double get_sequence_multiplier() const; 01104 01105 /** \brief Interface to control memory allocation for sequences 01106 * Provide a factor that controls the size of the sequence that gets allocated. 01107 * This is typically useful in the parallel setting when a-priori, the number of ghost entities 01108 * and the memory required for them within the same sequence as the owned entities are unknown. 01109 * The default factor is 1.0 but this can be appropriately updated at runtime so that we do not 01110 * have broken sequences. 01111 * 01112 * \param meshset User specified multiplier (should be greater than 1.0) 01113 */ 01114 virtual void set_sequence_multiplier( double factor ); 01115 01116 /**@}*/ 01117 01118 private: 01119 /**\brief Do not allow copying */ 01120 Core( const Core& copy ); 01121 /**\brief Do not allow copying */ 01122 Core& operator=( const Core& copy ); 01123 01124 void estimated_memory_use_internal( const Range* ents, unsigned long long* total_storage, 01125 unsigned long long* total_amortized_storage, unsigned long long* entity_storage, 01126 unsigned long long* amortized_entity_storage, 01127 unsigned long long* adjacency_storage, 01128 unsigned long long* amortized_adjacency_storage, const Tag* tag_array, 01129 unsigned num_tags, unsigned long long* tag_storage, 01130 unsigned long long* amortized_tag_storage ); 01131 01132 //! database init and de-init routines 01133 ErrorCode initialize(); 01134 void deinitialize(); 01135 01136 //! return the entity set representing the whole mesh 01137 EntityHandle get_root_set(); 01138 01139 //!\brief Clean up after a file reader returns failure. 01140 //! 01141 //! Delete all entities not contained in initial_entities 01142 //! and all tags not contained in initial_tags. 01143 void clean_up_failed_read( const Range& initial_entities, std::vector< Tag > initial_tags ); 01144 01145 // other interfaces for MB 01146 WriteUtil* mMBWriteUtil; 01147 ReadUtil* mMBReadUtil; 01148 ScdInterface* scdInterface; 01149 01150 //! store the total number of elements defined in this interface 01151 // int totalNumElements; 01152 01153 //! store the total number of nodes defined in this interface 01154 // int totalNumNodes; 01155 01156 //! the overall geometric dimension of this mesh 01157 int geometricDimension; 01158 01159 Tag materialTag; 01160 Tag neumannBCTag; 01161 Tag dirichletBCTag; 01162 Tag geomDimensionTag; 01163 Tag globalIdTag; 01164 01165 //! tag server for this interface 01166 std::list< TagInfo* > tagList; 01167 inline bool valid_tag_handle( const TagInfo* t ) const 01168 { 01169 return std::find( tagList.begin(), tagList.end(), t ) != tagList.end(); 01170 } 01171 01172 SequenceManager* sequenceManager; 01173 01174 AEntityFactory* aEntityFactory; 01175 01176 ReaderWriterSet* readerWriterSet; 01177 01178 Error* mError; 01179 bool mpiFinalize; 01180 int writeMPELog; 01181 bool initErrorHandlerInCore; 01182 01183 //! list of iterators 01184 std::vector< SetIterator* > setIterators; 01185 01186 #ifdef MOAB_HAVE_AHF 01187 HalfFacetRep* ahfRep; 01188 bool mesh_modified; 01189 #endif 01190 }; 01191 01192 } // namespace moab 01193 01194 #endif // MOAB_IMPL_GENERAL_HPP