Commit a6a22b19 authored by Iulian Grindeanu's avatar Iulian Grindeanu Committed by vijaysm
Browse files

- Add element ownership, connectivity, global IDs

- Add dirichlet set to the test file
- Also, report number of surface BCs, not neumann sets
   report number of vertices in diri sets, not number of diri sets
   (look at signature of GetPointerToSurfaceBC, etc)
parent 5a712939
......@@ -35,6 +35,7 @@ struct appData {
EntityHandle file_set;
Range all_verts;
Range primary_elems;
int dimension; // 2 or 3, dimension of primary elements (redundant?)
Range mat_sets;
std::map<int, int> matIndex; // map from global block id to index in mat_sets
Range neu_sets;
......@@ -448,31 +449,72 @@ ErrCode GetMeshInfo( iMOAB_AppID pid, int* num_visible_vertices, int* num_visibl
// this will include ghost elements
// we should keep a data structure with mesh, sets, etc, for each pid
//
EntityHandle fileSet=appDatas[*pid].file_set;
appData & data = appDatas[*pid];
EntityHandle fileSet=data.file_set;
ErrorCode rval = MBI->get_entities_by_type(fileSet, MBVERTEX, appDatas[*pid].all_verts, true); // recursive
if (MB_SUCCESS!=rval)
return 1;
*num_visible_vertices = (int) appDatas[*pid].all_verts.size();
*num_visible_vertices = (int) data.all_verts.size();
// is dimension 3?
rval = MBI->get_entities_by_dimension(fileSet, 3, appDatas[*pid].primary_elems, true); // recursive
rval = MBI->get_entities_by_dimension(fileSet, 3, data.primary_elems, true); // recursive
if (MB_SUCCESS!=rval)
return 1;
*num_visible_elements = (int) appDatas[*pid].primary_elems.size();
data.dimension = 3;
if (data.primary_elems.empty())
{
appDatas[*pid].dimension = 2;
rval = MBI->get_entities_by_dimension(fileSet, 2, data.primary_elems, true); // recursive
if (MB_SUCCESS!=rval)
return 1;
if (data.primary_elems.empty())
return 1; // no elements of dimension 2 or 3
}
*num_visible_elements = (int) data.primary_elems.size();
// get all blocks, BCs, etc
rval = MBI->get_entities_by_type_and_tag(fileSet, MBENTITYSET, &(gtags[0]), 0, 1, appDatas[*pid].mat_sets , Interface::UNION);
rval = MBI->get_entities_by_type_and_tag(fileSet, MBENTITYSET, &(gtags[0]), 0, 1, data.mat_sets , Interface::UNION);
if (MB_SUCCESS!=rval)
return 1;
*num_visible_blocks = (int)appDatas[*pid].mat_sets.size();
rval = MBI->get_entities_by_type_and_tag(fileSet, MBENTITYSET, &(gtags[1]), 0, 1, appDatas[*pid].neu_sets , Interface::UNION);
*num_visible_blocks = data.mat_sets.size();
rval = MBI->get_entities_by_type_and_tag(fileSet, MBENTITYSET, &(gtags[1]), 0, 1, data.neu_sets , Interface::UNION);
if (MB_SUCCESS!=rval)
return 1;
*num_visible_surfaceBC = (int)appDatas[*pid].neu_sets.size();
rval = MBI->get_entities_by_type_and_tag(fileSet, MBENTITYSET, &(gtags[2]), 0, 1, appDatas[*pid].diri_sets , Interface::UNION);
*num_visible_surfaceBC = 0;
// count how many faces are in each neu set, and how many regions are
// adjacent to them;
int numNeuSets = (int)data.neu_sets.size();
for (int i=0; i<numNeuSets; i++)
{
Range subents;
EntityHandle nset = data.neu_sets[i];
rval = MBI->get_entities_by_dimension(nset, data.dimension-1, subents);
if (MB_SUCCESS!=rval)
return 1;
for (Range::iterator it=subents.begin(); it!=subents.end(); ++it)
{
EntityHandle subent = *it;
Range adjPrimaryEnts;
rval = MBI->get_adjacencies(&subent, 1, data.dimension, false, adjPrimaryEnts);
if (MB_SUCCESS!=rval)
return 1;
*num_visible_surfaceBC += (int)adjPrimaryEnts.size();
}
}
rval = MBI->get_entities_by_type_and_tag(fileSet, MBENTITYSET, &(gtags[2]), 0, 1, data.diri_sets , Interface::UNION);
if (MB_SUCCESS!=rval)
return 1;
*num_visible_vertexBC= (int)appDatas[*pid].diri_sets.size();
*num_visible_vertexBC= 0;
int numDiriSets = (int)data.diri_sets.size();
for (int i=0; i<numDiriSets; i++)
{
Range verts;
EntityHandle diset = data.diri_sets[i];
rval = MBI->get_entities_by_dimension(diset, 0, verts);
if (MB_SUCCESS!=rval)
return 1;
*num_visible_vertexBC += (int)verts.size();
}
return 0;
}
......@@ -631,7 +673,7 @@ ErrCode GetBlockInfo(iMOAB_AppID pid, iMOAB_GlobalID global_block_ID,
return 0;
}
#if 0
/**
\fn ErrorCode GetElementConnectivity(iMOAB_AppID pid, iMOAB_GlobalID global_block_ID, int connectivity_length, int* element_connectivity)
\brief Get the connectivity for elements within a certain block, ordered based on global element IDs
......@@ -643,7 +685,35 @@ ErrCode GetBlockInfo(iMOAB_AppID pid, iMOAB_GlobalID global_block_ID,
\param[in] connectivity_length (int) The allocated size of array (typical <TT>size := vertices_per_element*num_visible_elements</TT>)
\param[out] element_connectivity (int*) The connectivity array to store element ordering in MOAB canonical numbering scheme (array allocated by client); array contains vertex identifiers with global ID numbering
*/
ErrorCode GetElementConnectivity(iMOAB_AppID pid, iMOAB_GlobalID global_block_ID, int connectivity_length, int* element_connectivity);
ErrCode GetElementConnectivity(iMOAB_AppID pid, iMOAB_GlobalID global_block_ID, int connectivity_length, int* element_connectivity)
{
std::map<int, int> & matMap = appDatas[*pid].matIndex;
std::map<int,int>::iterator it = matMap.find(global_block_ID);
if (it==matMap.end())
return 1; // error in finding block with id
int blockIndex = matMap[global_block_ID];
EntityHandle matMeshSet = appDatas[*pid].mat_sets[blockIndex];
std::vector<EntityHandle> elems;
ErrorCode rval = MBI-> get_entities_by_handle(matMeshSet, elems);
if (MB_SUCCESS!=rval || elems.empty() )
return 1;
std::vector<EntityHandle> vconnect;
rval = MBI->get_connectivity(&elems[0], elems.size(), vconnect);
if (MB_SUCCESS!=rval)
return 1;
if (connectivity_length!=(int)vconnect.size())
return 1; // mismatched sizes
//gtags[3] is global id tag
rval = MBI->tag_get_data(gtags[3], &vconnect[0], connectivity_length, element_connectivity);
if (MB_SUCCESS!=rval)
return 1;
return 0;
}
/**
\fn ErrorCode GetElementOwnership(iMOAB_AppID pid, iMOAB_GlobalID global_block_ID, int num_elements_in_block, int* element_ownership)
......@@ -656,7 +726,33 @@ ErrorCode GetElementConnectivity(iMOAB_AppID pid, iMOAB_GlobalID global_block_ID
\param[in] num_elements_in_block (int) The allocated size of ownership array, same as <TT>num_elements_in_block</TT> returned from GetBlockInfo()
\param[out] element_ownership (int*) The ownership array to store processor ID for all elements (array allocated by client)
*/
ErrorCode GetElementOwnership(iMOAB_AppID pid, iMOAB_GlobalID global_block_ID, int num_elements_in_block, int* element_ownership);
ErrCode GetElementOwnership(iMOAB_AppID pid, iMOAB_GlobalID global_block_ID, int num_elements_in_block, int* element_ownership)
{
std::map<int, int> & matMap = appDatas[*pid].matIndex;
ParallelComm * pco = pcomms[*pid];
std::map<int,int>::iterator it = matMap.find(global_block_ID);
if (it==matMap.end())
return 1; // error in finding block with id
int blockIndex = matMap[global_block_ID];
EntityHandle matMeshSet = appDatas[*pid].mat_sets[blockIndex];
Range elems;
ErrorCode rval = MBI-> get_entities_by_handle(matMeshSet, elems);
if (MB_SUCCESS!=rval || elems.empty() )
return 1;
if (num_elements_in_block!=(int)elems.size())
return 1; // bad memory allocation
int i=0;
for (Range::iterator vit=elems.begin(); vit!=elems.end(); vit++, i++)
{
rval = pco-> get_owner(*vit, element_ownership[i]);
if (MB_SUCCESS!=rval)
return 1;
}
return 0;
}
/**
\fn ErrorCode GetElementID(iMOAB_AppID pid, iMOAB_GlobalID global_block_ID, int num_elements_in_block, iMOAB_GlobalID* global_element_ID, iMOAB_LocalID* local_element_ID)
......@@ -670,8 +766,32 @@ ErrorCode GetElementOwnership(iMOAB_AppID pid, iMOAB_GlobalID global_block_ID, i
\param[out] global_element_ID (iMOAB_GlobalID*) The global IDs for all locally visible elements (array allocated by client)
\param[out] local_element_ID (iMOAB_LocalID*) (<I><TT>Optional</TT></I>) The local IDs for all locally visible elements (array allocated by client)
*/
ErrorCode GetElementID(iMOAB_AppID pid, iMOAB_GlobalID global_block_ID, int num_elements_in_block, iMOAB_GlobalID* global_element_ID, iMOAB_LocalID* local_element_ID);
ErrCode GetElementID(iMOAB_AppID pid, iMOAB_GlobalID global_block_ID, int num_elements_in_block, iMOAB_GlobalID* global_element_ID, iMOAB_LocalID* local_element_ID)
{
std::map<int, int> & matMap = appDatas[*pid].matIndex;
std::map<int,int>::iterator it = matMap.find(global_block_ID);
if (it==matMap.end())
return 1; // error in finding block with id
int blockIndex = matMap[global_block_ID];
EntityHandle matMeshSet = appDatas[*pid].mat_sets[blockIndex];
Range elems;
ErrorCode rval = MBI-> get_entities_by_handle(matMeshSet, elems);
if (MB_SUCCESS!=rval || elems.empty() )
return 1;
if (num_elements_in_block!=(int)elems.size())
return 1; // bad memory allocation
rval = MBI->tag_get_data(gtags[3], elems, global_element_ID);
if (MB_SUCCESS!=rval )
return 1;
for (int i=0; i<num_elements_in_block; i++)
local_element_ID[i]=i;
return 0;
}
#if 0
/**
\fn ErrorCode GetPointerToSurfaceBC(iMOAB_AppID pid, int surface_BC_length, iMOAB_GlobalID* global_element_ID, int* reference_surface_ID, int* boundary_condition_value)
\brief Get the surface boundary condition information
......
......@@ -53,8 +53,8 @@ int main(int argc, char * argv[])
" %d visible vertices\n"
" %d visible elements\n"
" %d visible blocks\n"
" %d visible neumann sets\n"
" %d visible dirichlet sets\n", rank, nverts, nelem, nblocks, nsbc, ndbc);
" %d visible neumann BCs\n"
" %d visible dirichlet BCs\n", rank, nverts, nelem, nblocks, nsbc, ndbc);
}
iMOAB_GlobalID * vGlobalID = (iMOAB_GlobalID*)malloc(nverts*sizeof(iMOAB_GlobalID)) ;
......@@ -95,6 +95,30 @@ int main(int argc, char * argv[])
rc = GetBlockInfo(pid, gbIDs[i] , &vertices_per_element, &num_elements_in_block);
CHECKRC(rc, "failed to elem block info");
printf(" has %4d elements with %d vertices per element\n", num_elements_in_block, vertices_per_element);
int size_conn= num_elements_in_block*vertices_per_element;
iMOAB_GlobalID * element_connectivity = (iMOAB_GlobalID*) malloc (sizeof(iMOAB_GlobalID)*size_conn);
rc = GetElementConnectivity(pid, gbIDs[i], size_conn, element_connectivity);
CHECKRC(rc, "failed to get block elem connectivity");
int * element_ownership = (int*) malloc (sizeof(int)*num_elements_in_block);
GetElementOwnership(pid, gbIDs[i], num_elements_in_block, element_ownership);
CHECKRC(rc, "failed to get block elem ownership");
iMOAB_GlobalID* global_element_ID = (iMOAB_GlobalID*)malloc(sizeof(iMOAB_GlobalID)*num_elements_in_block);
iMOAB_LocalID* local_element_ID =(iMOAB_LocalID*)malloc(sizeof(iMOAB_LocalID)*num_elements_in_block);
rc = GetElementID(pid, gbIDs[i], num_elements_in_block, global_element_ID, local_element_ID);
CHECKRC(rc, "failed to get block elem IDs");
for (int j=0; j< num_elements_in_block; j++)
{
printf(" elem %3d owned by %d gid: %4d -- ", j, element_ownership[j], global_element_ID[j]);
for (int k=0; k<vertices_per_element; k++)
printf( " %5d", element_connectivity[j*vertices_per_element+k]);
printf("\n");
}
free(global_element_ID);
free(local_element_ID);
free (element_connectivity);
free (element_ownership);
}
}
......
No preview for this file type
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment