Mesh Oriented datABase  (version 5.5.0)
An array-based unstructured mesh library
iMOAB.cpp File Reference
#include "moab/MOABConfig.h"
#include "moab/Core.hpp"
#include "DebugOutput.hpp"
#include "moab/iMOAB.h"
#include "moab/CartVect.hpp"
#include "MBTagConventions.hpp"
#include "moab/MeshTopoUtil.hpp"
#include "moab/ReadUtilIface.hpp"
#include "moab/MergeMesh.hpp"
#include <cassert>
#include <sstream>
#include <iostream>
+ Include dependency graph for iMOAB.cpp:

Go to the source code of this file.

Classes

struct  appData
 
struct  GlobalContext
 

Functions

ErrCode iMOAB_Initialize (int argc, iMOAB_String *argv)
 Initialize the iMOAB interface implementation. More...
 
ErrCode iMOAB_InitializeFortran ()
 Initialize the iMOAB interface implementation from Fortran driver. More...
 
ErrCode iMOAB_Finalize ()
 Finalize the iMOAB interface implementation. More...
 
ErrCode iMOAB_RegisterApplication (const iMOAB_String app_name, int *compid, iMOAB_AppID pid)
 Register application - Create a unique application ID and bootstrap interfaces for further queries. More...
 
ErrCode iMOAB_RegisterApplicationFortran (const iMOAB_String app_name, int *compid, iMOAB_AppID pid)
 Register a Fortran-based application - Create a unique application ID and bootstrap interfaces for further queries. More...
 
ErrCode iMOAB_DeregisterApplication (iMOAB_AppID pid)
 De-Register application: delete mesh (set) associated with the application ID. More...
 
ErrCode iMOAB_DeregisterApplicationFortran (iMOAB_AppID pid)
 De-Register the Fortran based application: delete mesh (set) associated with the application ID. More...
 
static void split_tag_names (std::string input_names, std::string &separator, std::vector< std::string > &list_tag_names)
 
ErrCode iMOAB_ReadHeaderInfo (const iMOAB_String filename, int *num_global_vertices, int *num_global_elements, int *num_dimension, int *num_parts)
 It should be called on master task only, and information obtained could be broadcasted by the user. It is a fast lookup in the header of the file. More...
 
ErrCode iMOAB_LoadMesh (iMOAB_AppID pid, const iMOAB_String filename, const iMOAB_String read_options, int *num_ghost_layers)
 Load a MOAB mesh file in parallel and exchange ghost layers as requested. More...
 
ErrCode iMOAB_WriteMesh (iMOAB_AppID pid, const iMOAB_String filename, const iMOAB_String write_options)
 Write a MOAB mesh along with the solution tags to a file. More...
 
ErrCode iMOAB_WriteLocalMesh (iMOAB_AppID pid, iMOAB_String prefix)
 Write a local MOAB mesh copy. More...
 
ErrCode iMOAB_UpdateMeshInfo (iMOAB_AppID pid)
 Update local mesh data structure, from file information. More...
 
ErrCode iMOAB_GetMeshInfo (iMOAB_AppID pid, int *num_visible_vertices, int *num_visible_elements, int *num_visible_blocks, int *num_visible_surfaceBC, int *num_visible_vertexBC)
 Retrieve all the important mesh information related to vertices, elements, vertex and surface boundary conditions. More...
 
ErrCode iMOAB_GetVertexID (iMOAB_AppID pid, int *vertices_length, iMOAB_GlobalID *global_vertex_ID)
 Get the global vertex IDs for all locally visible (owned and shared/ghosted) vertices. More...
 
ErrCode iMOAB_GetVertexOwnership (iMOAB_AppID pid, int *vertices_length, int *visible_global_rank_ID)
 Get vertex ownership information. More...
 
ErrCode iMOAB_GetVisibleVerticesCoordinates (iMOAB_AppID pid, int *coords_length, double *coordinates)
 Get vertex coordinates for all local (owned and ghosted) vertices. More...
 
ErrCode iMOAB_GetBlockID (iMOAB_AppID pid, int *block_length, iMOAB_GlobalID *global_block_IDs)
 Get the global block IDs for all locally visible (owned and shared/ghosted) blocks. More...
 
ErrCode iMOAB_GetBlockInfo (iMOAB_AppID pid, iMOAB_GlobalID *global_block_ID, int *vertices_per_element, int *num_elements_in_block)
 Get the global block information and number of visible elements of belonging to a block (MATERIAL SET). More...
 
ErrCode iMOAB_GetVisibleElementsInfo (iMOAB_AppID pid, int *num_visible_elements, iMOAB_GlobalID *element_global_IDs, int *ranks, iMOAB_GlobalID *block_IDs)
 Get the visible elements information. More...
 
ErrCode iMOAB_GetBlockElementConnectivities (iMOAB_AppID pid, iMOAB_GlobalID *global_block_ID, int *connectivity_length, int *element_connectivity)
 Get the connectivities for elements contained within a certain block. More...
 
ErrCode iMOAB_GetElementConnectivity (iMOAB_AppID pid, iMOAB_LocalID *elem_index, int *connectivity_length, int *element_connectivity)
 Get the connectivity for one element only. More...
 
ErrCode iMOAB_GetElementOwnership (iMOAB_AppID pid, iMOAB_GlobalID *global_block_ID, int *num_elements_in_block, int *element_ownership)
 Get the element ownership within a certain block i.e., processor ID of the element owner. More...
 
ErrCode iMOAB_GetElementID (iMOAB_AppID pid, iMOAB_GlobalID *global_block_ID, int *num_elements_in_block, iMOAB_GlobalID *global_element_ID, iMOAB_LocalID *local_element_ID)
 Get the global IDs for all locally visible elements belonging to a particular block. More...
 
ErrCode iMOAB_GetPointerToSurfaceBC (iMOAB_AppID pid, int *surface_BC_length, iMOAB_LocalID *local_element_ID, int *reference_surface_ID, int *boundary_condition_value)
 Get the surface boundary condition information. More...
 
ErrCode iMOAB_GetPointerToVertexBC (iMOAB_AppID pid, int *vertex_BC_length, iMOAB_LocalID *local_vertex_ID, int *boundary_condition_value)
 Get the vertex boundary condition information. More...
 
ErrCode iMOAB_DefineTagStorage (iMOAB_AppID pid, const iMOAB_String tag_storage_name, int *tag_type, int *components_per_entity, int *tag_index)
 Define a MOAB Tag corresponding to the application depending on requested types. More...
 
ErrCode iMOAB_SetIntTagStorage (iMOAB_AppID pid, const iMOAB_String tag_storage_name, int *num_tag_storage_length, int *ent_type, int *tag_storage_data)
 Store the specified values in a MOAB integer Tag. More...
 
ErrCode iMOAB_GetIntTagStorage (iMOAB_AppID pid, const iMOAB_String tag_storage_name, int *num_tag_storage_length, int *ent_type, int *tag_storage_data)
 Get the specified values in a MOAB integer Tag. More...
 
ErrCode iMOAB_SetDoubleTagStorage (iMOAB_AppID pid, const iMOAB_String tag_storage_names, int *num_tag_storage_length, int *ent_type, double *tag_storage_data)
 Store the specified values in a MOAB double Tag. More...
 
ErrCode iMOAB_SetDoubleTagStorageWithGid (iMOAB_AppID pid, const iMOAB_String tag_storage_names, int *num_tag_storage_length, int *ent_type, double *tag_storage_data, int *globalIds)
 Store the specified values in a MOAB double Tag, for given ids. More...
 
ErrCode iMOAB_GetDoubleTagStorage (iMOAB_AppID pid, const iMOAB_String tag_storage_names, int *num_tag_storage_length, int *ent_type, double *tag_storage_data)
 Retrieve the specified values in a MOAB double Tag. More...
 
ErrCode iMOAB_SynchronizeTags (iMOAB_AppID pid, int *num_tag, int *tag_indices, int *ent_type)
 Exchange tag values for the given tags across process boundaries. More...
 
ErrCode iMOAB_ReduceTagsMax (iMOAB_AppID pid, int *tag_index, int *ent_type)
 Perform global reductions with the processes in the current applications communicator. Specifically this routine performs reductions on the maximum value of the tag indicated by tag_index. More...
 
ErrCode iMOAB_GetNeighborElements (iMOAB_AppID pid, iMOAB_LocalID *local_index, int *num_adjacent_elements, iMOAB_LocalID *adjacent_element_IDs)
 Retrieve the adjacencies for the element entities. More...
 
ErrCode iMOAB_CreateVertices (iMOAB_AppID pid, int *coords_len, int *dim, double *coordinates)
 Create vertices for an app; it assumes no other vertices. More...
 
ErrCode iMOAB_CreateElements (iMOAB_AppID pid, int *num_elem, int *type, int *num_nodes_per_element, int *connectivity, int *block_ID)
 Create elements for an app; it assumes connectivity from local vertices, in order. More...
 
ErrCode iMOAB_SetGlobalInfo (iMOAB_AppID pid, int *num_global_verts, int *num_global_elems)
 Set global information for number of vertices and number of elements; It is usually available from the h5m file, or it can be computed with a MPI_Reduce. More...
 
ErrCode iMOAB_GetGlobalInfo (iMOAB_AppID pid, int *num_global_verts, int *num_global_elems)
 Get global information about number of vertices and number of elements. More...
 

Variables

static struct GlobalContext context
 

Function Documentation

◆ iMOAB_CreateElements()

ErrCode iMOAB_CreateElements ( iMOAB_AppID  pid,
int *  num_elem,
int *  type,
int *  num_nodes_per_element,
int *  connectivity,
int *  block_ID 
)

Create elements for an app; it assumes connectivity from local vertices, in order.

Note
Operations: Not collective
Parameters
[in]pid(iMOAB_AppID) The unique pointer to the application ID.
[in]num_elem(int*) Number of elements.
[in]type(int*) Type of element (moab type).
[in]num_nodes_per_element(int*) Number of nodes per element.
[in]connectivity(int *) Connectivity array, with respect to vertices (assumed contiguous).
[in]block_ID(int *) The local block identifier, which will now contain the elements.
Returns
ErrCode The error code indicating success or failure.

add the new ents to the clock set

Definition at line 2225 of file iMOAB.cpp.

2231 {
2232  // Create elements
2233  appData& data = context.appDatas[*pid];
2234 
2235  ReadUtilIface* read_iface;
2236  ErrorCode rval = context.MBI->query_interface( read_iface );MB_CHK_ERR( rval );
2237 
2238  EntityType mbtype = (EntityType)( *type );
2239  EntityHandle actual_start_handle;
2240  EntityHandle* array = NULL;
2241  rval = read_iface->get_element_connect( *num_elem, *num_nodes_per_element, mbtype, 1, actual_start_handle, array );MB_CHK_ERR( rval );
2242 
2243  // fill up with actual connectivity from input; assume the vertices are in order, and start
2244  // vertex is the first in the current data vertex range
2245  EntityHandle firstVertex = data.local_verts[0];
2246 
2247  for( int j = 0; j < *num_elem * ( *num_nodes_per_element ); j++ )
2248  {
2249  array[j] = connectivity[j] + firstVertex - 1;
2250  } // assumes connectivity uses 1 based array (from fortran, mostly)
2251 
2252  Range new_elems( actual_start_handle, actual_start_handle + *num_elem - 1 );
2253 
2254  rval = context.MBI->add_entities( data.file_set, new_elems );MB_CHK_ERR( rval );
2255 
2256  data.primary_elems.merge( new_elems );
2257 
2258  // add to adjacency
2259  rval = read_iface->update_adjacencies( actual_start_handle, *num_elem, *num_nodes_per_element, array );MB_CHK_ERR( rval );
2260  // organize all new elements in block, with the given block ID; if the block set is not
2261  // existing, create a new mesh set;
2262  Range sets;
2263  int set_no = *block_ID;
2264  const void* setno_ptr = &set_no;
2266  sets );
2267  EntityHandle block_set;
2268 
2269  if( MB_FAILURE == rval || sets.empty() )
2270  {
2271  // create a new set, with this block ID
2272  rval = context.MBI->create_meshset( MESHSET_SET, block_set );MB_CHK_ERR( rval );
2273 
2274  rval = context.MBI->tag_set_data( context.material_tag, &block_set, 1, &set_no );MB_CHK_ERR( rval );
2275 
2276  // add the material set to file set
2277  rval = context.MBI->add_entities( data.file_set, &block_set, 1 );MB_CHK_ERR( rval );
2278  }
2279  else
2280  {
2281  block_set = sets[0];
2282  } // first set is the one we want
2283 
2284  /// add the new ents to the clock set
2285  rval = context.MBI->add_entities( block_set, new_elems );MB_CHK_ERR( rval );
2286 
2287  return moab::MB_SUCCESS;
2288 }

References moab::Interface::add_entities(), GlobalContext::appDatas, context, moab::Interface::create_meshset(), moab::Range::empty(), ErrorCode, appData::file_set, moab::ReadUtilIface::get_element_connect(), moab::Interface::get_entities_by_type_and_tag(), appData::local_verts, GlobalContext::material_tag, MB_CHK_ERR, MB_SUCCESS, MBENTITYSET, GlobalContext::MBI, moab::Range::merge(), MESHSET_SET, appData::primary_elems, moab::Interface::query_interface(), moab::Interface::tag_set_data(), and moab::ReadUtilIface::update_adjacencies().

◆ iMOAB_CreateVertices()

ErrCode iMOAB_CreateVertices ( iMOAB_AppID  pid,
int *  coords_len,
int *  dim,
double *  coordinates 
)

Create vertices for an app; it assumes no other vertices.

Note
Operations: Not collective
Parameters
[in]pid(iMOAB_AppID) The unique pointer to the application ID.
[in]coords_len(int*) Size of the coordinates array (nverts * dim).
[in]dim(int*) Dimension of the vertex coordinates (usually 3).
[in]coordinates(double*) Coordinates of all vertices, interleaved.
Returns
ErrCode The error code indicating success or failure.

Definition at line 2204 of file iMOAB.cpp.

2205 {
2206  ErrorCode rval;
2207  appData& data = context.appDatas[*pid];
2208 
2209  if( !data.local_verts.empty() ) // we should have no vertices in the app
2210  {
2211  return moab::MB_FAILURE;
2212  }
2213 
2214  int nverts = *coords_len / *dim;
2215 
2216  rval = context.MBI->create_vertices( coordinates, nverts, data.local_verts );MB_CHK_ERR( rval );
2217 
2218  rval = context.MBI->add_entities( data.file_set, data.local_verts );MB_CHK_ERR( rval );
2219 
2220  // also add the vertices to the all_verts range
2221  data.all_verts.merge( data.local_verts );
2222  return moab::MB_SUCCESS;
2223 }

References moab::Interface::add_entities(), appData::all_verts, GlobalContext::appDatas, context, moab::Interface::create_vertices(), dim, moab::Range::empty(), ErrorCode, appData::file_set, appData::local_verts, MB_CHK_ERR, MB_SUCCESS, GlobalContext::MBI, and moab::Range::merge().

◆ iMOAB_DefineTagStorage()

ErrCode iMOAB_DefineTagStorage ( iMOAB_AppID  pid,
const iMOAB_String  tag_storage_name,
int *  tag_type,
int *  components_per_entity,
int *  tag_index 
)

Define a MOAB Tag corresponding to the application depending on requested types.

Note
In MOAB, for most solution vectors, we only need to create a "Dense", "Double" Tag. A sparse tag can be created too. If the tag is already existing in the file, it will not be created. If it is a new tag, memory will be allocated when setting the values. Default values are 0 for for integer tags, 0.0 for double tags, 0 for entity handle tags.

Operations: Collective

Parameters
[in]pid(iMOAB_AppID) The unique pointer to the application ID.
[in]tag_storage_name(iMOAB_String) The tag name to store/retrieve the data in MOAB.
[in]tag_type(int*) The type of MOAB tag (Dense/Sparse, Double/Int/EntityHandle), enum MOAB_TAG_TYPE.
[in]components_per_entity(int*) The total size of vector dimension per entity for the tag (e.g., number of doubles per entity).
[out]tag_index(int*) The tag index which can be used as identifier in synchronize methods.
Returns
ErrCode The error code indicating success or failure.

Definition at line 1461 of file iMOAB.cpp.

