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

Add routine for vertex bcs and surface bcs

- report from each task in driver
- minor interface changes
parent a6a22b19
......@@ -768,30 +768,39 @@ ErrCode GetElementOwnership(iMOAB_AppID pid, iMOAB_GlobalID global_block_ID, int
*/
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;
appData & data = appDatas[*pid];
std::map<int, int> & matMap = data.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];
EntityHandle matMeshSet = data.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;
// check that elems are among primary_elems in data
for (int i=0; i<num_elements_in_block; i++)
local_element_ID[i]=i;
{
local_element_ID[i]=data.primary_elems.index(elems[i]);
if (-1==local_element_ID[i])
return 1;// error, not in local primary elements
}
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
......@@ -804,7 +813,60 @@ ErrCode GetElementID(iMOAB_AppID pid, iMOAB_GlobalID global_block_ID, int num_el
\param[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)
\param[out] boundary_condition_value (int*) The boundary condition type as obtained from the mesh description (value of the NeumannSet defined on the element)
*/
ErrorCode GetPointerToSurfaceBC(iMOAB_AppID pid, int surface_BC_length, iMOAB_GlobalID* global_element_ID, int* reference_surface_ID, int* boundary_condition_value);
ErrCode GetPointerToSurfaceBC(iMOAB_AppID pid, int surface_BC_length, iMOAB_GlobalID* global_element_ID, int* reference_surface_ID, int* boundary_condition_value)
{
// we have to fill bc data for neumann sets;/
// it was counted above, in GetMeshInfo
appData & data = appDatas[*pid];
int numNeuSets = (int)data.neu_sets.size();
int index = 0; // index [0, surface_BC_length) for the arrays returned
for (int i=0; i<numNeuSets; i++)
{
Range subents;
EntityHandle nset = data.neu_sets[i];
ErrorCode rval = MBI->get_entities_by_dimension(nset, data.dimension-1, subents);
if (MB_SUCCESS!=rval)
return 1;
int neuVal ;
rval = MBI->tag_get_data(gtags[1], &nset, 1, &neuVal);
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;
// get global id of the primary ents, and side number of the quad/subentity
// this is moab ordering
for (Range::iterator pit=adjPrimaryEnts.begin(); pit!=adjPrimaryEnts.end(); pit++)
{
EntityHandle primaryEnt = *pit;
// get global id
int globalID;
rval = MBI->tag_get_data(gtags[3], &primaryEnt, 1, &globalID);
if (MB_SUCCESS!=rval)
return 1;
global_element_ID[index] = globalID;
int side_number, sense, offset;
rval = MBI->side_number(primaryEnt, subent, side_number, sense, offset);
if (MB_SUCCESS!=rval)
return 1;
reference_surface_ID[index] = side_number+1; // moab is from 0 to 5, it needs 1 to 6
boundary_condition_value[index] = neuVal;
index++;
}
}
}
if (index != surface_BC_length)
return 1; // error in array allocations
return 0;
}
/**
\fn ErrorCode GetPointerToVertexBC(iMOAB_AppID pid, int vertex_BC_length, iMOAB_GlobalID* global_vertext_ID, int* num_vertex_BC, int* boundary_condition_value)
......@@ -818,8 +880,44 @@ ErrorCode GetPointerToSurfaceBC(iMOAB_AppID pid, int surface_BC_length, iMOAB_Gl
\param[out] num_vertex_BC (int*) The allocated size of vertex boundary condition array, same as num_visible_vertexBC
\param[out] boundary_condition_value (int*) The boundary condition type as obtained from the mesh description (value of the DirichletSet defined on the vertex)
*/
ErrorCode GetPointerToVertexBC(iMOAB_AppID pid, int vertex_BC_length, iMOAB_GlobalID* global_vertext_ID, int* num_vertex_BC, int* boundary_condition_value);
ErrCode GetPointerToVertexBC(iMOAB_AppID pid, int vertex_BC_length,
iMOAB_GlobalID* global_vertext_ID, int* boundary_condition_value)
{
// it was counted above, in GetMeshInfo
appData & data = appDatas[*pid];
int numDiriSets = (int)data.diri_sets.size();
int index = 0; // index [0, vertex_BC_length) for the arrays returned
for (int i=0; i<numDiriSets; i++)
{
Range verts;
EntityHandle diset = data.diri_sets[i];
ErrorCode rval = MBI->get_entities_by_dimension(diset, 0, verts);
if (MB_SUCCESS!=rval)
return 1;
int diriVal;
rval = MBI->tag_get_data(gtags[2], &diset, 1, &diriVal);
if (MB_SUCCESS!=rval)
return 1;
for (Range::iterator vit=verts.begin(); vit!=verts.end(); ++vit)
{
EntityHandle vt =*vit;
int vgid;
rval = MBI->tag_get_data(gtags[3], &vt, 1, &vgid);
if (MB_SUCCESS!=rval)
return 1;
global_vertext_ID[index] = vgid;
boundary_condition_value[index] = diriVal;
index++;
}
}
if (vertex_BC_length!=index)
return 1; // array allocation issue
return 0;
}
#if 0
/**
\fn ErrorCode DefineTagStorage(iMOAB_AppID pid, iMOAB_String tag_storage_name, int* tag_type, int* components_per_entity, int tag_storage_name_length)
\brief Define a MOAB Tag corresponding to the application depending on requested types.
......
......@@ -305,10 +305,9 @@ ErrCode GetPointerToSurfaceBC(iMOAB_AppID pid, int surface_BC_length, iMOAB_Glob
\param[in] pid (iMOAB_AppID) The unique pointer to the application ID
\param[in] vertex_BC_length (int) The allocated size of vertex boundary condition array, same as <TT>num_visible_vertexBC</TT> returned by GetMeshInfo()
\param[out] global_vertext_ID (iMOAB_GlobalID*) The global vertex ID that has Dirichlet BC defined
\param[out] num_vertex_BC (int*) The allocated size of vertex boundary condition array, same as num_visible_vertexBC
\param[out] boundary_condition_value (int*) The boundary condition type as obtained from the mesh description (value of the DirichletSet defined on the vertex)
*/
ErrCode GetPointerToVertexBC(iMOAB_AppID pid, int vertex_BC_length, iMOAB_GlobalID* global_vertext_ID, int* num_vertex_BC, int* boundary_condition_value);
ErrCode GetPointerToVertexBC(iMOAB_AppID pid, int vertex_BC_length, iMOAB_GlobalID* global_vertext_ID, int* boundary_condition_value);
/**
\fn ErrCode DefineTagStorage(iMOAB_AppID pid, iMOAB_String tag_storage_name, int* tag_type, int* components_per_entity, int tag_storage_name_length)
......
......@@ -47,15 +47,6 @@ int main(int argc, char * argv[])
int nverts, nelem, nblocks, nsbc, ndbc;
rc = GetMeshInfo( pid, &nverts, &nelem, &nblocks, &nsbc, &ndbc);
CHECKRC(rc, "failed to get mesh info");
if (1==rank)
{
printf("on rank %d, there are \n"
" %d visible vertices\n"
" %d visible elements\n"
" %d visible blocks\n"
" %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)) ;
iMOAB_LocalID * vLocalID = (iMOAB_LocalID*)malloc(nverts*sizeof(iMOAB_GlobalID)) ;
......@@ -70,60 +61,90 @@ int main(int argc, char * argv[])
rc = GetVisibleVerticesCoordinates( pid, 3*nverts, coords);
CHECKRC(rc, "failed to get coordinates");
if (1==rank)
{
// printf some of the vertex id infos
int numToPrint = nverts;
printf("on rank %d some vertex info:\n", rank);
for (int i=0; i<numToPrint; i++)
printf(" vertex local id: %3d, rank ID:%d global ID: %3d coords: %g, %g, %g\n",vLocalID[i], vranks[i], vGlobalID[i],
coords[3*i], coords[3*i+1], coords[3*i+2]);
}
iMOAB_GlobalID * gbIDs = (iMOAB_GlobalID*) malloc(nblocks*sizeof(iMOAB_GlobalID));
iMOAB_LocalID * lbIDs = (iMOAB_LocalID*) malloc(nblocks*sizeof(iMOAB_LocalID));
rc = GetBlockID(pid, nblocks, gbIDs, lbIDs);
CHECKRC(rc, "failed to get block info");
if (1==rank)
for (int irank=0; irank<nprocs; irank++)
{
// printf some of the block info
printf("on rank %d some block info:\n", rank);
for (int i=0; i<nblocks; i++)
if (irank==rank)
{
printf(" block index: %3d, block ID: %3d \n", lbIDs[i], gbIDs[i] );
int vertices_per_element, num_elements_in_block;
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 some of the block info
printf("on rank %d, there are \n"
" %d visible vertices\n"
" %d visible elements\n"
" %d visible blocks\n"
" %d visible neumann BCs\n"
" %d visible dirichlet BCs\n", rank, nverts, nelem, nblocks, nsbc, ndbc);
// printf some of the vertex id infos
int numToPrint = nverts;
printf("on rank %d vertex info:\n", rank);
for (int i=0; i<numToPrint; i++)
printf(" vertex local id: %3d, rank ID:%d global ID: %3d coords: %g, %g, %g\n",vLocalID[i], vranks[i], vGlobalID[i],
coords[3*i], coords[3*i+1], coords[3*i+2]);
for (int i=0; i<nblocks; i++)
{
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");
printf(" block index: %3d, block ID: %3d \n", lbIDs[i], gbIDs[i] );
int vertices_per_element, num_elements_in_block;
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);
}
free(global_element_ID);
free(local_element_ID);
free (element_connectivity);
free (element_ownership);
// query surface BCs
iMOAB_GlobalID * surfBC_ID = (iMOAB_GlobalID*) malloc (sizeof(iMOAB_GlobalID)*nsbc);
int * ref_surf = (int *) malloc (sizeof(int)*nsbc);
int * bc_value = (int *) malloc (sizeof(int)*nsbc);
rc = GetPointerToSurfaceBC(pid, nsbc, surfBC_ID, ref_surf, bc_value);
CHECKRC(rc, "failed to get surf boundary conditions");
printf(" Surface boundary conditions:\n");
for (int i=0; i<nsbc; i++)
{
printf(" elemID %4d side:%d BC:%2d\n",surfBC_ID[i] ,ref_surf[i], bc_value[i] );
}
free(surfBC_ID);
free(ref_surf);
free(bc_value);
// query vertex BCs
iMOAB_GlobalID * vertBC_ID = (iMOAB_GlobalID*) malloc (sizeof(iMOAB_GlobalID)*ndbc);
int * vertBC_value = (int *) malloc (sizeof(int)*ndbc);
rc = GetPointerToVertexBC(pid, ndbc, vertBC_ID, vertBC_value);
CHECKRC(rc, "failed to get vertex boundary conditions");
printf(" Vertex boundary conditions:\n");
for (int i=0; i<ndbc; i++)
{
printf(" vertex %4d BC:%2d\n",vertBC_ID[i], vertBC_value[i] );
}
free(vertBC_ID);
free(vertBC_value);
}
MPI_Barrier(comm); // to avoid printing problems, as we write all procs data
}
// free allocated data
free(coords);
free (vGlobalID);
......
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