45 #include <PanzerAdaptersSTK_config.hpp>
46 #include <stk_mesh/base/FEMHelpers.hpp>
51 namespace panzer_stk {
66 PANZER_FUNC_TIME_MONITOR(
"panzer::CubeHexMeshFactory::buildMesh()");
72 mesh->initialize(parallelMach);
82 PANZER_FUNC_TIME_MONITOR(
"panzer::CubeHexMeshFactory::buildUncomittedMesh()");
86 machRank_ = stk::parallel_machine_rank(parallelMach);
87 machSize_ = stk::parallel_machine_size(parallelMach);
99 const int maxFactor = 50;
102 int factors[maxFactor];
103 for (
int jj = 0; jj < maxFactor; jj++) factors[jj] = 0;
104 for (
int jj = 2; jj < maxFactor; jj++) {
107 int temp = ProcTemp/jj;
108 if (temp*jj == ProcTemp) {
118 for (
int jj = maxFactor-1; jj > 0; jj--) {
119 while (factors[jj] != 0) {
135 "Cannot build CubeHexMeshFactory, the product of \"X Procs\", \"Y Procs\", and \"Z Procs\""
136 " must equal the number of processors.");
150 PANZER_FUNC_TIME_MONITOR(
"panzer::CubeHexMeshFactory::completeMeshConstruction()");
166 out <<
"CubeHexMesh: Building sub cells" << std::endl;
171 out <<
"CubeHexMesh: NOT building sub cells" << std::endl;
207 x0_ = paramList->
get<
double>(
"X0");
208 y0_ = paramList->
get<
double>(
"Y0");
209 z0_ = paramList->
get<
double>(
"Z0");
211 xf_ = paramList->
get<
double>(
"Xf");
212 yf_ = paramList->
get<
double>(
"Yf");
213 zf_ = paramList->
get<
double>(
"Zf");
233 if (not buildSubcells_ && createEdgeBlocks_) {
238 out <<
"CubeHexMesh: NOT creating edge blocks because building sub cells disabled" << std::endl;
239 createEdgeBlocks_ =
false;
241 if (not buildSubcells_ && createFaceBlocks_) {
246 out <<
"CubeHexMesh: NOT creating face blocks because building sub cells disabled" << std::endl;
247 createFaceBlocks_ =
false;
260 if(defaultParams == Teuchos::null) {
263 defaultParams->
set<
double>(
"X0",0.0);
264 defaultParams->
set<
double>(
"Y0",0.0);
265 defaultParams->
set<
double>(
"Z0",0.0);
267 defaultParams->
set<
double>(
"Xf",1.0);
268 defaultParams->
set<
double>(
"Yf",1.0);
269 defaultParams->
set<
double>(
"Zf",1.0);
271 defaultParams->
set<
int>(
"X Blocks",1);
272 defaultParams->
set<
int>(
"Y Blocks",1);
273 defaultParams->
set<
int>(
"Z Blocks",1);
275 defaultParams->
set<
int>(
"X Procs",-1);
276 defaultParams->
set<
int>(
"Y Procs",1);
277 defaultParams->
set<
int>(
"Z Procs",1);
279 defaultParams->
set<
int>(
"X Elements",5);
280 defaultParams->
set<
int>(
"Y Elements",5);
281 defaultParams->
set<
int>(
"Z Elements",5);
283 defaultParams->
set<
bool>(
"Build Interface Sidesets",
false);
285 defaultParams->
set<
bool>(
"Build Subcells",
true);
288 defaultParams->
set<
bool>(
"Create Edge Blocks",
false,
"Create edge blocks in the mesh");
289 defaultParams->
set<
bool>(
"Create Face Blocks",
false,
"Create face blocks in the mesh");
292 bcs.
set<
int>(
"Count",0);
295 return defaultParams;
317 typedef shards::Hexahedron<8> HexTopo;
318 const CellTopologyData * ctd = shards::getCellTopologyData<HexTopo>();
319 const CellTopologyData * side_ctd = shards::CellTopology(ctd).getBaseCellTopologyData(2,0);
320 const CellTopologyData * edge_ctd = shards::CellTopology(ctd).getBaseCellTopologyData(1,0);
321 const CellTopologyData * face_ctd = shards::CellTopology(ctd).getBaseCellTopologyData(2,0);
329 std::stringstream ebPostfix;
330 ebPostfix <<
"-" << bx <<
"_" << by <<
"_" << bz;
359 std::stringstream ss;
360 ss <<
"vertical_" << bx-1;
364 std::stringstream ss;
365 ss <<
"horizontal_" << by-1;
369 std::stringstream ss;
370 ss <<
"transverse_" << bz-1;
383 for(
int xBlock=0;xBlock<
xBlocks_;xBlock++) {
384 for(
int yBlock=0;yBlock<
yBlocks_;yBlock++) {
385 for(
int zBlock=0;zBlock<
zBlocks_;zBlock++) {
386 buildBlock(parallelMach,xBlock,yBlock,zBlock,mesh);
400 panzer::GlobalOrdinal myXElems_start = sizeAndStartX.first;
401 panzer::GlobalOrdinal myXElems_end = myXElems_start+sizeAndStartX.second;
402 panzer::GlobalOrdinal myYElems_start = sizeAndStartY.first;
403 panzer::GlobalOrdinal myYElems_end = myYElems_start+sizeAndStartY.second;
404 panzer::GlobalOrdinal myZElems_start = sizeAndStartZ.first;
405 panzer::GlobalOrdinal myZElems_end = myZElems_start+sizeAndStartZ.second;
411 double deltaX = (
xf_-
x0_)/
double(totalXElems);
412 double deltaY = (
yf_-
y0_)/
double(totalYElems);
413 double deltaZ = (
zf_-
z0_)/
double(totalZElems);
415 std::vector<double> coord(3,0.0);
418 for(panzer::GlobalOrdinal nx=myXElems_start;nx<myXElems_end+1;++nx) {
420 for(panzer::GlobalOrdinal ny=myYElems_start;ny<myYElems_end+1;++ny) {
422 for(panzer::GlobalOrdinal nz=myZElems_start;nz<myZElems_end+1;++nz) {
425 mesh.
addNode(nz*(totalYElems+1)*(totalXElems+1)+ny*(totalXElems+1)+nx+1,coord);
430 std::stringstream blockName;
431 blockName <<
"eblock-" << xBlock <<
"_" << yBlock <<
"_" << zBlock;
435 for(panzer::GlobalOrdinal nx=myXElems_start;nx<myXElems_end;++nx) {
436 for(panzer::GlobalOrdinal ny=myYElems_start;ny<myYElems_end;++ny) {
437 for(panzer::GlobalOrdinal nz=myZElems_start;nz<myZElems_end;++nz) {
438 stk::mesh::EntityId gid = totalXElems*totalYElems*nz+totalXElems*ny+nx+1;
439 std::vector<stk::mesh::EntityId> nodes(8);
440 nodes[0] = nx+1+ny*(totalXElems+1) +nz*(totalYElems+1)*(totalXElems+1);
441 nodes[1] = nodes[0]+1;
442 nodes[2] = nodes[1]+(totalXElems+1);
443 nodes[3] = nodes[2]-1;
444 nodes[4] = nodes[0]+(totalYElems+1)*(totalXElems+1);
445 nodes[5] = nodes[1]+(totalYElems+1)*(totalXElems+1);
446 nodes[6] = nodes[2]+(totalYElems+1)*(totalXElems+1);
447 nodes[7] = nodes[3]+(totalYElems+1)*(totalXElems+1);
459 panzer::GlobalOrdinal minElements =
nXElems_/size;
460 panzer::GlobalOrdinal extra =
nXElems_ - minElements*size;
466 panzer::GlobalOrdinal nume=0, start=0;
467 if(panzer::GlobalOrdinal(xProcLoc)<extra) {
468 nume = minElements+1;
469 start = xProcLoc*(minElements+1);
473 start = extra*(minElements+1)+(xProcLoc-extra)*minElements;
476 return std::make_pair(start+
nXElems_*xBlock,nume);
485 panzer::GlobalOrdinal minElements =
nYElems_/size;
486 panzer::GlobalOrdinal extra =
nYElems_ - minElements*size;
492 panzer::GlobalOrdinal nume=0, start=0;
493 if(panzer::GlobalOrdinal(yProcLoc)<extra) {
494 nume = minElements+1;
495 start = yProcLoc*(minElements+1);
499 start = extra*(minElements+1)+(yProcLoc-extra)*minElements;
502 return std::make_pair(start+
nYElems_*yBlock,nume);
510 panzer::GlobalOrdinal minElements =
nZElems_/size;
511 panzer::GlobalOrdinal extra =
nZElems_ - minElements*size;
517 panzer::GlobalOrdinal nume=0, start=0;
518 if(zProcLoc<Teuchos::as<std::size_t>(extra)) {
519 nume = minElements+1;
520 start = zProcLoc*(minElements+1);
524 start = extra*(minElements+1)+(zProcLoc-extra)*minElements;
527 return std::make_pair(start+
nZElems_*zBlock,nume);
539 std::vector<stk::mesh::Entity> localElmts;
542 stk::mesh::EntityId offset[6];
544 offset[1] = offset[0] + totalXElems*totalZElems;
545 offset[2] = offset[1] + totalYElems*totalZElems;
546 offset[3] = offset[2] + totalXElems*totalZElems;
547 offset[4] = offset[3] + totalYElems*totalZElems;
548 offset[5] = offset[4] + totalXElems*totalYElems;
553 std::vector<stk::mesh::Entity>::const_iterator itr;
554 for(itr=localElmts.begin();itr!=localElmts.end();++itr) {
555 stk::mesh::Entity element = (*itr);
558 std::size_t nx,ny,nz;
559 nz = (gid-1) / (totalXElems*totalYElems);
560 gid = (gid-1)-nz*(totalXElems*totalYElems);
561 ny = gid / totalXElems;
562 nx = gid-ny*totalXElems;
564 std::vector<stk::mesh::Part*> parts;
568 mesh.
getBulkData()->declare_element_side(element, 4, parts);
570 if(nz+1==totalZElems) {
572 mesh.
getBulkData()->declare_element_side(element, 5, parts);
577 mesh.
getBulkData()->declare_element_side(element, 0, parts);
579 if(ny+1==totalYElems) {
581 mesh.
getBulkData()->declare_element_side(element, 2, parts);
586 mesh.
getBulkData()->declare_element_side(element, 3, parts);
588 if(nx+1==totalXElems) {
590 mesh.
getBulkData()->declare_element_side(element, 1, parts);
601 const stk::mesh::EntityRank side_rank = mesh.
getSideRank();
608 stk::mesh::Part * left = mesh.
getSideset(
"left");
609 stk::mesh::Part * right = mesh.
getSideset(
"right");
610 stk::mesh::Part * top = mesh.
getSideset(
"top");
611 stk::mesh::Part * bottom = mesh.
getSideset(
"bottom");
612 stk::mesh::Part * front = mesh.
getSideset(
"front");
613 stk::mesh::Part * back = mesh.
getSideset(
"back");
615 std::vector<stk::mesh::Part*> vertical;
616 std::vector<stk::mesh::Part*> horizontal;
617 std::vector<stk::mesh::Part*> transverse;
621 std::stringstream ss;
622 ss <<
"vertical_" << bx-1;
623 vertical.push_back(mesh.
getSideset(ss.str()));
626 std::stringstream ss;
627 ss <<
"horizontal_" << by-1;
628 horizontal.push_back(mesh.
getSideset(ss.str()));
631 std::stringstream ss;
632 ss <<
"transverse_" << bz-1;
633 transverse.push_back(mesh.
getSideset(ss.str()));
637 std::vector<stk::mesh::Entity> localElmts;
643 std::vector<stk::mesh::Entity>::const_iterator itr;
644 for(itr=localElmts.begin();itr!=localElmts.end();++itr) {
645 stk::mesh::Entity element = (*itr);
648 std::size_t nx,ny,nz;
649 nz = (gid-1) / (totalXElems*totalYElems);
650 gid = (gid-1)-nz*(totalXElems*totalYElems);
651 ny = gid / totalXElems;
652 nx = gid-ny*totalXElems;
674 if(nz+1==totalZElems) {
705 if(ny+1==totalYElems) {
736 if(nx+1==totalXElems) {
754 stk::mesh::Part * origin = mesh.
getNodeset(
"origin");
760 stk::mesh::Entity node = bulkData->get_entity(mesh.
getNodeRank(),1);
774 stk::mesh::Selector owned_block = metaData->locally_owned_part();
776 std::vector<stk::mesh::Entity> edges;
777 bulkData->get_entities(mesh.
getEdgeRank(), owned_block, edges);
790 stk::mesh::Selector owned_block = metaData->locally_owned_part();
792 std::vector<stk::mesh::Entity> faces;
793 bulkData->get_entities(mesh.
getFaceRank(), owned_block, faces);
800 std::size_t i=0,j=0,k=0;
806 return Teuchos::tuple(i,j,k);
void addSides(STK_Interface &mesh) const
void initializeWithDefaults()
virtual ~CubeHexMeshFactory()
Destructor.
void addNodeset(const std::string &name)
unsigned entityOwnerRank(stk::mesh::Entity entity) const
void addEntityToNodeset(stk::mesh::Entity entity, stk::mesh::Part *nodeset)
void addEdgeBlocks(STK_Interface &mesh) const
static void parsePeriodicBCList(const Teuchos::RCP< Teuchos::ParameterList > &pl, std::vector< Teuchos::RCP< const PeriodicBC_MatcherBase > > &periodicBC, bool &useBBoxSearch)
void addFaceBlocks(STK_Interface &mesh) const
void buildMetaData(stk::ParallelMachine parallelMach, STK_Interface &mesh) const
basic_FancyOStream & setShowProcRank(const bool showProcRank)
T & get(const std::string &name, T def_value)
stk::mesh::Part * getElementBlockPart(const std::string &name) const
get the block part
void buildElements(stk::ParallelMachine parallelMach, STK_Interface &mesh) const
panzer::GlobalOrdinal nXElems_
ParameterList & set(std::string const &name, T const &value, std::string const &docString="", RCP< const ParameterEntryValidator > const &validator=null)
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
void addEntityToSideset(stk::mesh::Entity entity, stk::mesh::Part *sideset)
void addSideSets(STK_Interface &mesh) const
stk::mesh::Part * getNodeset(const std::string &name) const
stk::mesh::Part * getFaceBlock(const std::string &name) const
get the block part
void addEdgeBlock(const std::string &elemBlockName, const std::string &edgeBlockName, const stk::topology &topology)
stk::mesh::Part * getSideset(const std::string &name) const
void addElement(const Teuchos::RCP< ElementDescriptor > &ed, stk::mesh::Part *block)
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
stk::mesh::EntityRank getSideRank() const
void initialize(stk::ParallelMachine parallelMach, bool setupIO=true, const bool buildRefinementSupport=false)
void rebalance(STK_Interface &mesh) const
Teuchos::RCP< stk::mesh::BulkData > getBulkData() const
void buildBlock(stk::ParallelMachine machRank, int xBlock, int yBlock, int zBlock, STK_Interface &mesh) const
void addNodeSets(STK_Interface &mesh) const
stk::mesh::EntityId elementGlobalId(std::size_t lid) const
void setMyParamList(const RCP< ParameterList > ¶mList)
void validateParametersAndSetDefaults(ParameterList const &validParamList, int const depth=1000)
static const std::string edgeBlockString
basic_FancyOStream & setOutputToRootOnly(const int rootRank)
void setParameterList(const Teuchos::RCP< Teuchos::ParameterList > ¶mList)
From ParameterListAcceptor.
std::pair< panzer::GlobalOrdinal, panzer::GlobalOrdinal > determineZElemSizeAndStart(int zBlock, unsigned int size, unsigned int rank) const
stk::mesh::EntityRank getFaceRank() const
void addNode(stk::mesh::EntityId gid, const std::vector< double > &coord)
void addEntitiesToEdgeBlock(std::vector< stk::mesh::Entity > entities, stk::mesh::Part *edgeblock)
void buildSubcells()
force the mesh to build subcells: edges and faces
Teuchos::Tuple< std::size_t, 3 > procTuple_
Teuchos::RCP< STK_Interface > buildMesh(stk::ParallelMachine parallelMach) const
Build the mesh object.
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const
From ParameterListAcceptor.
CubeHexMeshFactory()
Constructor.
Teuchos::RCP< stk::mesh::MetaData > getMetaData() const
bool isInitialized() const
Has initialize been called on this mesh object?
stk::mesh::Entity findConnectivityById(stk::mesh::Entity src, stk::mesh::EntityRank tgt_rank, unsigned rel_id) const
void addSideset(const std::string &name, const CellTopologyData *ctData)
static const std::string faceBlockString
double getMeshCoord(const int nx, const double deltaX, const double x0) const
std::pair< panzer::GlobalOrdinal, panzer::GlobalOrdinal > determineXElemSizeAndStart(int xBlock, unsigned int size, unsigned int rank) const
stk::mesh::Part * getEdgeBlock(const std::string &name) const
get the block part
std::vector< Teuchos::RCP< const PeriodicBC_MatcherBase > > periodicBCVec_
stk::mesh::EntityRank getEdgeRank() const
void addEntitiesToFaceBlock(std::vector< stk::mesh::Entity > entities, stk::mesh::Part *faceblock)
bool buildInterfaceSidesets_
void buildLocalElementIDs()
void addFaceBlock(const std::string &elemBlockName, const std::string &faceBlockName, const stk::topology &topology)
ParameterList & sublist(const std::string &name, bool mustAlreadyExist=false, const std::string &docString="")
std::pair< panzer::GlobalOrdinal, panzer::GlobalOrdinal > determineYElemSizeAndStart(int yBlock, unsigned int size, unsigned int rank) const
panzer::GlobalOrdinal nYElems_
void getMyElements(std::vector< stk::mesh::Entity > &elements) const
#define TEUCHOS_ASSERT(assertion_test)
stk::mesh::EntityRank getNodeRank() const
virtual Teuchos::RCP< STK_Interface > buildUncommitedMesh(stk::ParallelMachine parallelMach) const
panzer::GlobalOrdinal nZElems_
void addElementBlock(const std::string &name, const CellTopologyData *ctData)
std::string edgeBlockName_
std::string faceBlockName_
virtual void completeMeshConstruction(STK_Interface &mesh, stk::ParallelMachine parallelMach) const
Teuchos::Tuple< std::size_t, 3 > procRankToProcTuple(std::size_t procRank) const
what is the 3D tuple describe this processor distribution