1466 {
1467  // see if the tag is already existing, and if yes, check the type, length
1468  if( *tag_type < 0 || *tag_type > 5 )
1469  {
1470  return moab::MB_FAILURE;
1471  } // we have 6 types of tags supported so far
1472 
1473  DataType tagDataType;
1474  TagType tagType;
1475  void* defaultVal = NULL;
1476  int* defInt = new int[*components_per_entity];
1477  double* defDouble = new double[*components_per_entity];
1478  EntityHandle* defHandle = new EntityHandle[*components_per_entity];
1479 
1480  for( int i = 0; i < *components_per_entity; i++ )
1481  {
1482  defInt[i] = 0;
1483  defDouble[i] = -1e+10;
1484  defHandle[i] = (EntityHandle)0;
1485  }
1486 
1487  switch( *tag_type )
1488  {
1489  case 0:
1490  tagDataType = MB_TYPE_INTEGER;
1491  tagType = MB_TAG_DENSE;
1492  defaultVal = defInt;
1493  break;
1494 
1495  case 1:
1496  tagDataType = MB_TYPE_DOUBLE;
1497  tagType = MB_TAG_DENSE;
1498  defaultVal = defDouble;
1499  break;
1500 
1501  case 2:
1502  tagDataType = MB_TYPE_HANDLE;
1503  tagType = MB_TAG_DENSE;
1504  defaultVal = defHandle;
1505  break;
1506 
1507  case 3:
1508  tagDataType = MB_TYPE_INTEGER;
1509  tagType = MB_TAG_SPARSE;
1510  defaultVal = defInt;
1511  break;
1512 
1513  case 4:
1514  tagDataType = MB_TYPE_DOUBLE;
1515  tagType = MB_TAG_SPARSE;
1516  defaultVal = defDouble;
1517  break;
1518 
1519  case 5:
1520  tagDataType = MB_TYPE_HANDLE;
1521  tagType = MB_TAG_SPARSE;
1522  defaultVal = defHandle;
1523  break;
1524 
1525  default: {
1526  delete[] defInt;
1527  delete[] defDouble;
1528  delete[] defHandle;
1529  return moab::MB_FAILURE;
1530  } // error
1531  }
1532 
1533  Tag tagHandle;
1534  // split storage names if separated list
1535 
1536  std::string tag_name( tag_storage_name );
1537 
1538  // first separate the names of the tags
1539  // we assume that there are separators ":" between the tag names
1540  std::vector< std::string > tagNames;
1541  std::string separator( ":" );
1542  split_tag_names( tag_name, separator, tagNames );
1543 
1544  ErrorCode rval = moab::MB_SUCCESS; // assume success already :)
1545  appData& data = context.appDatas[*pid];
1546  int already_defined_tags = 0;
1547 
1548  for( size_t i = 0; i < tagNames.size(); i++ )
1549  {
1550  rval = context.MBI->tag_get_handle( tagNames[i].c_str(), *components_per_entity, tagDataType, tagHandle,
1551  tagType | MB_TAG_EXCL | MB_TAG_CREAT, defaultVal );
1552 
1553  if( MB_ALREADY_ALLOCATED == rval )
1554  {
1555  std::map< std::string, Tag >& mTags = data.tagMap;
1556  std::map< std::string, Tag >::iterator mit = mTags.find( tagNames[i].c_str() );
1557 
1558  if( mit == mTags.end() )
1559  {
1560  // add it to the map
1561  mTags[tagNames[i]] = tagHandle;
1562  // push it to the list of tags, too
1563  *tag_index = (int)data.tagList.size();
1564  data.tagList.push_back( tagHandle );
1565  }
1566  rval = MB_SUCCESS;
1567  already_defined_tags++;
1568  }
1569  else if( MB_SUCCESS == rval )
1570  {
1571  data.tagMap[tagNames[i]] = tagHandle;
1572  *tag_index = (int)data.tagList.size();
1573  data.tagList.push_back( tagHandle );
1574  }
1575  else
1576  {
1577  rval = moab::MB_FAILURE; // some tags were not created
1578  }
1579  }
1580  // we don't need default values anymore, avoid leaks
1581  int rankHere = 0;
1582 #ifdef MOAB_HAVE_MPI
1583  ParallelComm* pco = context.pcomms[*pid];
1584  rankHere = pco->rank();
1585 #endif
1586  if( !rankHere && already_defined_tags )
1587  std::cout << " application with ID: " << *pid << " global id: " << data.global_id << " name: " << data.name
1588  << " has " << already_defined_tags << " already defined tags out of " << tagNames.size()
1589  << " tags \n";
1590  delete[] defInt;
1591  delete[] defDouble;
1592  delete[] defHandle;
1593  return rval;
1594 }

References GlobalContext::appDatas, context, ErrorCode, appData::global_id, MB_ALREADY_ALLOCATED, MB_SUCCESS, MB_TAG_CREAT, MB_TAG_DENSE, MB_TAG_EXCL, MB_TAG_SPARSE, MB_TYPE_DOUBLE, MB_TYPE_HANDLE, MB_TYPE_INTEGER, GlobalContext::MBI, appData::name, moab::ParallelComm::rank(), split_tag_names(), moab::Interface::tag_get_handle(), appData::tagList, appData::tagMap, and TagType.

◆ iMOAB_DeregisterApplication()

ErrCode iMOAB_DeregisterApplication ( iMOAB_AppID  pid)

De-Register application: delete mesh (set) associated with the application ID.

Note
The associated communicator will be released, and all associated mesh entities and sets will be deleted from the mesh data structure. Associated tag storage data will be freed too.

Operations: Collective

Parameters
[in]pid(iMOAB_AppID) The unique pointer to the application ID
Returns
ErrCode The error code indicating success or failure.

Definition at line 329 of file iMOAB.cpp.

330 {
331  // the file set , parallel comm are all in vectors indexed by *pid
332  // assume we did not delete anything yet
333  // *pid will not be reused if we register another application
334  appData& data = context.appDatas[*pid];
335  int rankHere = 0;
336 #ifdef MOAB_HAVE_MPI
337  ParallelComm* pco = context.pcomms[*pid];
338  rankHere = pco->rank();
339 #endif
340  if( !rankHere )
341  std::cout << " application with ID: " << *pid << " global id: " << data.global_id << " name: " << data.name
342  << " is de-registered now \n";
343 
344  EntityHandle fileSet = data.file_set;
345  // get all entities part of the file set
346  Range fileents;
347  ErrorCode rval = context.MBI->get_entities_by_handle( fileSet, fileents, /*recursive */ true );MB_CHK_ERR( rval );
348 
349  fileents.insert( fileSet );
350 
351  rval = context.MBI->get_entities_by_type( fileSet, MBENTITYSET, fileents );MB_CHK_ERR( rval ); // append all mesh sets
352 
353 #ifdef MOAB_HAVE_TEMPESTREMAP
354  if( data.tempestData.remapper ) delete data.tempestData.remapper;
355  if( data.tempestData.weightMaps.size() ) data.tempestData.weightMaps.clear();
356 #endif
357 
358 #ifdef MOAB_HAVE_MPI
359 
360  // we could get the pco also with
361  // ParallelComm * pcomm = ParallelComm::get_pcomm(context.MBI, *pid);
362 
363  std::map< int, ParCommGraph* >& pargs = data.pgraph;
364 
365  // free the parallel comm graphs associated with this app
366  for( std::map< int, ParCommGraph* >::iterator mt = pargs.begin(); mt != pargs.end(); mt++ )
367  {
368  ParCommGraph* pgr = mt->second;
369  delete pgr;
370  pgr = NULL;
371  }
372  if( pco ) delete pco;
373 #endif
374 
375  // delete first all except vertices
376  Range vertices = fileents.subset_by_type( MBVERTEX );
377  Range noverts = subtract( fileents, vertices );
378 
379  rval = context.MBI->delete_entities( noverts );MB_CHK_ERR( rval );
380  // now retrieve connected elements that still exist (maybe in other sets, pids?)
381  Range adj_ents_left;
382  rval = context.MBI->get_adjacencies( vertices, 1, false, adj_ents_left, Interface::UNION );MB_CHK_ERR( rval );
383  rval = context.MBI->get_adjacencies( vertices, 2, false, adj_ents_left, Interface::UNION );MB_CHK_ERR( rval );
384  rval = context.MBI->get_adjacencies( vertices, 3, false, adj_ents_left, Interface::UNION );MB_CHK_ERR( rval );
385 
386  if( !adj_ents_left.empty() )
387  {
388  Range conn_verts;
389  rval = context.MBI->get_connectivity( adj_ents_left, conn_verts );MB_CHK_ERR( rval );
390  vertices = subtract( vertices, conn_verts );
391  }
392 
393  rval = context.MBI->delete_entities( vertices );MB_CHK_ERR( rval );
394 
395  std::map< std::string, int >::iterator mit;
396 
397  for( mit = context.appIdMap.begin(); mit != context.appIdMap.end(); mit++ )
398  {
399  int pidx = mit->second;
400 
401  if( *pid == pidx )
402  {
403  break;
404  }
405  }
406 
407  context.appIdMap.erase( mit );
408  std::map< int, int >::iterator mit1;
409 
410  for( mit1 = context.appIdCompMap.begin(); mit1 != context.appIdCompMap.end(); mit1++ )
411  {
412  int pidx = mit1->second;
413 
414  if( *pid == pidx )
415  {
416  break;
417  }
418  }
419 
420  context.appIdCompMap.erase( mit1 );
421 
422  context.unused_pid--; // we have to go backwards always TODO
423  context.appDatas.pop_back();
424 #ifdef MOAB_HAVE_MPI
425  context.pcomms.pop_back();
426 #endif
427  return moab::MB_SUCCESS;
428 }

References GlobalContext::appDatas, GlobalContext::appIdCompMap, GlobalContext::appIdMap, context, moab::Interface::delete_entities(), moab::Range::empty(), ErrorCode, appData::file_set, moab::Interface::get_adjacencies(), moab::Interface::get_connectivity(), moab::Interface::get_entities_by_handle(), moab::Interface::get_entities_by_type(), appData::global_id, moab::Range::insert(), MB_CHK_ERR, MB_SUCCESS, MBENTITYSET, GlobalContext::MBI, MBVERTEX, appData::name, moab::ParallelComm::rank(), moab::Range::subset_by_type(), moab::subtract(), moab::Interface::UNION, and GlobalContext::unused_pid.

Referenced by iMOAB_DeregisterApplicationFortran().

◆ iMOAB_DeregisterApplicationFortran()

ErrCode iMOAB_DeregisterApplicationFortran ( iMOAB_AppID  pid)

De-Register the Fortran based application: delete mesh (set) associated with the application ID.

Note
The associated communicator will be released, and all associated mesh entities and sets will be deleted from the mesh data structure. Associated tag storage data will be freed too.

Operations: Collective

Parameters
[in]pid(iMOAB_AppID) The unique pointer to the application ID
Returns
ErrCode The error code indicating success or failure.

Definition at line 430 of file iMOAB.cpp.

431 {
432  // release any Fortran specific allocations here before we pass it on
433  context.appDatas[*pid].is_fortran = false;
434 
435  // release all datastructure allocations
436  return iMOAB_DeregisterApplication( pid );
437 }

References GlobalContext::appDatas, context, and iMOAB_DeregisterApplication().

◆ iMOAB_Finalize()

ErrCode iMOAB_Finalize ( )

Finalize the iMOAB interface implementation.

It will delete the internally reference counted MOAB instance, if the reference count reaches 0.

Operations: Collective

Returns
ErrCode The error code indicating success or failure.

Definition at line 186 of file iMOAB.cpp.

187 {
189 
190  if( 0 == context.refCountMB )
191  {
192  delete context.MBI;
193  }
194 
195  return MB_SUCCESS;
196 }

References context, MB_SUCCESS, GlobalContext::MBI, and GlobalContext::refCountMB.

◆ iMOAB_GetBlockElementConnectivities()

ErrCode iMOAB_GetBlockElementConnectivities ( iMOAB_AppID  pid,
iMOAB_GlobalID global_block_ID,
int *  connectivity_length,
int *  element_connectivity 
)

Get the connectivities for elements contained within a certain block.

Note
input is the block ID. Should we change to visible block local index?

Operations: Not Collective

Parameters
[in]pid(iMOAB_AppID) The unique pointer to the application ID.
[in]global_block_ID(iMOAB_GlobalID*) The global block ID of the set being queried.
[in]connectivity_length(int*) The allocated size of array. (typical size := vertices_per_element*num_elements_in_block)
[out]element_connectivity(int*) The connectivity array to store element ordering in MOAB canonical numbering scheme (array allocated by user); array contains vertex indices in the local numbering order for vertices elements are in the same order as provided by GetElementOwnership and GetElementID.
Returns
ErrCode The error code indicating success or failure.

Definition at line 1153 of file iMOAB.cpp.

1157 {
1158  assert( global_block_ID ); // ensure global block ID argument is not null
1159  assert( connectivity_length ); // ensure connectivity length argument is not null
1160 
1161  appData& data = context.appDatas[*pid];
1162  std::map< int, int >& matMap = data.matIndex;
1163  std::map< int, int >::iterator it = matMap.find( *global_block_ID );
1164 
1165  if( it == matMap.end() )
1166  {
1167  return moab::MB_FAILURE;
1168  } // error in finding block with id
1169 
1170  int blockIndex = matMap[*global_block_ID];
1171  EntityHandle matMeshSet = data.mat_sets[blockIndex];
1172  std::vector< EntityHandle > elems;
1173 
1174  ErrorCode rval = context.MBI->get_entities_by_handle( matMeshSet, elems );MB_CHK_ERR( rval );
1175 
1176  if( elems.empty() )
1177  {
1178  return moab::MB_FAILURE;
1179  }
1180 
1181  std::vector< EntityHandle > vconnect;
1182  rval = context.MBI->get_connectivity( &elems[0], elems.size(), vconnect );MB_CHK_ERR( rval );
1183 
1184  if( *connectivity_length != (int)vconnect.size() )
1185  {
1186  return moab::MB_FAILURE;
1187  } // mismatched sizes
1188 
1189  for( int i = 0; i < *connectivity_length; i++ )
1190  {
1191  int inx = data.all_verts.index( vconnect[i] );
1192 
1193  if( -1 == inx )
1194  {
1195  return moab::MB_FAILURE;
1196  } // error, vertex not in local range
1197 
1198  element_connectivity[i] = inx;
1199  }
1200 
1201  return moab::MB_SUCCESS;
1202 }

References appData::all_verts, GlobalContext::appDatas, context, ErrorCode, moab::Interface::get_connectivity(), moab::Interface::get_entities_by_handle(), moab::Range::index(), appData::mat_sets, appData::matIndex, MB_CHK_ERR, MB_SUCCESS, and GlobalContext::MBI.

◆ iMOAB_GetBlockID()

ErrCode iMOAB_GetBlockID ( iMOAB_AppID  pid,
int *  block_length,
iMOAB_GlobalID global_block_IDs 
)

Get the global block IDs for all locally visible (owned and shared/ghosted) blocks.

Note
Block IDs are corresponding to MATERIAL_SET tags for material sets. Usually the block ID is exported from Cubit as a unique integer value. First blocks are local, and next blocks are fully ghosted. First blocks have at least one owned cell/element, ghost blocks have only ghost cells. Internally, a block corresponds to a mesh set with a MATERIAL_SET tag value equal to the block ID.

Operations: Not Collective

Parameters
[in]pid(iMOAB_AppID) The unique pointer to the application ID.
[in]block_length(int*) The allocated size of array (typical size := num_visible_blocks).
[out]global_block_IDs(iMOAB_GlobalID*) The global IDs for all locally visible blocks (array allocated by user).
Returns
ErrCode The error code indicating success or failure.

Definition at line 1030 of file iMOAB.cpp.

1031 {
1032  // local id blocks? they are counted from 0 to number of visible blocks ...
1033  // will actually return material set tag value for global
1034  Range& matSets = context.appDatas[*pid].mat_sets;
1035 
1036  if( *block_length != (int)matSets.size() )
1037  {
1038  return moab::MB_FAILURE;
1039  }
1040 
1041  // return material set tag gtags[0 is material set tag
1042  ErrorCode rval = context.MBI->tag_get_data( context.material_tag, matSets, global_block_IDs );MB_CHK_ERR( rval );
1043 
1044  // populate map with index
1045  std::map< int, int >& matIdx = context.appDatas[*pid].matIndex;
1046  for( unsigned i = 0; i < matSets.size(); i++ )
1047  {
1048  matIdx[global_block_IDs[i]] = i;
1049  }
1050 
1051  return moab::MB_SUCCESS;
1052 }

References GlobalContext::appDatas, context, ErrorCode, GlobalContext::material_tag, MB_CHK_ERR, MB_SUCCESS, GlobalContext::MBI, moab::Range::size(), and moab::Interface::tag_get_data().

◆ iMOAB_GetBlockInfo()

ErrCode iMOAB_GetBlockInfo ( iMOAB_AppID  pid,
iMOAB_GlobalID global_block_ID,
int *  vertices_per_element,
int *  num_elements_in_block 
)

Get the global block information and number of visible elements of belonging to a block (MATERIAL SET).

Note
A block has to be homogeneous, it can contain elements of a single type.

Operations: Not Collective

Parameters
[in]pid(iMOAB_AppID) The unique pointer to the application ID.
[in]global_block_ID(iMOAB_GlobalID) The global block ID of the set to be queried.
[out]vertices_per_element(int*) The number of vertices per element.
[out]num_elements_in_block(int*) The number of elements in block.
Returns
ErrCode The error code indicating success or failure.

Definition at line 1054 of file iMOAB.cpp.

1058 {
1059  assert( global_block_ID );
1060 
1061  std::map< int, int >& matMap = context.appDatas[*pid].matIndex;
1062  std::map< int, int >::iterator it = matMap.find( *global_block_ID );
1063 
1064  if( it == matMap.end() )
1065  {
1066  return moab::MB_FAILURE;
1067  } // error in finding block with id
1068 
1069  int blockIndex = matMap[*global_block_ID];
1070  EntityHandle matMeshSet = context.appDatas[*pid].mat_sets[blockIndex];
1071  Range blo_elems;
1072  ErrorCode rval = context.MBI->get_entities_by_handle( matMeshSet, blo_elems );
1073 
1074  if( MB_SUCCESS != rval || blo_elems.empty() )
1075  {
1076  return moab::MB_FAILURE;
1077  }
1078 
1079  EntityType type = context.MBI->type_from_handle( blo_elems[0] );
1080 
1081  if( !blo_elems.all_of_type( type ) )
1082  {
1083  return moab::MB_FAILURE;
1084  } // not all of same type
1085 
1086  const EntityHandle* conn = NULL;
1087  int num_verts = 0;
1088  rval = context.MBI->get_connectivity( blo_elems[0], conn, num_verts );MB_CHK_ERR( rval );
1089 
1090  *vertices_per_element = num_verts;
1091  *num_elements_in_block = (int)blo_elems.size();
1092 
1093  return moab::MB_SUCCESS;
1094 }

References moab::Range::all_of_type(), GlobalContext::appDatas, context, moab::Range::empty(), ErrorCode, moab::Interface::get_connectivity(), moab::Interface::get_entities_by_handle(), MB_CHK_ERR, MB_SUCCESS, GlobalContext::MBI, moab::Range::size(), and moab::Interface::type_from_handle().

◆ iMOAB_GetDoubleTagStorage()

ErrCode iMOAB_GetDoubleTagStorage ( iMOAB_AppID  pid,
const iMOAB_String  tag_storage_name,
int *  num_tag_storage_length,
int *  entity_type,
double *  tag_storage_data 
)

Retrieve the specified values in a MOAB double Tag.

Note
Operations: Collective
Parameters
[in]pid(iMOAB_AppID) The unique pointer to the application ID.
[in]tag_storage_name(iMOAB_String) The tag names to retrieve the data in MOAB.
[in]num_tag_storage_length(int) The size of total tag storage data (e.g., num_visible_vertices*components_per_entity*num_tags or num_visible_elements*components_per_entity*num_tags)
[in]entity_type(int*) Type=0 for vertices, and Type=1 for primary elements.
[out]tag_storage_data(double*) The array data of type double to replace the internal tag memory; The data is assumed to be contiguous over the local set of visible entities (either vertices or elements). unrolled by tags
Returns
ErrCode The error code indicating success or failure.

Definition at line 2021 of file iMOAB.cpp.

2026 {
2027  ErrorCode rval;
2028  // exactly the same code, except tag type check
2029  std::string tag_names( tag_storage_names );
2030  // exactly the same code as for int tag :) maybe should check the type of tag too
2031  std::vector< std::string > tagNames;
2032  std::vector< Tag > tagHandles;
2033  std::string separator( ":" );
2034  split_tag_names( tag_names, separator, tagNames );
2035 
2036  appData& data = context.appDatas[*pid];
2037 
2038  // set it on a subset of entities, based on type and length
2039  Range* ents_to_get = NULL;
2040 
2041  if( *ent_type == 0 ) // vertices
2042  {
2043  ents_to_get = &data.all_verts;
2044  }
2045  else if( *ent_type == 1 )
2046  {
2047  ents_to_get = &data.primary_elems;
2048  }
2049  int nents_to_get = (int)ents_to_get->size();
2050  int position = 0;
2051  for( size_t i = 0; i < tagNames.size(); i++ )
2052  {
2053  if( data.tagMap.find( tagNames[i] ) == data.tagMap.end() )
2054  {
2055  return moab::MB_FAILURE;
2056  } // tag not defined
2057 
2058  Tag tag = data.tagMap[tagNames[i]];
2059 
2060  int tagLength = 0;
2061  rval = context.MBI->tag_get_length( tag, tagLength );MB_CHK_ERR( rval );
2062 
2063  DataType dtype;
2064  rval = context.MBI->tag_get_data_type( tag, dtype );MB_CHK_ERR( rval );
2065 
2066  if( dtype != MB_TYPE_DOUBLE )
2067  {
2068  return moab::MB_FAILURE;
2069  }
2070 
2071  if( position + nents_to_get * tagLength > *num_tag_storage_length )
2072  return moab::MB_FAILURE; // too many entity values to get
2073 
2074  rval = context.MBI->tag_get_data( tag, *ents_to_get, &tag_storage_data[position] );MB_CHK_ERR( rval );
2075  position = position + nents_to_get * tagLength;
2076  }
2077 
2078  return moab::MB_SUCCESS; // no error
2079 }

References appData::all_verts, GlobalContext::appDatas, context, ErrorCode, MB_CHK_ERR, MB_SUCCESS, MB_TYPE_DOUBLE, GlobalContext::MBI, appData::primary_elems, moab::Range::size(), split_tag_names(), moab::Interface::tag_get_data(), moab::Interface::tag_get_data_type(), moab::Interface::tag_get_length(), and appData::tagMap.

◆ iMOAB_GetElementConnectivity()

ErrCode iMOAB_GetElementConnectivity ( iMOAB_AppID  pid,
iMOAB_LocalID elem_index,
int *  connectivity_length,
int *  element_connectivity 
)

Get the connectivity for one element only.

Note
This is a convenience method and in general should not be used unless necessary. You should use colaesced calls for multiple elements for efficiency.

Operations: Not Collective

Parameters
[in]pid(iMOAB_AppID) The unique pointer to the application ID.
[in]elem_index(iMOAB_LocalID *) Local element index.
[in,out]connectivity_length(int *) On input, maximum length of connectivity. On output, actual length.
[out]element_connectivity(int*) The connectivity array to store connectivity in MOAB canonical numbering scheme. Array contains vertex indices in the local numbering order for vertices.
Returns
ErrCode The error code indicating success or failure.

Definition at line 1204 of file iMOAB.cpp.

1208 {
1209  assert( elem_index ); // ensure element index argument is not null
1210  assert( connectivity_length ); // ensure connectivity length argument is not null
1211 
1212  appData& data = context.appDatas[*pid];
1213  assert( ( *elem_index >= 0 ) && ( *elem_index < (int)data.primary_elems.size() ) );
1214 
1215  int num_nodes;
1216  const EntityHandle* conn;
1217 
1218  EntityHandle eh = data.primary_elems[*elem_index];
1219 
1220  ErrorCode rval = context.MBI->get_connectivity( eh, conn, num_nodes );MB_CHK_ERR( rval );
1221 
1222  if( *connectivity_length < num_nodes )
1223  {
1224  return moab::MB_FAILURE;
1225  } // wrong number of vertices
1226 
1227  for( int i = 0; i < num_nodes; i++ )
1228  {
1229  int index = data.all_verts.index( conn[i] );
1230 
1231  if( -1 == index )
1232  {
1233  return moab::MB_FAILURE;
1234  }
1235 
1236  element_connectivity[i] = index;
1237  }
1238 
1239  *connectivity_length = num_nodes;
1240 
1241  return moab::MB_SUCCESS;
1242 }

References appData::all_verts, GlobalContext::appDatas, context, ErrorCode, moab::Interface::get_connectivity(), moab::Range::index(), MB_CHK_ERR, MB_SUCCESS, GlobalContext::MBI, appData::primary_elems, and moab::Range::size().

◆ iMOAB_GetElementID()

ErrCode iMOAB_GetElementID ( iMOAB_AppID  pid,
iMOAB_GlobalID global_block_ID,
int *  num_elements_in_block,
iMOAB_GlobalID global_element_ID,
iMOAB_LocalID local_element_ID 
)

Get the global IDs for all locally visible elements belonging to a particular block.

Note
The method will return also the local index of each element, in the local range that contains all visible elements.

Operations: Not Collective

Parameters
[in]pid(iMOAB_AppID) The unique pointer to the application ID.
[in]global_block_ID(iMOAB_GlobalID*) The global block ID of the set being queried.
[in]num_elements_in_block(int*) The allocated size of global element ID array, same as num_elements_in_block returned from GetBlockInfo().
[out]global_element_ID(iMOAB_GlobalID*) The global IDs for all locally visible elements (array allocated by user).
[out]local_element_ID(iMOAB_LocalID*) (Optional) The local IDs for all locally visible elements (index in the range of all primary elements in the rank)
Returns
ErrCode The error code indicating success or failure.

Definition at line 1294 of file iMOAB.cpp.

1299 {
1300  assert( global_block_ID ); // ensure global block ID argument is not null
1301  assert( num_elements_in_block ); // ensure number of elements in block argument is not null
1302 
1303  appData& data = context.appDatas[*pid];
1304  std::map< int, int >& matMap = data.matIndex;
1305 
1306  std::map< int, int >::iterator it = matMap.find( *global_block_ID );
1307 
1308  if( it == matMap.end() )
1309  {
1310  return moab::MB_FAILURE;
1311  } // error in finding block with id
1312 
1313  int blockIndex = matMap[*global_block_ID];
1314  EntityHandle matMeshSet = data.mat_sets[blockIndex];
1315  Range elems;
1316  ErrorCode rval = context.MBI->get_entities_by_handle( matMeshSet, elems );MB_CHK_ERR( rval );
1317 
1318  if( elems.empty() )
1319  {
1320  return moab::MB_FAILURE;
1321  }
1322 
1323  if( *num_elements_in_block != (int)elems.size() )
1324  {
1325  return moab::MB_FAILURE;
1326  } // bad memory allocation
1327 
1328  rval = context.MBI->tag_get_data( context.globalID_tag, elems, global_element_ID );MB_CHK_ERR( rval );
1329 
1330  // check that elems are among primary_elems in data
1331  for( int i = 0; i < *num_elements_in_block; i++ )
1332  {
1333  local_element_ID[i] = data.primary_elems.index( elems[i] );
1334 
1335  if( -1 == local_element_ID[i] )
1336  {
1337  return moab::MB_FAILURE;
1338  } // error, not in local primary elements
1339  }
1340 
1341  return moab::MB_SUCCESS;
1342 }

References GlobalContext::appDatas, context, moab::Range::empty(), ErrorCode, moab::Interface::get_entities_by_handle(), GlobalContext::globalID_tag, moab::Range::index(), appData::mat_sets, appData::matIndex, MB_CHK_ERR, MB_SUCCESS, GlobalContext::MBI, appData::primary_elems, moab::Range::size(), and moab::Interface::tag_get_data().

◆ iMOAB_GetElementOwnership()

ErrCode iMOAB_GetElementOwnership ( iMOAB_AppID  pid,
iMOAB_GlobalID global_block_ID,
int *  num_elements_in_block,
int *  element_ownership 
)

Get the element ownership within a certain block i.e., processor ID of the element owner.

Note
: Should we use local block index for input, instead of the global block ID ?

Operations: Not Collective

Parameters
[in]pid(iMOAB_AppID) The unique pointer to the application ID.
[in]global_block_ID(iMOAB_GlobalID) The global block ID of the set being queried.
[in]num_elements_in_block(int*) The allocated size of ownership array, same as num_elements_in_block returned from GetBlockInfo().
[out]element_ownership(int*) The ownership array to store processor ID for all elements (array allocated by user).
Returns
ErrCode The error code indicating success or failure.

Definition at line 1244 of file iMOAB.cpp.

1248 {
1249  assert( global_block_ID ); // ensure global block ID argument is not null
1250  assert( num_elements_in_block ); // ensure number of elements in block argument is not null
1251 
1252  std::map< int, int >& matMap = context.appDatas[*pid].matIndex;
1253 
1254  std::map< int, int >::iterator it = matMap.find( *global_block_ID );
1255 
1256  if( it == matMap.end() )
1257  {
1258  return moab::MB_FAILURE;
1259  } // error in finding block with id
1260 
1261  int blockIndex = matMap[*global_block_ID];
1262  EntityHandle matMeshSet = context.appDatas[*pid].mat_sets[blockIndex];
1263  Range elems;
1264 
1265  ErrorCode rval = context.MBI->get_entities_by_handle( matMeshSet, elems );MB_CHK_ERR( rval );
1266 
1267  if( elems.empty() )
1268  {
1269  return moab::MB_FAILURE;
1270  }
1271 
1272  if( *num_elements_in_block != (int)elems.size() )
1273  {
1274  return moab::MB_FAILURE;
1275  } // bad memory allocation
1276 
1277  int i = 0;
1278 #ifdef MOAB_HAVE_MPI
1279  ParallelComm* pco = context.pcomms[*pid];
1280 #endif
1281 
1282  for( Range::iterator vit = elems.begin(); vit != elems.end(); vit++, i++ )
1283  {
1284 #ifdef MOAB_HAVE_MPI
1285  rval = pco->get_owner( *vit, element_ownership[i] );MB_CHK_ERR( rval );
1286 #else
1287  element_ownership[i] = 0; /* owned by 0 */
1288 #endif
1289  }
1290 
1291  return moab::MB_SUCCESS;
1292 }

References GlobalContext::appDatas, moab::Range::begin(), context, moab::Range::empty(), moab::Range::end(), ErrorCode, moab::Interface::get_entities_by_handle(), moab::ParallelComm::get_owner(), MB_CHK_ERR, MB_SUCCESS, GlobalContext::MBI, and moab::Range::size().

◆ iMOAB_GetGlobalInfo()

ErrCode iMOAB_GetGlobalInfo ( iMOAB_AppID  pid,
int *  num_global_verts,
int *  num_global_elems 
)

Get global information about number of vertices and number of elements.

Parameters
[in]pid(iMOAB_AppID) The unique pointer to the application ID.
[in]num_global_verts(int*) The number of total vertices in the mesh.
[in]num_global_elems(MPI_Comm) The number of total elements in the mesh.
Returns
ErrCode The error code indicating success or failure.

Definition at line 2298 of file iMOAB.cpp.

2299 {
2300  appData& data = context.appDatas[*pid];
2301  if( NULL != num_global_verts )
2302  {
2303  *num_global_verts = data.num_global_vertices;
2304  }
2305  if( NULL != num_global_elems )
2306  {
2307  *num_global_elems = data.num_global_elements;
2308  }
2309 
2310  return moab::MB_SUCCESS;
2311 }

References GlobalContext::appDatas, context, MB_SUCCESS, appData::num_global_elements, and appData::num_global_vertices.

◆ iMOAB_GetIntTagStorage()

ErrCode iMOAB_GetIntTagStorage ( iMOAB_AppID  pid,
const iMOAB_String  tag_storage_name,
int *  num_tag_storage_length,
int *  entity_type,
int *  tag_storage_data 
)

Get the specified values in a MOAB integer Tag.

Note
Operations: Collective
Parameters
[in]pid(iMOAB_AppID) The unique pointer to the application ID.
[in]tag_storage_name(iMOAB_String) The tag name to store/retreive the data in MOAB.
[in]num_tag_storage_length(int*) The size of tag storage data (e.g., num_visible_vertices*components_per_entity or num_visible_elements*components_per_entity).
[in]entity_type(int*) Type=0 for vertices, and Type=1 for primary elements.
[out]tag_storage_data(int*) The array data of type int to be copied from the internal tag memory; The data is assumed to be contiguous over the local set of visible entities (either vertices or elements).
Returns
ErrCode The error code indicating success or failure.

Definition at line 1653 of file iMOAB.cpp.

1658 {
1659  ErrorCode rval;
1660  std::string tag_name( tag_storage_name );
1661 
1662  appData& data = context.appDatas[*pid];
1663 
1664  if( data.tagMap.find( tag_name ) == data.tagMap.end() )
1665  {
1666  return moab::MB_FAILURE;
1667  } // tag not defined
1668 
1669  Tag tag = data.tagMap[tag_name];
1670 
1671  int tagLength = 0;
1672  rval = context.MBI->tag_get_length( tag, tagLength );MB_CHK_ERR( rval );
1673 
1674  DataType dtype;
1675  rval = context.MBI->tag_get_data_type( tag, dtype );MB_CHK_ERR( rval );
1676 
1677  if( dtype != MB_TYPE_INTEGER )
1678  {
1679  return moab::MB_FAILURE;
1680  }
1681 
1682  // set it on a subset of entities, based on type and length
1683  Range* ents_to_get;
1684 
1685  if( *ent_type == 0 ) // vertices
1686  {
1687  ents_to_get = &data.all_verts;
1688  }
1689  else // if (*ent_type == 1)
1690  {
1691  ents_to_get = &data.primary_elems;
1692  }
1693 
1694  int nents_to_get = *num_tag_storage_length / tagLength;
1695 
1696  if( nents_to_get > (int)ents_to_get->size() || nents_to_get < 1 )
1697  {
1698  return moab::MB_FAILURE;
1699  } // to many entities to get, or too little
1700 
1701  // restrict the range; everything is contiguous; or not?
1702  // Range contig_range ( * ( ents_to_get->begin() ), * ( ents_to_get->begin() + nents_to_get - 1
1703  // ) );
1704 
1705  rval = context.MBI->tag_get_data( tag, *ents_to_get, tag_storage_data );MB_CHK_ERR( rval );
1706 
1707  return moab::MB_SUCCESS; // no error
1708 }

References appData::all_verts, GlobalContext::appDatas, context, ErrorCode, MB_CHK_ERR, MB_SUCCESS, MB_TYPE_INTEGER, GlobalContext::MBI, appData::primary_elems, moab::Range::size(), moab::Interface::tag_get_data(), moab::Interface::tag_get_data_type(), moab::Interface::tag_get_length(), and appData::tagMap.

◆ iMOAB_GetMeshInfo()

ErrCode iMOAB_GetMeshInfo ( iMOAB_AppID  pid,
int *  num_visible_vertices,
int *  num_visible_elements,
int *  num_visible_blocks,
int *  num_visible_surfaceBC,
int *  num_visible_vertexBC 
)

Retrieve all the important mesh information related to vertices, elements, vertex and surface boundary conditions.

Note
Number of visible vertices and cells include ghost entities. All arrays returned have size 3. Local entities are first, then ghost entities are next. Shared vertices can be owned in MOAB sense by different tasks. Ghost vertices and cells are always owned by other tasks.

Operations: Not Collective

Parameters
[in]pid(iMOAB_AppID) The unique pointer to the application ID.
[out]num_visible_vertices(int*) The number of vertices in the current partition/process arranged. as: owned/shared only, ghosted, total_visible (array allocated by user, size := 3).
[out]num_visible_elements(int*) The number of elements in current partition/process arranged as: owned only, ghosted/shared, total_visible (array allocated by user, size := 3).
[out]num_visible_blocks(int*) The number of material sets in local mesh in current partition/process arranged as: owned only, ghosted/shared, total_visible (array allocated by user, size := 3).
[out]num_visible_surfaceBC(int*) The number of mesh surfaces that have a NEUMANN_SET BC defined in local mesh in current partition/process arranged as: owned only, ghosted/shared, total_visible (array allocated by user, size := 3).
[out]num_visible_vertexBC(int*) The number of vertices that have a DIRICHLET_SET BC defined in local mesh in current partition/process arranged as: owned only, ghosted/shared, total_visible (array allocated by user, size := 3).
Returns
ErrCode The error code indicating success or failure.

Definition at line 872 of file iMOAB.cpp.

878 {
879  ErrorCode rval;
880  appData& data = context.appDatas[*pid];
881  EntityHandle fileSet = data.file_set;
882 
883  // this will include ghost elements
884  // first clear all data ranges; this can be called after ghosting
885  if( num_visible_elements )
886  {
887  num_visible_elements[2] = static_cast< int >( data.primary_elems.size() );
888  // separate ghost and local/owned primary elements
889  num_visible_elements[0] = static_cast< int >( data.owned_elems.size() );
890  num_visible_elements[1] = static_cast< int >( data.ghost_elems.size() );
891  }
892  if( num_visible_vertices )
893  {
894  num_visible_vertices[2] = static_cast< int >( data.all_verts.size() );
895  num_visible_vertices[1] = static_cast< int >( data.ghost_vertices.size() );
896  // local are those that are not ghosts; they may include shared, not owned vertices
897  num_visible_vertices[0] = num_visible_vertices[2] - num_visible_vertices[1];
898  }
899 
900  if( num_visible_blocks )
901  {
903  data.mat_sets, Interface::UNION );MB_CHK_ERR( rval );
904 
905  num_visible_blocks[2] = data.mat_sets.size();
906  num_visible_blocks[0] = num_visible_blocks[2];
907  num_visible_blocks[1] = 0;
908  }
909 
910  if( num_visible_surfaceBC )
911  {
913  data.neu_sets, Interface::UNION );MB_CHK_ERR( rval );
914 
915  num_visible_surfaceBC[2] = 0;
916  // count how many faces are in each neu set, and how many regions are
917  // adjacent to them;
918  int numNeuSets = (int)data.neu_sets.size();
919 
920  for( int i = 0; i < numNeuSets; i++ )
921  {
922  Range subents;
923  EntityHandle nset = data.neu_sets[i];
924  rval = context.MBI->get_entities_by_dimension( nset, data.dimension - 1, subents );MB_CHK_ERR( rval );
925 
926  for( Range::iterator it = subents.begin(); it != subents.end(); ++it )
927  {
928  EntityHandle subent = *it;
929  Range adjPrimaryEnts;
930  rval = context.MBI->get_adjacencies( &subent, 1, data.dimension, false, adjPrimaryEnts );MB_CHK_ERR( rval );
931 
932  num_visible_surfaceBC[2] += (int)adjPrimaryEnts.size();
933  }
934  }
935 
936  num_visible_surfaceBC[0] = num_visible_surfaceBC[2];
937  num_visible_surfaceBC[1] = 0;
938  }
939 
940  if( num_visible_vertexBC )
941  {
943  data.diri_sets, Interface::UNION );MB_CHK_ERR( rval );
944 
945  num_visible_vertexBC[2] = 0;
946  int numDiriSets = (int)data.diri_sets.size();
947 
948  for( int i = 0; i < numDiriSets; i++ )
949  {
950  Range verts;
951  EntityHandle diset = data.diri_sets[i];
952  rval = context.MBI->get_entities_by_dimension( diset, 0, verts );MB_CHK_ERR( rval );
953 
954  num_visible_vertexBC[2] += (int)verts.size();
955  }
956 
957  num_visible_vertexBC[0] = num_visible_vertexBC[2];
958  num_visible_vertexBC[1] = 0;
959  }
960 
961  return moab::MB_SUCCESS;
962 }

References appData::all_verts, GlobalContext::appDatas, moab::Range::begin(), context, appData::dimension, appData::diri_sets, GlobalContext::dirichlet_tag, moab::Range::end(), ErrorCode, appData::file_set, moab::Interface::get_adjacencies(), moab::Interface::get_entities_by_dimension(), moab::Interface::get_entities_by_type_and_tag(), appData::ghost_elems, appData::ghost_vertices, appData::mat_sets, GlobalContext::material_tag, MB_CHK_ERR, MB_SUCCESS, MBENTITYSET, GlobalContext::MBI, appData::neu_sets, GlobalContext::neumann_tag, appData::owned_elems, appData::primary_elems, moab::Range::size(), and moab::Interface::UNION.

◆ iMOAB_GetNeighborElements()

ErrCode iMOAB_GetNeighborElements ( iMOAB_AppID  pid,
iMOAB_LocalID local_index,
int *  num_adjacent_elements,
iMOAB_LocalID adjacent_element_IDs 
)

Retrieve the adjacencies for the element entities.

Note
Operations: Not Collective
Parameters
[in]pid(iMOAB_AppID) The unique pointer to the application ID.
[in]local_index(iMOAB_LocalID*) The local element ID for which adjacency information is needed.
[out]num_adjacent_elements(int*) The total number of adjacent elements.
[out]adjacent_element_IDs(iMOAB_LocalID*) The local element IDs of all adjacent elements to the current one (typically, num_total_sides for internal elements or num_total_sides-num_sides_on_boundary for boundary elements).
Returns
ErrCode The error code indicating success or failure.

Definition at line 2166 of file iMOAB.cpp.

2170 {
2171  ErrorCode rval;
2172 
2173  // one neighbor for each subentity of dimension-1
2174  MeshTopoUtil mtu( context.MBI );
2175  appData& data = context.appDatas[*pid];
2176  EntityHandle eh = data.primary_elems[*local_index];
2177  Range adjs;
2178  rval = mtu.get_bridge_adjacencies( eh, data.dimension - 1, data.dimension, adjs );MB_CHK_ERR( rval );
2179 
2180  if( *num_adjacent_elements < (int)adjs.size() )
2181  {
2182  return moab::MB_FAILURE;
2183  } // not dimensioned correctly
2184 
2185  *num_adjacent_elements = (int)adjs.size();
2186 
2187  for( int i = 0; i < *num_adjacent_elements; i++ )
2188  {
2189  adjacent_element_IDs[i] = data.primary_elems.index( adjs[i] );
2190  }
2191 
2192  return moab::MB_SUCCESS;
2193 }

References GlobalContext::appDatas, context, appData::dimension, ErrorCode, moab::MeshTopoUtil::get_bridge_adjacencies(), moab::Range::index(), MB_CHK_ERR, MB_SUCCESS, GlobalContext::MBI, appData::primary_elems, and moab::Range::size().

◆ iMOAB_GetPointerToSurfaceBC()

ErrCode iMOAB_GetPointerToSurfaceBC ( iMOAB_AppID  pid,
int *  surface_BC_length,
iMOAB_LocalID local_element_ID,
int *  reference_surface_ID,
int *  boundary_condition_value 
)

Get the surface boundary condition information.

Note
Operations: Not Collective
Parameters
[in]pid(iMOAB_AppID) The unique pointer to the application ID.
[in]surface_BC_length(int*) The allocated size of surface boundary condition array, same as num_visible_surfaceBC returned by GetMeshInfo().
[out]local_element_ID(iMOAB_LocalID*) The local element IDs that contains the side with the surface BC.
[out]reference_surface_ID(int*) The surface number with the BC in the canonical reference element (e.g., 1 to 6 for HEX, 1-4 for TET).
[out]boundary_condition_value(int*) The boundary condition type as obtained from the mesh description (value of the NeumannSet defined on the element).
Returns
ErrCode The error code indicating success or failure.

Definition at line 1344 of file iMOAB.cpp.

1349 {
1350  // we have to fill bc data for neumann sets;/
1351  ErrorCode rval;
1352 
1353  // it was counted above, in GetMeshInfo
1354  appData& data = context.appDatas[*pid];
1355  int numNeuSets = (int)data.neu_sets.size();
1356 
1357  int index = 0; // index [0, surface_BC_length) for the arrays returned
1358 
1359  for( int i = 0; i < numNeuSets; i++ )
1360  {
1361  Range subents;
1362  EntityHandle nset = data.neu_sets[i];
1363  rval = context.MBI->get_entities_by_dimension( nset, data.dimension - 1, subents );MB_CHK_ERR( rval );
1364 
1365  int neuVal;
1366  rval = context.MBI->tag_get_data( context.neumann_tag, &nset, 1, &neuVal );MB_CHK_ERR( rval );
1367 
1368  for( Range::iterator it = subents.begin(); it != subents.end(); ++it )
1369  {
1370  EntityHandle subent = *it;
1371  Range adjPrimaryEnts;
1372  rval = context.MBI->get_adjacencies( &subent, 1, data.dimension, false, adjPrimaryEnts );MB_CHK_ERR( rval );
1373 
1374  // get global id of the primary ents, and side number of the quad/subentity
1375  // this is moab ordering
1376  for( Range::iterator pit = adjPrimaryEnts.begin(); pit != adjPrimaryEnts.end(); pit++ )
1377  {
1378  EntityHandle primaryEnt = *pit;
1379  // get global id
1380  /*int globalID;
1381  rval = context.MBI->tag_get_data(gtags[3], &primaryEnt, 1, &globalID);
1382  if (MB_SUCCESS!=rval)
1383  return moab::MB_FAILURE;
1384  global_element_ID[index] = globalID;*/
1385 
1386  // get local element id
1387  local_element_ID[index] = data.primary_elems.index( primaryEnt );
1388 
1389  if( -1 == local_element_ID[index] )
1390  {
1391  return moab::MB_FAILURE;
1392  } // did not find the element locally
1393 
1394  int side_number, sense, offset;
1395  rval = context.MBI->side_number( primaryEnt, subent, side_number, sense, offset );MB_CHK_ERR( rval );
1396 
1397  reference_surface_ID[index] = side_number + 1; // moab is from 0 to 5, it needs 1 to 6
1398  boundary_condition_value[index] = neuVal;
1399  index++;
1400  }
1401  }
1402  }
1403 
1404  if( index != *surface_BC_length )
1405  {
1406  return moab::MB_FAILURE;
1407  } // error in array allocations
1408 
1409  return moab::MB_SUCCESS;
1410 }

References GlobalContext::appDatas, moab::Range::begin(), context, appData::dimension, moab::Range::end(), ErrorCode, moab::Interface::get_adjacencies(), moab::Interface::get_entities_by_dimension(), moab::Range::index(), MB_CHK_ERR, MB_SUCCESS, GlobalContext::MBI, appData::neu_sets, GlobalContext::neumann_tag, appData::primary_elems, moab::Interface::side_number(), moab::side_number(), moab::Range::size(), and moab::Interface::tag_get_data().

◆ iMOAB_GetPointerToVertexBC()

ErrCode iMOAB_GetPointerToVertexBC ( iMOAB_AppID  pid,
int *  vertex_BC_length,
iMOAB_LocalID local_vertex_ID,
int *  boundary_condition_value 
)

Get the vertex boundary condition information.

Note
Operations: Not Collective
Parameters
[in]pid(iMOAB_AppID) The unique pointer to the application ID.
[in]vertex_BC_length(int) The allocated size of vertex boundary condition array, same as num_visible_vertexBC returned by GetMeshInfo().
[out]local_vertex_ID(iMOAB_LocalID*) The local vertex ID that has Dirichlet BC defined.
[out]boundary_condition_value(int*) The boundary condition type as obtained from the mesh description (value of the Dirichlet_Set tag defined on the vertex).
Returns
ErrCode The error code indicating success or failure.

Definition at line 1412 of file iMOAB.cpp.

1416 {
1417  ErrorCode rval;
1418 
1419  // it was counted above, in GetMeshInfo
1420  appData& data = context.appDatas[*pid];
1421  int numDiriSets = (int)data.diri_sets.size();
1422  int index = 0; // index [0, *vertex_BC_length) for the arrays returned
1423 
1424  for( int i = 0; i < numDiriSets; i++ )
1425  {
1426  Range verts;
1427  EntityHandle diset = data.diri_sets[i];
1428  rval = context.MBI->get_entities_by_dimension( diset, 0, verts );MB_CHK_ERR( rval );
1429 
1430  int diriVal;
1431  rval = context.MBI->tag_get_data( context.dirichlet_tag, &diset, 1, &diriVal );MB_CHK_ERR( rval );
1432 
1433  for( Range::iterator vit = verts.begin(); vit != verts.end(); ++vit )
1434  {
1435  EntityHandle vt = *vit;
1436  /*int vgid;
1437  rval = context.MBI->tag_get_data(gtags[3], &vt, 1, &vgid);
1438  if (MB_SUCCESS!=rval)
1439  return moab::MB_FAILURE;
1440  global_vertext_ID[index] = vgid;*/
1441  local_vertex_ID[index] = data.all_verts.index( vt );
1442 
1443  if( -1 == local_vertex_ID[index] )
1444  {
1445  return moab::MB_FAILURE;
1446  } // vertex was not found
1447 
1448  boundary_condition_value[index] = diriVal;
1449  index++;
1450  }
1451  }
1452 
1453  if( *vertex_BC_length != index )
1454  {
1455  return moab::MB_FAILURE;
1456  } // array allocation issue
1457 
1458  return moab::MB_SUCCESS;
1459 }

References appData::all_verts, GlobalContext::appDatas, moab::Range::begin(), context, appData::diri_sets, GlobalContext::dirichlet_tag, moab::Range::end(), ErrorCode, moab::Interface::get_entities_by_dimension(), moab::Range::index(), MB_CHK_ERR, MB_SUCCESS, GlobalContext::MBI, moab::Range::size(), and moab::Interface::tag_get_data().

◆ iMOAB_GetVertexID()

ErrCode iMOAB_GetVertexID ( iMOAB_AppID  pid,
int *  vertices_length,
iMOAB_GlobalID global_vertex_ID 
)

Get the global vertex IDs for all locally visible (owned and shared/ghosted) vertices.

Note
The array should be allocated by the user, sized with the total number of visible vertices from iMOAB_GetMeshInfo method.

Operations: Not collective

Parameters
[in]pid(iMOAB_AppID) The unique pointer to the application ID.
[in]vertices_length(int*) The allocated size of array (typical size := num_visible_vertices).
[out]global_vertex_ID(iMOAB_GlobalID*) The global IDs for all locally visible vertices (array allocated by user).
Returns
ErrCode The error code indicating success or failure.

Definition at line 964 of file iMOAB.cpp.

965 {
966  IMOAB_CHECKPOINTER( vertices_length, 2 );
967  IMOAB_CHECKPOINTER( global_vertex_ID, 3 );
968 
969  const Range& verts = context.appDatas[*pid].all_verts;
970  // check for problems with array length
971  IMOAB_ASSERT( *vertices_length == static_cast< int >( verts.size() ), "Invalid vertices length provided" );
972 
973  // global id tag is context.globalID_tag
974  return context.MBI->tag_get_data( context.globalID_tag, verts, global_vertex_ID );
975 }

References GlobalContext::appDatas, context, GlobalContext::globalID_tag, IMOAB_ASSERT, IMOAB_CHECKPOINTER, GlobalContext::MBI, moab::Range::size(), and moab::Interface::tag_get_data().

◆ iMOAB_GetVertexOwnership()

ErrCode iMOAB_GetVertexOwnership ( iMOAB_AppID  pid,
int *  vertices_length,
int *  visible_global_rank_ID 
)

Get vertex ownership information.

For each vertex based on the local ID, return the process that owns the vertex (local, shared or ghost)

Note
Shared vertices could be owned by different tasks. Local and shared vertices are first, ghost vertices are next in the array. Ghost vertices are always owned by a different process ID. Array allocated by the user with total size of visible vertices.

Operations: Not Collective

Parameters
[in]pid(iMOAB_AppID) The unique pointer to the application ID.
[in]vertices_length(int*) The allocated size of array (typically size := num_visible_vertices).
[out]visible_global_rank_ID(int*) The processor rank owning each of the local vertices.
Returns
ErrCode The error code indicating success or failure.

Definition at line 977 of file iMOAB.cpp.

978 {
979  assert( vertices_length && *vertices_length );
980  assert( visible_global_rank_ID );
981 
982  Range& verts = context.appDatas[*pid].all_verts;
983  int i = 0;
984 #ifdef MOAB_HAVE_MPI
985  ParallelComm* pco = context.pcomms[*pid];
986 
987  for( Range::iterator vit = verts.begin(); vit != verts.end(); vit++, i++ )
988  {
989  ErrorCode rval = pco->get_owner( *vit, visible_global_rank_ID[i] );MB_CHK_ERR( rval );
990  }
991 
992  if( i != *vertices_length )
993  {
994  return moab::MB_FAILURE;
995  } // warning array allocation problem
996 
997 #else
998 
999  /* everything owned by proc 0 */
1000  if( (int)verts.size() != *vertices_length )
1001  {
1002  return moab::MB_FAILURE;
1003  } // warning array allocation problem
1004 
1005  for( Range::iterator vit = verts.begin(); vit != verts.end(); vit++, i++ )
1006  {
1007  visible_global_rank_ID[i] = 0;
1008  } // all vertices are owned by processor 0, as this is serial run
1009 
1010 #endif
1011 
1012  return moab::MB_SUCCESS;
1013 }

References GlobalContext::appDatas, moab::Range::begin(), context, moab::Range::end(), ErrorCode, moab::ParallelComm::get_owner(), MB_CHK_ERR, MB_SUCCESS, and moab::Range::size().

◆ iMOAB_GetVisibleElementsInfo()

ErrCode iMOAB_GetVisibleElementsInfo ( iMOAB_AppID  pid,
int *  num_visible_elements,
iMOAB_GlobalID element_global_IDs,
int *  ranks,
iMOAB_GlobalID block_IDs 
)

Get the visible elements information.

Return for all visible elements the global IDs, ranks they belong to, block ids they belong to.

Note
Operations: Not Collective
Parameters
[in]pid(iMOAB_AppID) The unique pointer to the application ID.
[in]num_visible_elements(int*) The number of visible elements (returned by GetMeshInfo).
[out]element_global_IDs(iMOAB_GlobalID*) Element global ids to be added to the block.
[out]ranks(int*) The owning ranks of elements.
[out]block_IDs(iMOAB_GlobalID*) The block ids containing the elements.
Returns
ErrCode The error code indicating success or failure.

Definition at line 1096 of file iMOAB.cpp.

1101 {
1102  appData& data = context.appDatas[*pid];
1103 #ifdef MOAB_HAVE_MPI
1104  ParallelComm* pco = context.pcomms[*pid];
1105 #endif
1106 
1107  ErrorCode rval = context.MBI->tag_get_data( context.globalID_tag, data.primary_elems, element_global_IDs );MB_CHK_ERR( rval );
1108 
1109  int i = 0;
1110 
1111  for( Range::iterator eit = data.primary_elems.begin(); eit != data.primary_elems.end(); ++eit, ++i )
1112  {
1113 #ifdef MOAB_HAVE_MPI
1114  rval = pco->get_owner( *eit, ranks[i] );MB_CHK_ERR( rval );
1115 
1116 #else
1117  /* everything owned by task 0 */
1118  ranks[i] = 0;
1119 #endif
1120  }
1121 
1122  for( Range::iterator mit = data.mat_sets.begin(); mit != data.mat_sets.end(); ++mit )
1123  {
1124  EntityHandle matMeshSet = *mit;
1125  Range elems;
1126  rval = context.MBI->get_entities_by_handle( matMeshSet, elems );MB_CHK_ERR( rval );
1127 
1128  int valMatTag;
1129  rval = context.MBI->tag_get_data( context.material_tag, &matMeshSet, 1, &valMatTag );MB_CHK_ERR( rval );
1130 
1131  for( Range::iterator eit = elems.begin(); eit != elems.end(); ++eit )
1132  {
1133  EntityHandle eh = *eit;
1134  int index = data.primary_elems.index( eh );
1135 
1136  if( -1 == index )
1137  {
1138  return moab::MB_FAILURE;
1139  }
1140 
1141  if( -1 >= *num_visible_elements )
1142  {
1143  return moab::MB_FAILURE;
1144  }
1145 
1146  block_IDs[index] = valMatTag;
1147  }
1148  }
1149 
1150  return moab::MB_SUCCESS;
1151 }

References GlobalContext::appDatas, moab::Range::begin(), context, moab::Range::end(), ErrorCode, moab::Interface::get_entities_by_handle(), moab::ParallelComm::get_owner(), GlobalContext::globalID_tag, moab::Range::index(), appData::mat_sets, GlobalContext::material_tag, MB_CHK_ERR, MB_SUCCESS, GlobalContext::MBI, appData::primary_elems, and moab::Interface::tag_get_data().

◆ iMOAB_GetVisibleVerticesCoordinates()

ErrCode iMOAB_GetVisibleVerticesCoordinates ( iMOAB_AppID  pid,
int *  coords_length,
double *  coordinates 
)

Get vertex coordinates for all local (owned and ghosted) vertices.

Note
coordinates are returned in an array allocated by user, interleaved. (do need an option for blocked coordinates ?) size of the array is dimension times number of visible vertices. The local ordering is implicit, owned/shared vertices are first, then ghosts.

Operations: Not Collective

Parameters
[in]pid(iMOAB_AppID) The unique pointer to the application ID.
[in]coords_length(int*) The size of the allocated coordinate array (array allocated by user, size := 3*num_visible_vertices).
[out]coordinates(double*) The pointer to user allocated memory that will be filled with interleaved coordinates.
Returns
ErrCode The error code indicating success or failure.

Definition at line 1015 of file iMOAB.cpp.

1016 {
1017  Range& verts = context.appDatas[*pid].all_verts;
1018 
1019  // interleaved coordinates, so that means deep copy anyway
1020  if( *coords_length != 3 * (int)verts.size() )
1021  {
1022  return moab::MB_FAILURE;
1023  }
1024 
1025  ErrorCode rval = context.MBI->get_coords( verts, coordinates );MB_CHK_ERR( rval );
1026 
1027  return moab::MB_SUCCESS;
1028 }

References GlobalContext::appDatas, context, ErrorCode, moab::Interface::get_coords(), MB_CHK_ERR, MB_SUCCESS, GlobalContext::MBI, and moab::Range::size().

◆ iMOAB_Initialize()

ErrCode iMOAB_Initialize ( int  argc,
iMOAB_String argv 
)

Initialize the iMOAB interface implementation.

Note
Will create the MOAB instance, if not created already (reference counted).

Operations: Collective

Parameters
[in]argc(int) Number of command line arguments
[in]argv(iMOAB_String*) Command line arguments
Returns
ErrCode

Definition at line 135 of file iMOAB.cpp.

136 {
137  if( argc ) IMOAB_CHECKPOINTER( argv, 1 );
138 
139  context.iArgc = argc;
140  context.iArgv = argv; // shallow copy
141 
142  if( 0 == context.refCountMB )
143  {
144  context.MBI = new( std::nothrow ) moab::Core;
145  // retrieve the default tags
146  const char* const shared_set_tag_names[] = { MATERIAL_SET_TAG_NAME, NEUMANN_SET_TAG_NAME,
148  // blocks, visible surfaceBC(neumann), vertexBC (Dirichlet), global id, parallel partition
149  Tag gtags[4];
150 
151  for( int i = 0; i < 4; i++ )
152  {
153 
154  ErrorCode rval =
155  context.MBI->tag_get_handle( shared_set_tag_names[i], 1, MB_TYPE_INTEGER, gtags[i], MB_TAG_ANY );MB_CHK_ERR( rval );
156  }
157 
158  context.material_tag = gtags[0];
159  context.neumann_tag = gtags[1];
160  context.dirichlet_tag = gtags[2];
161  context.globalID_tag = gtags[3];
162  }
163 
164  context.MPI_initialized = false;
165 #ifdef MOAB_HAVE_MPI
166  int flagInit;
167  MPI_Initialized( &flagInit );
168 
169  if( flagInit && !context.MPI_initialized )
170  {
171  MPI_Comm_size( MPI_COMM_WORLD, &context.worldprocs );
172  MPI_Comm_rank( MPI_COMM_WORLD, &context.globalrank );
173  context.MPI_initialized = true;
174  }
175 #endif
176 
178  return moab::MB_SUCCESS;
179 }

References context, DIRICHLET_SET_TAG_NAME, GlobalContext::dirichlet_tag, ErrorCode, GLOBAL_ID_TAG_NAME, GlobalContext::globalID_tag, GlobalContext::globalrank, GlobalContext::iArgc, GlobalContext::iArgv, IMOAB_CHECKPOINTER, MATERIAL_SET_TAG_NAME, GlobalContext::material_tag, MB_CHK_ERR, MB_SUCCESS, MB_TAG_ANY, MB_TYPE_INTEGER, GlobalContext::MBI, GlobalContext::MPI_initialized, NEUMANN_SET_TAG_NAME, GlobalContext::neumann_tag, GlobalContext::refCountMB, moab::Interface::tag_get_handle(), and GlobalContext::worldprocs.

Referenced by iMOAB_InitializeFortran().

◆ iMOAB_InitializeFortran()

ErrCode iMOAB_InitializeFortran ( )

Initialize the iMOAB interface implementation from Fortran driver.

It will create the MOAB instance, if not created already (reference counted).

Operations: Collective

Returns
ErrCode The error code indicating success or failure.

Definition at line 181 of file iMOAB.cpp.

182 {
183  return iMOAB_Initialize( 0, 0 );
184 }

References iMOAB_Initialize().

◆ iMOAB_LoadMesh()

ErrCode iMOAB_LoadMesh ( iMOAB_AppID  pid,
const iMOAB_String  filename,
const iMOAB_String  read_options,
int *  num_ghost_layers 
)

Load a MOAB mesh file in parallel and exchange ghost layers as requested.

Note
All communication is MPI-based, and read options include parallel loading information, resolving shared entities. Local MOAB instance is populated with mesh cells and vertices in the corresponding local partitions.

This will also exchange ghost cells and vertices, as requested. The default bridge dimension is 0 (vertices), and all additional lower dimensional sub-entities are exchanged (mesh edges and faces). The tags in the file are not exchanged by default. Default tag information for GLOBAL_ID, MATERIAL_SET, NEUMANN_SET and DIRICHLET_SET is exchanged. Global ID tag is exchanged for all cells and vertices. Material sets, Neumann sets and Dirichlet sets are all augmented with the ghost entities.

Operations: Collective

Parameters
[in]pid(iMOAB_AppID) The unique pointer to the application ID.
[in]filename(iMOAB_String) The MOAB mesh file (H5M) to load onto the internal application mesh set.
[in]read_options(iMOAB_String) Additional options for reading the MOAB mesh file in parallel.
[in]num_ghost_layers(int*) The total number of ghost layers to exchange during mesh loading.
Returns
ErrCode The error code indicating success or failure.

Definition at line 605 of file iMOAB.cpp.

609 {
610  IMOAB_CHECKPOINTER( filename, 2 );
611  IMOAB_ASSERT( strlen( filename ), "Invalid filename length." );
612 
613  // make sure we use the file set and pcomm associated with the *pid
614  std::ostringstream newopts;
615  if( read_options ) newopts << read_options;
616 
617 #ifdef MOAB_HAVE_MPI
618 
620  {
621  if( context.worldprocs > 1 )
622  {
623  std::string opts( ( read_options ? read_options : "" ) );
624  std::string pcid( "PARALLEL_COMM=" );
625  std::size_t found = opts.find( pcid );
626 
627  if( found != std::string::npos )
628  {
629  std::cerr << " cannot specify PARALLEL_COMM option, it is implicit \n";
630  return moab::MB_FAILURE;
631  }
632 
633  // in serial, apply PARALLEL_COMM option only for h5m files; it does not work for .g
634  // files (used in test_remapping)
635  std::string filen( filename );
636  std::string::size_type idx = filen.rfind( '.' );
637 
638  if( idx != std::string::npos )
639  {
640  std::string extension = filen.substr( idx + 1 );
641  if( ( extension == std::string( "h5m" ) ) || ( extension == std::string( "nc" ) ) )
642  newopts << ";;PARALLEL_COMM=" << *pid;
643  }
644 
645  if( *num_ghost_layers >= 1 )
646  {
647  // if we want ghosts, we will want additional entities, the last .1
648  // because the addl ents can be edges, faces that are part of the neumann sets
649  std::string pcid2( "PARALLEL_GHOSTS=" );
650  std::size_t found2 = opts.find( pcid2 );
651 
652  if( found2 != std::string::npos )
653  {
654  std::cout << " PARALLEL_GHOSTS option is already specified, ignore passed "
655  "number of layers \n";
656  }
657  else
658  {
659  // dimension of primary entities is 3 here, but it could be 2 for climate
660  // meshes; we would need to pass PARALLEL_GHOSTS explicitly for 2d meshes, for
661  // example: ";PARALLEL_GHOSTS=2.0.1"
662  newopts << ";PARALLEL_GHOSTS=3.0." << *num_ghost_layers << ".3";
663  }
664  }
665  }
666  }
667 #else
668  IMOAB_ASSERT( *num_ghost_layers == 0, "Cannot provide ghost layers in serial." );
669 #endif
670 
671  // Now let us actually load the MOAB file with the appropriate read options
672  ErrorCode rval = context.MBI->load_file( filename, &context.appDatas[*pid].file_set, newopts.str().c_str() );MB_CHK_ERR( rval );
673 
674 #ifdef VERBOSE
675  // some debugging stuff
676  std::ostringstream outfile;
677 #ifdef MOAB_HAVE_MPI
678  int rank = context.pcomms[*pid]->rank();
679  int nprocs = context.pcomms[*pid]->size();
680  outfile << "TaskMesh_n" << nprocs << "." << rank << ".h5m";
681 #else
682  outfile << "TaskMesh_n1.0.h5m";
683 #endif
684  // the mesh contains ghosts too, but they are not part of mat/neumann set
685  // write in serial the file, to see what tags are missing
686  rval = context.MBI->write_file( outfile.str().c_str() );MB_CHK_ERR( rval ); // everything on current task, written in serial
687 #endif
688 
689  // Update mesh information
690  return iMOAB_UpdateMeshInfo( pid );
691 }

References GlobalContext::appDatas, context, ErrorCode, IMOAB_ASSERT, IMOAB_CHECKPOINTER, iMOAB_UpdateMeshInfo(), moab::Interface::load_file(), MB_CHK_ERR, GlobalContext::MBI, GlobalContext::MPI_initialized, GlobalContext::worldprocs, and moab::Interface::write_file().

◆ iMOAB_ReadHeaderInfo()

ErrCode iMOAB_ReadHeaderInfo ( const iMOAB_String  filename,
int *  num_global_vertices,
int *  num_global_elements,
int *  num_dimension,
int *  num_parts 
)

It should be called on master task only, and information obtained could be broadcasted by the user. It is a fast lookup in the header of the file.

Note
Operations: Not collective
Parameters
[in]filename(iMOAB_String) The MOAB mesh file (H5M) to probe for header information.
[out]num_global_vertices(int*) The total number of vertices in the mesh file.
[out]num_global_elements(int*) The total number of elements (of highest dimension only).
[out]num_dimension(int*) The highest dimension of elements in the mesh (Edge=1, Tri/Quad=2, Tet/Hex/Prism/Pyramid=3).
[out]num_parts(int*) The total number of partitions available in the mesh file, typically partitioned with mbpart during pre-processing.
Returns
ErrCode The error code indicating success or failure.

Definition at line 461 of file iMOAB.cpp.

466 {
467  IMOAB_CHECKPOINTER( filename, 1 );
468  IMOAB_ASSERT( strlen( filename ), "Invalid filename length." );
469 
470 #ifdef MOAB_HAVE_HDF5
471  std::string filen( filename );
472 
473  int edges = 0;
474  int faces = 0;
475  int regions = 0;
476  if( num_global_vertices ) *num_global_vertices = 0;
477  if( num_global_elements ) *num_global_elements = 0;
478  if( num_dimension ) *num_dimension = 0;
479  if( num_parts ) *num_parts = 0;
480 
481  mhdf_FileHandle file;
482  mhdf_Status status;
483  unsigned long max_id;
484  struct mhdf_FileDesc* data;
485 
486  file = mhdf_openFile( filen.c_str(), 0, &max_id, -1, &status );
487 
488  if( mhdf_isError( &status ) )
489  {
490  fprintf( stderr, "%s: %s\n", filename, mhdf_message( &status ) );
491  return moab::MB_FAILURE;
492  }
493 
494  data = mhdf_getFileSummary( file, H5T_NATIVE_ULONG, &status,
495  1 ); // will use extra set info; will get parallel partition tag info too!
496 
497  if( mhdf_isError( &status ) )
498  {
499  fprintf( stderr, "%s: %s\n", filename, mhdf_message( &status ) );
500  return moab::MB_FAILURE;
501  }
502 
503  if( num_dimension ) *num_dimension = data->nodes.vals_per_ent;
504  if( num_global_vertices ) *num_global_vertices = (int)data->nodes.count;
505 
506  for( int i = 0; i < data->num_elem_desc; i++ )
507  {
508  struct mhdf_ElemDesc* el_desc = &( data->elems[i] );
509  struct mhdf_EntDesc* ent_d = &( el_desc->desc );
510 
511  if( 0 == strcmp( el_desc->type, mhdf_EDGE_TYPE_NAME ) )
512  {
513  edges += ent_d->count;
514  }
515 
516  if( 0 == strcmp( el_desc->type, mhdf_TRI_TYPE_NAME ) )
517  {
518  faces += ent_d->count;
519  }
520 
521  if( 0 == strcmp( el_desc->type, mhdf_QUAD_TYPE_NAME ) )
522  {
523  faces += ent_d->count;
524  }
525 
526  if( 0 == strcmp( el_desc->type, mhdf_POLYGON_TYPE_NAME ) )
527  {
528  faces += ent_d->count;
529  }
530 
531  if( 0 == strcmp( el_desc->type, mhdf_TET_TYPE_NAME ) )
532  {
533  regions += ent_d->count;
534  }
535 
536  if( 0 == strcmp( el_desc->type, mhdf_PYRAMID_TYPE_NAME ) )
537  {
538  regions += ent_d->count;
539  }
540 
541  if( 0 == strcmp( el_desc->type, mhdf_PRISM_TYPE_NAME ) )
542  {
543  regions += ent_d->count;
544  }
545 
546  if( 0 == strcmp( el_desc->type, mdhf_KNIFE_TYPE_NAME ) )
547  {
548  regions += ent_d->count;
549  }
550 
551  if( 0 == strcmp( el_desc->type, mdhf_HEX_TYPE_NAME ) )
552  {
553  regions += ent_d->count;
554  }
555 
556  if( 0 == strcmp( el_desc->type, mhdf_POLYHEDRON_TYPE_NAME ) )
557  {
558  regions += ent_d->count;
559  }
560 
561  if( 0 == strcmp( el_desc->type, mhdf_SEPTAHEDRON_TYPE_NAME ) )
562  {
563  regions += ent_d->count;
564  }
565  }
566 
567  if( num_parts ) *num_parts = data->numEntSets[0];
568 
569  // is this required?
570  if( edges > 0 )
571  {
572  if( num_dimension ) *num_dimension = 1; // I don't think it will ever return 1
573  if( num_global_elements ) *num_global_elements = edges;
574  }
575 
576  if( faces > 0 )
577  {
578  if( num_dimension ) *num_dimension = 2;
579  if( num_global_elements ) *num_global_elements = faces;
580  }
581 
582  if( regions > 0 )
583  {
584  if( num_dimension ) *num_dimension = 3;
585  if( num_global_elements ) *num_global_elements = regions;
586  }
587 
588  mhdf_closeFile( file, &status );
589 
590  free( data );
591 
592 #else
593  std::cout << filename
594  << ": Please reconfigure with HDF5. Cannot retrieve header information for file "
595  "formats other than a h5m file.\n";
596  if( num_global_vertices ) *num_global_vertices = 0;
597  if( num_global_elements ) *num_global_elements = 0;
598  if( num_dimension ) *num_dimension = 0;
599  if( num_parts ) *num_parts = 0;
600 #endif
601 
602  return moab::MB_SUCCESS;
603 }

References mhdf_EntDesc::count, mhdf_ElemDesc::desc, mhdf_FileDesc::elems, IMOAB_ASSERT, IMOAB_CHECKPOINTER, MB_SUCCESS, mdhf_HEX_TYPE_NAME, mdhf_KNIFE_TYPE_NAME, mhdf_closeFile(), mhdf_EDGE_TYPE_NAME, mhdf_getFileSummary(), mhdf_isError(), mhdf_message(), mhdf_openFile(), mhdf_POLYGON_TYPE_NAME, mhdf_POLYHEDRON_TYPE_NAME, mhdf_PRISM_TYPE_NAME, mhdf_PYRAMID_TYPE_NAME, mhdf_QUAD_TYPE_NAME, mhdf_SEPTAHEDRON_TYPE_NAME, mhdf_TET_TYPE_NAME, mhdf_TRI_TYPE_NAME, mhdf_FileDesc::nodes, mhdf_FileDesc::num_elem_desc, mhdf_FileDesc::numEntSets, mhdf_ElemDesc::type, and mhdf_EntDesc::vals_per_ent.

◆ iMOAB_ReduceTagsMax()

ErrCode iMOAB_ReduceTagsMax ( iMOAB_AppID  pid,
int *  tag_index,
int *  ent_type 
)

Perform global reductions with the processes in the current applications communicator. Specifically this routine performs reductions on the maximum value of the tag indicated by tag_index.

Note
Operations: Collective
Parameters
[in]pid(iMOAB_AppID) The unique pointer to the application ID.
[in]tag_index(int*) The tag index of interest.
[in]ent_type(int*) The type of entity for tag reduction operation (vertices = 0, elements = 1)
Returns
ErrCode The error code indicating success or failure.

Definition at line 2126 of file iMOAB.cpp.

2127 {
2128 #ifdef MOAB_HAVE_MPI
2129  appData& data = context.appDatas[*pid];
2130  Range ent_exchange;
2131 
2132  if( *tag_index < 0 || *tag_index >= (int)data.tagList.size() )
2133  {
2134  return moab::MB_FAILURE;
2135  } // error in tag index
2136 
2137  Tag tagh = data.tagList[*tag_index];
2138 
2139  if( *ent_type == 0 )
2140  {
2141  ent_exchange = data.all_verts;
2142  }
2143  else if( *ent_type == 1 )
2144  {
2145  ent_exchange = data.primary_elems;
2146  }
2147  else
2148  {
2149  return moab::MB_FAILURE;
2150  } // unexpected type
2151 
2152  ParallelComm* pco = context.pcomms[*pid];
2153  // we could do different MPI_Op; do not bother now, we will call from fortran
2154  ErrorCode rval = pco->reduce_tags( tagh, MPI_MAX, ent_exchange );MB_CHK_ERR( rval );
2155 
2156 #else
2157  /* do nothing if serial */
2158  // just silence the warning
2159  // do not call sync tags in serial!
2160  int k = *pid + *tag_index + *ent_type;
2161  k++; // just do junk, to avoid complaints
2162 #endif
2163  return moab::MB_SUCCESS;
2164 }

References appData::all_verts, GlobalContext::appDatas, context, ErrorCode, MB_CHK_ERR, MB_SUCCESS, appData::primary_elems, moab::ParallelComm::reduce_tags(), and appData::tagList.

◆ iMOAB_RegisterApplication()

ErrCode iMOAB_RegisterApplication ( const iMOAB_String  app_name,
int *  compid,
iMOAB_AppID  pid 
)

Register application - Create a unique application ID and bootstrap interfaces for further queries.

Note
Internally, a mesh set will be associated with the application ID and all subsequent queries on the MOAB instance will be directed to this mesh/file set.

Operations: Collective

Parameters
[in]app_name(iMOAB_String) Application name (PROTEUS, NEK5000, etc)
[in]comm(MPI_Comm*) MPI communicator to be used for all mesh-related queries originating from this application
[in]compid(int*) The unique external component identifier
[out]pid(iMOAB_AppID) The unique pointer to the application ID
Returns
ErrCode The error code indicating success or failure.

Definition at line 198 of file iMOAB.cpp.

204 {
205  IMOAB_CHECKPOINTER( app_name, 1 );
206 #ifdef MOAB_HAVE_MPI
207  IMOAB_CHECKPOINTER( comm, 2 );
208  IMOAB_CHECKPOINTER( compid, 3 );
209 #else
210  IMOAB_CHECKPOINTER( compid, 2 );
211 #endif
212 
213  // will create a parallel comm for this application too, so there will be a
214  // mapping from *pid to file set and to parallel comm instances
215  std::string name( app_name );
216 
217  if( context.appIdMap.find( name ) != context.appIdMap.end() )
218  {
219  std::cout << " application " << name << " already registered \n";
220  return moab::MB_FAILURE;
221  }
222 
223  *pid = context.unused_pid++;
224  context.appIdMap[name] = *pid;
225  int rankHere = 0;
226 #ifdef MOAB_HAVE_MPI
227  MPI_Comm_rank( *comm, &rankHere );
228 #endif
229  if( !rankHere ) std::cout << " application " << name << " with ID = " << *pid << " is registered now \n";
230  if( *compid <= 0 )
231  {
232  std::cout << " convention for external application is to have its id positive \n";
233  return moab::MB_FAILURE;
234  }
235 
236  if( context.appIdCompMap.find( *compid ) != context.appIdCompMap.end() )
237  {
238  std::cout << " external application with comp id " << *compid << " is already registered\n";
239  return moab::MB_FAILURE;
240  }
241 
242  context.appIdCompMap[*compid] = *pid;
243 
244  // now create ParallelComm and a file set for this application
245 #ifdef MOAB_HAVE_MPI
246  if( *comm )
247  {
248  ParallelComm* pco = new ParallelComm( context.MBI, *comm );
249 
250 #ifndef NDEBUG
251  int index = pco->get_id(); // it could be useful to get app id from pcomm instance ...
252  assert( index == *pid );
253  // here, we assert the the pid is the same as the id of the ParallelComm instance
254  // useful for writing in parallel
255 #endif
256  context.pcomms.push_back( pco );
257  }
258  else
259  {
260  context.pcomms.push_back( 0 );
261  }
262 #endif
263 
264  // create now the file set that will be used for loading the model in
265  EntityHandle file_set;
266  ErrorCode rval = context.MBI->create_meshset( MESHSET_SET, file_set );MB_CHK_ERR( rval );
267 
268  appData app_data;
269  app_data.file_set = file_set;
270  app_data.global_id = *compid; // will be used mostly for par comm graph
271  app_data.name = name; // save the name of application
272 
273 #ifdef MOAB_HAVE_TEMPESTREMAP
274  app_data.tempestData.remapper = NULL; // Only allocate as needed
275 #endif
276 
277  app_data.point_cloud = false;
278  app_data.is_fortran = false;
279 
280  context.appDatas.push_back(
281  app_data ); // it will correspond to app_FileSets[*pid] will be the file set of interest
282  return moab::MB_SUCCESS;
283 }

References GlobalContext::appDatas, GlobalContext::appIdCompMap, GlobalContext::appIdMap, context, moab::Interface::create_meshset(), ErrorCode, appData::file_set, moab::ParallelComm::get_id(), appData::global_id, IMOAB_CHECKPOINTER, appData::is_fortran, MB_CHK_ERR, MB_SUCCESS, GlobalContext::MBI, MESHSET_SET, appData::name, appData::point_cloud, and GlobalContext::unused_pid.

Referenced by iMOAB_RegisterApplicationFortran().

◆ iMOAB_RegisterApplicationFortran()

ErrCode iMOAB_RegisterApplicationFortran ( const iMOAB_String  app_name,
int *  compid,
iMOAB_AppID  pid 
)

Register a Fortran-based application - Create a unique application ID and bootstrap interfaces for further queries.

Note
Internally, the comm object, usually a 32 bit integer in Fortran, will be converted using MPI_Comm_f2c and stored as MPI_Comm. A mesh set will be associated with the application ID and all subsequent queries on the MOAB instance will be directed to this mesh/file set.

Operations: Collective

Parameters
[in]app_name(iMOAB_String) Application name ("MySolver", "E3SM", "DMMoab", "PROTEUS", "NEK5000", etc)
[in]communicator(int*) Fortran MPI communicator to be used for all mesh-related queries originating from this application.
[in]compid(int*) External component id, (A const unique identifier).
[out]pid(iMOAB_AppID) The unique pointer to the application ID.
Returns
ErrCode The error code indicating success or failure.

Definition at line 285 of file iMOAB.cpp.

291 {
292  IMOAB_CHECKPOINTER( app_name, 1 );
293 #ifdef MOAB_HAVE_MPI
294  IMOAB_CHECKPOINTER( comm, 2 );
295  IMOAB_CHECKPOINTER( compid, 3 );
296 #else
297  IMOAB_CHECKPOINTER( compid, 2 );
298 #endif
299 
300  ErrCode err;
301  assert( app_name != nullptr );
302  std::string name( app_name );
303 
304 #ifdef MOAB_HAVE_MPI
305  MPI_Comm ccomm;
306  if( comm )
307  {
308  // convert from Fortran communicator to a C communicator
309  // see transfer of handles
310  // http://www.mpi-forum.org/docs/mpi-2.2/mpi22-report/node361.htm
311  ccomm = MPI_Comm_f2c( (MPI_Fint)*comm );
312  }
313 #endif
314 
315  // now call C style registration function:
316  err = iMOAB_RegisterApplication( app_name,
317 #ifdef MOAB_HAVE_MPI
318  &ccomm,
319 #endif
320  compid, pid );
321 
322  // Now that we have created the application context, store that
323  // the application being registered is from a Fortran context
324  context.appDatas[*pid].is_fortran = true;
325 
326  return err;
327 }

References GlobalContext::appDatas, context, ErrCode, IMOAB_CHECKPOINTER, and iMOAB_RegisterApplication().

◆ iMOAB_SetDoubleTagStorage()

ErrCode iMOAB_SetDoubleTagStorage ( iMOAB_AppID  pid,
const iMOAB_String  tag_storage_name,
int *  num_tag_storage_length,
int *  entity_type,
double *  tag_storage_data 
)

Store the specified values in a MOAB double Tag.

Note
Operations: Collective
Parameters
[in]pid(iMOAB_AppID) The unique pointer to the application ID.
[in]tag_storage_name(iMOAB_String) The tag names, separated, to store the data.
[in]num_tag_storage_length(int*) The size of total tag storage data (e.g., num_visible_vertices*components_per_entity or num_visible_elements*components_per_entity*num_tags).
[in]entity_type(int*) Type=0 for vertices, and Type=1 for primary elements.
[out]tag_storage_data(double*) The array data of type double to replace the internal tag memory; The data is assumed to be contiguous over the local set of visible entities (either vertices or elements). unrolled by tags
Returns
ErrCode The error code indicating success or failure.

Definition at line 1710 of file iMOAB.cpp.

1715 {
1716  ErrorCode rval;
1717  std::string tag_names( tag_storage_names );
1718  // exactly the same code as for int tag :) maybe should check the type of tag too
1719  std::vector< std::string > tagNames;
1720  std::vector< Tag > tagHandles;
1721  std::string separator( ":" );
1722  split_tag_names( tag_names, separator, tagNames );
1723 
1724  appData& data = context.appDatas[*pid];
1725  Range* ents_to_set = NULL;
1726 
1727  if( *ent_type == 0 ) // vertices
1728  {
1729  ents_to_set = &data.all_verts;
1730  }
1731  else if( *ent_type == 1 )
1732  {
1733  ents_to_set = &data.primary_elems;
1734  }
1735 
1736  int nents_to_be_set = (int)( *ents_to_set ).size();
1737  int position = 0;
1738 
1739  for( size_t i = 0; i < tagNames.size(); i++ )
1740  {
1741  if( data.tagMap.find( tagNames[i] ) == data.tagMap.end() )
1742  {
1743  return moab::MB_FAILURE;
1744  } // some tag not defined yet in the app
1745 
1746  Tag tag = data.tagMap[tagNames[i]];
1747 
1748  int tagLength = 0;
1749  rval = context.MBI->tag_get_length( tag, tagLength );MB_CHK_ERR( rval );
1750 
1751  DataType dtype;
1752  rval = context.MBI->tag_get_data_type( tag, dtype );MB_CHK_ERR( rval );
1753 
1754  if( dtype != MB_TYPE_DOUBLE )
1755  {
1756  return moab::MB_FAILURE;
1757  }
1758 
1759  // set it on the subset of entities, based on type and length
1760  if( position + tagLength * nents_to_be_set > *num_tag_storage_length )
1761  return moab::MB_FAILURE; // too many entity values to be set
1762 
1763  rval = context.MBI->tag_set_data( tag, *ents_to_set, &tag_storage_data[position] );MB_CHK_ERR( rval );
1764  // increment position to next tag
1765  position = position + tagLength * nents_to_be_set;
1766  }
1767  return moab::MB_SUCCESS; // no error
1768 }

References appData::all_verts, GlobalContext::appDatas, context, ErrorCode, MB_CHK_ERR, MB_SUCCESS, MB_TYPE_DOUBLE, GlobalContext::MBI, appData::primary_elems, moab::Range::size(), split_tag_names(), moab::Interface::tag_get_data_type(), moab::Interface::tag_get_length(), moab::Interface::tag_set_data(), and appData::tagMap.

◆ iMOAB_SetDoubleTagStorageWithGid()

ErrCode iMOAB_SetDoubleTagStorageWithGid ( iMOAB_AppID  pid,
const iMOAB_String  tag_storage_names,
int *  num_tag_storage_length,
int *  entity_type,
double *  tag_storage_data,
int *  globalIds 
)

Store the specified values in a MOAB double Tag, for given ids.

Note
Operations: Collective
Parameters
[in]pid(iMOAB_AppID) The unique pointer to the application ID.
[in]tag_storage_name(iMOAB_String) The tag names, separated, to store the data.
[in]num_tag_storage_length(int*) The size of total tag storage data (e.g., num_visible_vertices*components_per_entity or num_visible_elements*components_per_entity*num_tags).
[in]entity_type(int*) Type=0 for vertices, and Type=1 for primary elements.
[in]tag_storage_data(double*) The array data of type double to replace the internal tag memory; The data is assumed to be permuted over the local set of visible entities (either vertices or elements). unrolled by tags in parallel, might be different order, or different entities on the task
[in]globalIdsglobal ids of the cells to be set;
Returns
ErrCode The error code indicating success or failure.

Definition at line 1770 of file iMOAB.cpp.

1776 {
1777  ErrorCode rval;
1778  std::string tag_names( tag_storage_names );
1779  // exactly the same code as for int tag :) maybe should check the type of tag too
1780  std::vector< std::string > tagNames;
1781  std::vector< Tag > tagHandles;
1782  std::string separator( ":" );
1783  split_tag_names( tag_names, separator, tagNames );
1784 
1785  appData& data = context.appDatas[*pid];
1786  Range* ents_to_set = NULL;
1787 
1788  if( *ent_type == 0 ) // vertices
1789  {
1790  ents_to_set = &data.all_verts;
1791  }
1792  else if( *ent_type == 1 )
1793  {
1794  ents_to_set = &data.primary_elems;
1795  }
1796 
1797  int nents_to_be_set = (int)( *ents_to_set ).size();
1798 
1799  Tag gidTag = context.MBI->globalId_tag();
1800  std::vector< int > gids;
1801  gids.resize( nents_to_be_set );
1802  rval = context.MBI->tag_get_data( gidTag, *ents_to_set, &gids[0] );MB_CHK_ERR( rval );
1803 
1804  // so we will need to set the tags according to the global id passed;
1805  // so the order in tag_storage_data is the same as the order in globalIds, but the order
1806  // in local range is gids
1807  std::map< int, EntityHandle > eh_by_gid;
1808  int i = 0;
1809  for( Range::iterator it = ents_to_set->begin(); it != ents_to_set->end(); ++it, ++i )
1810  {
1811  eh_by_gid[gids[i]] = *it;
1812  }
1813 
1814  std::vector< int > tagLengths( tagNames.size() );
1815  std::vector< Tag > tagList;
1816  size_t total_tag_len = 0;
1817  for( size_t i = 0; i < tagNames.size(); i++ )
1818  {
1819  if( data.tagMap.find( tagNames[i] ) == data.tagMap.end() )
1820  {
1821  MB_SET_ERR( moab::MB_FAILURE, "tag missing" );
1822  } // some tag not defined yet in the app
1823 
1824  Tag tag = data.tagMap[tagNames[i]];
1825  tagList.push_back( tag );
1826 
1827  int tagLength = 0;
1828  rval = context.MBI->tag_get_length( tag, tagLength );MB_CHK_ERR( rval );
1829 
1830  total_tag_len += tagLength;
1831  tagLengths[i] = tagLength;
1832  DataType dtype;
1833  rval = context.MBI->tag_get_data_type( tag, dtype );MB_CHK_ERR( rval );
1834 
1835  if( dtype != MB_TYPE_DOUBLE )
1836  {
1837  MB_SET_ERR( moab::MB_FAILURE, "tag not double type" );
1838  }
1839  }
1840  bool serial = true;
1841 #ifdef MOAB_HAVE_MPI
1842  ParallelComm* pco = context.pcomms[*pid];
1843  unsigned num_procs = pco->size();
1844  if( num_procs > 1 ) serial = false;
1845 #endif
1846 
1847  if( serial )
1848  {
1849  assert( total_tag_len * nents_to_be_set - *num_tag_storage_length == 0 );
1850  // tags are unrolled, we loop over global ids first, then careful about tags
1851  for( int i = 0; i < nents_to_be_set; i++ )
1852  {
1853  int gid = globalIds[i];
1854  EntityHandle eh = eh_by_gid[gid];
1855  // now loop over tags
1856  int indexInTagValues = 0; //
1857  for( size_t j = 0; j < tagList.size(); j++ )
1858  {
1859  indexInTagValues += i * tagLengths[j];
1860  rval = context.MBI->tag_set_data( tagList[j], &eh, 1, &tag_storage_data[indexInTagValues] );MB_CHK_ERR( rval );
1861  // advance the pointer/index
1862  indexInTagValues += ( nents_to_be_set - i ) * tagLengths[j]; // at the end of tag data
1863  }
1864  }
1865  }
1866 #ifdef MOAB_HAVE_MPI
1867  else // it can be not serial only if pco->size() > 1, parallel
1868  {
1869  // in this case, we have to use 2 crystal routers, to send data to the processor that needs it
1870  // we will create first a tuple to rendevous points, then from there send to the processor that requested it
1871  // it is a 2-hop global gather scatter
1872  int nbLocalVals = *num_tag_storage_length / ( (int)tagNames.size() );
1873  assert( nbLocalVals * tagNames.size() - *num_tag_storage_length == 0 );
1874  TupleList TLsend;
1875  TLsend.initialize( 2, 0, 0, total_tag_len, nbLocalVals ); // to proc, marker(gid), total_tag_len doubles
1876  TLsend.enableWriteAccess();
1877  // the processor id that processes global_id is global_id / num_ents_per_proc
1878 
1879  int indexInRealLocal = 0;
1880  for( int i = 0; i < nbLocalVals; i++ )
1881  {
1882  // to proc, marker, element local index, index in el
1883  int marker = globalIds[i];
1884  int to_proc = marker % num_procs;
1885  int n = TLsend.get_n();
1886  TLsend.vi_wr[2 * n] = to_proc; // send to processor
1887  TLsend.vi_wr[2 * n + 1] = marker;
1888  int indexInTagValues = 0;
1889  // tag data collect by number of tags
1890  for( size_t j = 0; j < tagList.size(); j++ )
1891  {
1892  indexInTagValues += i * tagLengths[j];
1893  for( int k = 0; k < tagLengths[j]; k++ )
1894  {
1895  TLsend.vr_wr[indexInRealLocal++] = tag_storage_data[indexInTagValues + k];
1896  }
1897  indexInTagValues += ( nbLocalVals - i ) * tagLengths[j];
1898  }
1899  TLsend.inc_n();
1900  }
1901 
1902  assert( nbLocalVals * total_tag_len - indexInRealLocal == 0 );
1903  // send now requests, basically inform the rendez-vous point who needs a particular global id
1904  // send the data to the other processors:
1905  ( pco->proc_config().crystal_router() )->gs_transfer( 1, TLsend, 0 );
1906  TupleList TLreq;
1907  TLreq.initialize( 2, 0, 0, 0, nents_to_be_set );
1908  TLreq.enableWriteAccess();
1909  for( int i = 0; i < nents_to_be_set; i++ )
1910  {
1911  // to proc, marker
1912  int marker = gids[i];
1913  int to_proc = marker % num_procs;
1914  int n = TLreq.get_n();
1915  TLreq.vi_wr[2 * n] = to_proc; // send to processor
1916  TLreq.vi_wr[2 * n + 1] = marker;
1917  // tag data collect by number of tags
1918  TLreq.inc_n();
1919  }
1920  ( pco->proc_config().crystal_router() )->gs_transfer( 1, TLreq, 0 );
1921 
1922  // we know now that process TLreq.vi_wr[2 * n] needs tags for gid TLreq.vi_wr[2 * n + 1]
1923  // we should first order by global id, and then build the new TL with send to proc, global id and
1924  // tags for it
1925  // sort by global ids the tuple lists
1926  moab::TupleList::buffer sort_buffer;
1927  sort_buffer.buffer_init( TLreq.get_n() );
1928  TLreq.sort( 1, &sort_buffer );
1929  sort_buffer.reset();
1930  sort_buffer.buffer_init( TLsend.get_n() );
1931  TLsend.sort( 1, &sort_buffer );
1932  sort_buffer.reset();
1933  // now send the tag values to the proc that requested it
1934  // in theory, for a full partition, TLreq and TLsend should have the same size, and
1935  // each dof should have exactly one target proc. Is that true or not in general ?
1936  // how do we plan to use this? Is it better to store the comm graph for future
1937 
1938  // start copy from comm graph settle
1939  TupleList TLBack;
1940  TLBack.initialize( 3, 0, 0, total_tag_len, 0 ); // to proc, marker, tag from proc , tag values
1941  TLBack.enableWriteAccess();
1942 
1943  int n1 = TLreq.get_n();
1944  int n2 = TLsend.get_n();
1945 
1946  int indexInTLreq = 0;
1947  int indexInTLsend = 0; // advance both, according to the marker
1948  if( n1 > 0 && n2 > 0 )
1949  {
1950 
1951  while( indexInTLreq < n1 && indexInTLsend < n2 ) // if any is over, we are done
1952  {
1953  int currentValue1 = TLreq.vi_rd[2 * indexInTLreq + 1];
1954  int currentValue2 = TLsend.vi_rd[2 * indexInTLsend + 1];
1955  if( currentValue1 < currentValue2 )
1956  {
1957  // we have a big problem; basically, we are saying that
1958  // dof currentValue is on one model and not on the other
1959  // std::cout << " currentValue1:" << currentValue1 << " missing in comp2" << "\n";
1960  indexInTLreq++;
1961  continue;
1962  }
1963  if( currentValue1 > currentValue2 )
1964  {
1965  // std::cout << " currentValue2:" << currentValue2 << " missing in comp1" << "\n";
1966  indexInTLsend++;
1967  continue;
1968  }
1969  int size1 = 1;
1970  int size2 = 1;
1971  while( indexInTLreq + size1 < n1 && currentValue1 == TLreq.vi_rd[2 * ( indexInTLreq + size1 ) + 1] )
1972  size1++;
1973  while( indexInTLsend + size2 < n2 && currentValue2 == TLsend.vi_rd[2 * ( indexInTLsend + size2 ) + 1] )
1974  size2++;
1975  // must be found in both lists, find the start and end indices
1976  for( int i1 = 0; i1 < size1; i1++ )
1977  {
1978  for( int i2 = 0; i2 < size2; i2++ )
1979  {
1980  // send the info back to components
1981  int n = TLBack.get_n();
1982  TLBack.reserve();
1983  TLBack.vi_wr[3 * n] = TLreq.vi_rd[2 * ( indexInTLreq + i1 )]; // send back to the proc marker
1984  // came from, info from comp2
1985  TLBack.vi_wr[3 * n + 1] = currentValue1; // initial value (resend, just for verif ?)
1986  TLBack.vi_wr[3 * n + 2] = TLsend.vi_rd[2 * ( indexInTLsend + i2 )]; // from proc on comp2
1987  // also fill tag values
1988  for( size_t k = 0; k < total_tag_len; k++ )
1989  {
1990  TLBack.vr_rd[total_tag_len * n + k] =
1991  TLsend.vr_rd[total_tag_len * indexInTLsend + k]; // deep copy of tag values
1992  }
1993  }
1994  }
1995  indexInTLreq += size1;
1996  indexInTLsend += size2;
1997  }
1998  }
1999  ( pco->proc_config().crystal_router() )->gs_transfer( 1, TLBack, 0 );
2000  // end copy from comm graph
2001  // after we are done sending, we need to set those tag values, in a reverse process compared to send
2002  n1 = TLBack.get_n();
2003  double* ptrVal = &TLBack.vr_rd[0]; //
2004  for( int i = 0; i < n1; i++ )
2005  {
2006  int gid = TLBack.vi_rd[3 * i + 1]; // marker
2007  EntityHandle eh = eh_by_gid[gid];
2008  // now loop over tags
2009 
2010  for( size_t j = 0; j < tagList.size(); j++ )
2011  {
2012  rval = context.MBI->tag_set_data( tagList[j], &eh, 1, (void*)ptrVal );MB_CHK_ERR( rval );
2013  // advance the pointer/index
2014  ptrVal += tagLengths[j]; // at the end of tag data per call
2015  }
2016  }
2017  }
2018 #endif
2019  return MB_SUCCESS;
2020 }

References appData::all_verts, GlobalContext::appDatas, moab::Range::begin(), context, moab::ProcConfig::crystal_router(), moab::TupleList::enableWriteAccess(), moab::Range::end(), ErrorCode, moab::TupleList::get_n(), moab::Interface::globalId_tag(), moab::TupleList::inc_n(), moab::TupleList::initialize(), MB_CHK_ERR, MB_SET_ERR, MB_SUCCESS, MB_TYPE_DOUBLE, GlobalContext::MBI, appData::primary_elems, moab::ParallelComm::proc_config(), moab::TupleList::reserve(), moab::TupleList::buffer::reset(), moab::Range::size(), moab::ParallelComm::size(), moab::TupleList::sort(), split_tag_names(), moab::Interface::tag_get_data(), moab::Interface::tag_get_data_type(), moab::Interface::tag_get_length(), moab::Interface::tag_set_data(), appData::tagMap, moab::TupleList::vi_rd, moab::TupleList::vi_wr, moab::TupleList::vr_rd, and moab::TupleList::vr_wr.

◆ iMOAB_SetGlobalInfo()

ErrCode iMOAB_SetGlobalInfo ( iMOAB_AppID  pid,
int *  num_global_verts,
int *  num_global_elems 
)

Set global information for number of vertices and number of elements; It is usually available from the h5m file, or it can be computed with a MPI_Reduce.

Parameters
[in]pid(iMOAB_AppID) The unique pointer to the application ID.
[in]num_global_verts(int*) The number of total vertices in the mesh.
[in]num_global_elems(MPI_Comm) The number of total elements in the mesh.
Returns
ErrCode The error code indicating success or failure.

Definition at line 2290 of file iMOAB.cpp.

2291 {
2292  appData& data = context.appDatas[*pid];
2293  data.num_global_vertices = *num_global_verts;
2294  data.num_global_elements = *num_global_elems;
2295  return moab::MB_SUCCESS;
2296 }

References GlobalContext::appDatas, context, MB_SUCCESS, appData::num_global_elements, and appData::num_global_vertices.

Referenced by iMOAB_UpdateMeshInfo().

◆ iMOAB_SetIntTagStorage()

ErrCode iMOAB_SetIntTagStorage ( iMOAB_AppID  pid,
const iMOAB_String  tag_storage_name,
int *  num_tag_storage_length,
int *  entity_type,
int *  tag_storage_data 
)

Store the specified values in a MOAB integer Tag.

Values are set on vertices or elements, depending on entity_type

Note
we could change the api to accept as input the tag index, as returned by iMOAB_DefineTagStorage.

Operations: Collective

Parameters
[in]pid(iMOAB_AppID) The unique pointer to the application ID.
[in]tag_storage_name(iMOAB_String) The tag name to store/retreive the data in MOAB.
[in]num_tag_storage_length(int*) The size of tag storage data (e.g., num_visible_vertices*components_per_entity or num_visible_elements*components_per_entity).
[in]entity_type(int*) Type=0 for vertices, and Type=1 for primary elements.
[out]tag_storage_data(int*) The array data of type int to replace the internal tag memory; The data is assumed to be contiguous over the local set of visible entities (either vertices or elements).
Returns
ErrCode The error code indicating success or failure.

Definition at line 1596 of file iMOAB.cpp.

1601 {
1602  std::string tag_name( tag_storage_name );
1603 
1604  // Get the application data
1605  appData& data = context.appDatas[*pid];
1606 
1607  if( data.tagMap.find( tag_name ) == data.tagMap.end() )
1608  {
1609  return moab::MB_FAILURE;
1610  } // tag not defined
1611 
1612  Tag tag = data.tagMap[tag_name];
1613 
1614  int tagLength = 0;
1615  ErrorCode rval = context.MBI->tag_get_length( tag, tagLength );MB_CHK_ERR( rval );
1616 
1617  DataType dtype;
1618  rval = context.MBI->tag_get_data_type( tag, dtype );
1619 
1620  if( MB_SUCCESS != rval || dtype != MB_TYPE_INTEGER )
1621  {
1622  return moab::MB_FAILURE;
1623  }
1624 
1625  // set it on a subset of entities, based on type and length
1626  Range* ents_to_set;
1627 
1628  if( *ent_type == 0 ) // vertices
1629  {
1630  ents_to_set = &data.all_verts;
1631  }
1632  else // if (*ent_type == 1) // *ent_type can be 0 (vertices) or 1 (elements)
1633  {
1634  ents_to_set = &data.primary_elems;
1635  }
1636 
1637  int nents_to_be_set = *num_tag_storage_length / tagLength;
1638 
1639  if( nents_to_be_set > (int)ents_to_set->size() || nents_to_be_set < 1 )
1640  {
1641  return moab::MB_FAILURE;
1642  } // to many entities to be set or too few
1643 
1644  // restrict the range; everything is contiguous; or not?
1645 
1646  // Range contig_range ( * ( ents_to_set->begin() ), * ( ents_to_set->begin() + nents_to_be_set -
1647  // 1 ) );
1648  rval = context.MBI->tag_set_data( tag, *ents_to_set, tag_storage_data );MB_CHK_ERR( rval );
1649 
1650  return moab::MB_SUCCESS; // no error
1651 }

References appData::all_verts, GlobalContext::appDatas, context, ErrorCode, MB_CHK_ERR, MB_SUCCESS, MB_TYPE_INTEGER, GlobalContext::MBI, appData::primary_elems, moab::Range::size(), moab::Interface::tag_get_data_type(), moab::Interface::tag_get_length(), moab::Interface::tag_set_data(), and appData::tagMap.

◆ iMOAB_SynchronizeTags()

ErrCode iMOAB_SynchronizeTags ( iMOAB_AppID  pid,
int *  num_tag,
int *  tag_indices,
int *  ent_type 
)

Exchange tag values for the given tags across process boundaries.

Note
Operations: Collective
Parameters
[in]pid(iMOAB_AppID) The unique pointer to the application ID.
[in]num_tags(int*) Number of tags to exchange.
[in]tag_indices(int*) Array with tag indices of interest (size = *num_tags).
[in]ent_type(int*) The type of entity for tag exchange.
Returns
ErrCode The error code indicating success or failure.

Definition at line 2081 of file iMOAB.cpp.

2082 {
2083 #ifdef MOAB_HAVE_MPI
2084  appData& data = context.appDatas[*pid];
2085  Range ent_exchange;
2086  std::vector< Tag > tags;
2087 
2088  for( int i = 0; i < *num_tag; i++ )
2089  {
2090  if( tag_indices[i] < 0 || tag_indices[i] >= (int)data.tagList.size() )
2091  {
2092  return moab::MB_FAILURE;
2093  } // error in tag index
2094 
2095  tags.push_back( data.tagList[tag_indices[i]] );
2096  }
2097 
2098  if( *ent_type == 0 )
2099  {
2100  ent_exchange = data.all_verts;
2101  }
2102  else if( *ent_type == 1 )
2103  {
2104  ent_exchange = data.primary_elems;
2105  }
2106  else
2107  {
2108  return moab::MB_FAILURE;
2109  } // unexpected type
2110 
2111  ParallelComm* pco = context.pcomms[*pid];
2112 
2113  ErrorCode rval = pco->exchange_tags( tags, tags, ent_exchange );MB_CHK_ERR( rval );
2114 
2115 #else
2116  /* do nothing if serial */
2117  // just silence the warning
2118  // do not call sync tags in serial!
2119  int k = *pid + *num_tag + *tag_indices + *ent_type;
2120  k++;
2121 #endif
2122 
2123  return moab::MB_SUCCESS;
2124 }

References appData::all_verts, GlobalContext::appDatas, context, ErrorCode, moab::ParallelComm::exchange_tags(), MB_CHK_ERR, MB_SUCCESS, appData::primary_elems, and appData::tagList.

◆ iMOAB_UpdateMeshInfo()

ErrCode iMOAB_UpdateMeshInfo ( iMOAB_AppID  pid)

Update local mesh data structure, from file information.

Note
The method should be called after mesh modifications, for example reading a file or creating mesh in memory

Operations: Collective

Parameters
[in]pid(iMOAB_AppID) The unique pointer to the application ID.
Returns
ErrCode The error code indicating success or failure.

Definition at line 772 of file iMOAB.cpp.

773 {
774  // this will include ghost elements info
775  appData& data = context.appDatas[*pid];
776  EntityHandle fileSet = data.file_set;
777 
778  // first clear all data ranges; this can be called after ghosting
779  data.all_verts.clear();
780  data.primary_elems.clear();
781  data.local_verts.clear();
782  data.owned_verts.clear();
783  data.ghost_vertices.clear();
784  data.owned_elems.clear();
785  data.ghost_elems.clear();
786  data.mat_sets.clear();
787  data.neu_sets.clear();
788  data.diri_sets.clear();
789 
790  // Let us get all the vertex entities
791  ErrorCode rval = context.MBI->get_entities_by_type( fileSet, MBVERTEX, data.all_verts, true );MB_CHK_ERR( rval ); // recursive
792 
793  // Let us check first entities of dimension = 3
794  data.dimension = 3;
795  rval = context.MBI->get_entities_by_dimension( fileSet, data.dimension, data.primary_elems, true );MB_CHK_ERR( rval ); // recursive
796 
797  if( data.primary_elems.empty() )
798  {
799  // Now 3-D elements. Let us check entities of dimension = 2
800  data.dimension = 2;
801  rval = context.MBI->get_entities_by_dimension( fileSet, data.dimension, data.primary_elems, true );MB_CHK_ERR( rval ); // recursive
802 
803  if( data.primary_elems.empty() )
804  {
805  // Now 3-D/2-D elements. Let us check entities of dimension = 1
806  data.dimension = 1;
807  rval = context.MBI->get_entities_by_dimension( fileSet, data.dimension, data.primary_elems, true );MB_CHK_ERR( rval ); // recursive
808 
809  if( data.primary_elems.empty() )
810  {
811  // no elements of dimension 1 or 2 or 3; it could happen for point clouds
812  data.dimension = 0;
813  }
814  }
815  }
816 
817  // check if the current mesh is just a point cloud
818  data.point_cloud = ( ( data.primary_elems.size() == 0 && data.all_verts.size() > 0 ) || data.dimension == 0 );
819 
820 #ifdef MOAB_HAVE_MPI
821 
823  {
824  ParallelComm* pco = context.pcomms[*pid];
825 
826  // filter ghost vertices, from local
827  rval = pco->filter_pstatus( data.all_verts, PSTATUS_GHOST, PSTATUS_NOT, -1, &data.local_verts );MB_CHK_ERR( rval );
828 
829  // Store handles for all ghosted entities
830  data.ghost_vertices = subtract( data.all_verts, data.local_verts );
831 
832  // filter ghost elements, from local
833  rval = pco->filter_pstatus( data.primary_elems, PSTATUS_GHOST, PSTATUS_NOT, -1, &data.owned_elems );MB_CHK_ERR( rval );
834 
835  data.ghost_elems = subtract( data.primary_elems, data.owned_elems );
836  // now update global number of primary cells and global number of vertices
837  // determine first number of owned vertices
838  // Get local owned vertices
839  rval = pco->filter_pstatus( data.all_verts, PSTATUS_NOT_OWNED, PSTATUS_NOT, -1, &data.owned_verts );MB_CHK_ERR( rval );
840  int local[2], global[2];
841  local[0] = data.owned_verts.size();
842  local[1] = data.owned_elems.size();
843  MPI_Allreduce( local, global, 2, MPI_INT, MPI_SUM, pco->comm() );
844  rval = iMOAB_SetGlobalInfo( pid, &( global[0] ), &( global[1] ) );MB_CHK_ERR( rval );
845  }
846  else
847  {
848  data.local_verts = data.all_verts;
849  data.owned_elems = data.primary_elems;
850  }
851 
852 #else
853 
854  data.local_verts = data.all_verts;
855  data.owned_elems = data.primary_elems;
856 
857 #endif
858 
859  // Get the references for some standard internal tags such as material blocks, BCs, etc
861  data.mat_sets, Interface::UNION );MB_CHK_ERR( rval );
862 
864  data.neu_sets, Interface::UNION );MB_CHK_ERR( rval );
865 
867  data.diri_sets, Interface::UNION );MB_CHK_ERR( rval );
868 
869  return moab::MB_SUCCESS;
870 }

References appData::all_verts, GlobalContext::appDatas, moab::Range::clear(), moab::ParallelComm::comm(), context, appData::dimension, appData::diri_sets, GlobalContext::dirichlet_tag, moab::Range::empty(), ErrorCode, appData::file_set, moab::ParallelComm::filter_pstatus(), moab::Interface::get_entities_by_dimension(), moab::Interface::get_entities_by_type(), moab::Interface::get_entities_by_type_and_tag(), appData::ghost_elems, appData::ghost_vertices, iMOAB_SetGlobalInfo(), appData::local_verts, appData::mat_sets, GlobalContext::material_tag, MB_CHK_ERR, MB_SUCCESS, MBENTITYSET, GlobalContext::MBI, MBVERTEX, GlobalContext::MPI_initialized, appData::neu_sets, GlobalContext::neumann_tag, appData::owned_elems, appData::owned_verts, appData::point_cloud, appData::primary_elems, PSTATUS_GHOST, PSTATUS_NOT, PSTATUS_NOT_OWNED, moab::Range::size(), moab::subtract(), and moab::Interface::UNION.

Referenced by iMOAB_LoadMesh().

◆ iMOAB_WriteLocalMesh()

ErrCode iMOAB_WriteLocalMesh ( iMOAB_AppID  pid,
iMOAB_String  prefix 
)

Write a local MOAB mesh copy.

Note
The interface will write one file per task. Only the local mesh set and solution data associated to the application will be written to the file. Very convenient method used for debugging.

Operations: Not Collective (independent operation).

Parameters
[in]pid(iMOAB_AppID) The unique pointer to the application ID.
[in]prefix(iMOAB_String) The MOAB file prefix that will be used with task ID in parallel.
Returns
ErrCode The error code indicating success or failure.

Definition at line 754 of file iMOAB.cpp.

755 {
756  IMOAB_CHECKPOINTER( prefix, 2 );
757  IMOAB_ASSERT( strlen( prefix ), "Invalid prefix string length." );
758 
759  std::ostringstream file_name;
760  int rank = 0, size = 1;
761 #ifdef MOAB_HAVE_MPI
762  rank = context.pcomms[*pid]->rank();
763  size = context.pcomms[*pid]->size();
764 #endif
765  file_name << prefix << "_" << size << "_" << rank << ".h5m";
766  // Now let us actually write the file to disk with appropriate options
767  ErrorCode rval = context.MBI->write_file( file_name.str().c_str(), 0, 0, &context.appDatas[*pid].file_set, 1 );MB_CHK_ERR( rval );
768 
769  return moab::MB_SUCCESS;
770 }

References GlobalContext::appDatas, context, ErrorCode, IMOAB_ASSERT, IMOAB_CHECKPOINTER, MB_CHK_ERR, MB_SUCCESS, GlobalContext::MBI, size, and moab::Interface::write_file().

◆ iMOAB_WriteMesh()

ErrCode iMOAB_WriteMesh ( iMOAB_AppID  pid,
const iMOAB_String  filename,
const iMOAB_String  write_options 
)

Write a MOAB mesh along with the solution tags to a file.

Note
The interface will write one single file (H5M) and for serial files (VTK/Exodus), it will write one file per task. Write options include parallel write options, if needed. Only the mesh set and solution data associated to the application will be written to the file.

Operations: Collective for parallel write, non collective for serial write.

Parameters
[in]pid(iMOAB_AppID) The unique pointer to the application ID.
[in]filename(iMOAB_String) The MOAB mesh file (H5M) to write all the entities contained in the internal application mesh set.
[in]write_options(iMOAB_String) Additional options for writing the MOAB mesh in parallel.
Returns
ErrCode The error code indicating success or failure.

Definition at line 693 of file iMOAB.cpp.

694 {
695  IMOAB_CHECKPOINTER( filename, 2 );
696  IMOAB_ASSERT( strlen( filename ), "Invalid filename length." );
697 
698  appData& data = context.appDatas[*pid];
699  EntityHandle fileSet = data.file_set;
700 
701  std::ostringstream newopts;
702 #ifdef MOAB_HAVE_MPI
703  std::string write_opts( ( write_options ? write_options : "" ) );
704  std::string pcid( "PARALLEL_COMM=" );
705  std::size_t found = write_opts.find( pcid );
706 
707  if( found != std::string::npos )
708  {
709  std::cerr << " cannot specify PARALLEL_COMM option, it is implicit \n";
710  return moab::MB_FAILURE;
711  }
712 
713  // if write in parallel, add pc option, to be sure about which ParallelComm instance is used
714  std::string pw( "PARALLEL=WRITE_PART" );
715  found = write_opts.find( pw );
716 
717  if( found != std::string::npos )
718  {
719  newopts << "PARALLEL_COMM=" << *pid << ";";
720  }
721 
722 #endif
723 
724  if( write_options ) newopts << write_options;
725 
726  std::vector< Tag > copyTagList = data.tagList;
727  // append Global ID and Parallel Partition
728  std::string gid_name_tag( "GLOBAL_ID" );
729 
730  // export global id tag, we need it always
731  if( data.tagMap.find( gid_name_tag ) == data.tagMap.end() )
732  {
733  Tag gid = context.MBI->globalId_tag();
734  copyTagList.push_back( gid );
735  }
736  // also Parallel_Partition PARALLEL_PARTITION
737  std::string pp_name_tag( "PARALLEL_PARTITION" );
738 
739  // write parallel part tag too, if it exists
740  if( data.tagMap.find( pp_name_tag ) == data.tagMap.end() )
741  {
742  Tag ptag;
743  context.MBI->tag_get_handle( pp_name_tag.c_str(), ptag );
744  if( ptag ) copyTagList.push_back( ptag );
745  }
746 
747  // Now let us actually write the file to disk with appropriate options
748  ErrorCode rval = context.MBI->write_file( filename, 0, newopts.str().c_str(), &fileSet, 1, &copyTagList[0],
749  (int)copyTagList.size() );MB_CHK_ERR( rval );
750 
751  return moab::MB_SUCCESS;
752 }

References GlobalContext::appDatas, context, ErrorCode, appData::file_set, moab::Interface::globalId_tag(), IMOAB_ASSERT, IMOAB_CHECKPOINTER, MB_CHK_ERR, MB_SUCCESS, GlobalContext::MBI, moab::Interface::tag_get_handle(), appData::tagList, appData::tagMap, and moab::Interface::write_file().

◆ split_tag_names()

static void split_tag_names ( std::string  input_names,
std::string &  separator,
std::vector< std::string > &  list_tag_names 
)
static

Definition at line 440 of file iMOAB.cpp.

443 {
444  size_t pos = 0;
445  std::string token;
446  while( ( pos = input_names.find( separator ) ) != std::string::npos )
447  {
448  token = input_names.substr( 0, pos );
449  if( !token.empty() ) list_tag_names.push_back( token );
450  // std::cout << token << std::endl;
451  input_names.erase( 0, pos + separator.length() );
452  }
453  if( !input_names.empty() )
454  {
455  // if leftover something, or if not ended with delimiter
456  list_tag_names.push_back( input_names );
457  }
458  return;
459 }

Referenced by iMOAB_DefineTagStorage(), iMOAB_GetDoubleTagStorage(), iMOAB_SetDoubleTagStorage(), and iMOAB_SetDoubleTagStorageWithGid().

Variable Documentation

◆ context