45 #include <PanzerAdaptersSTK_config.hpp>
52 namespace panzer_stk {
67 PANZER_FUNC_TIME_MONITOR(
"panzer::SquareQuadMeshFactory::buildMesh()");
73 mesh->initialize(parallelMach);
83 PANZER_FUNC_TIME_MONITOR(
"panzer::SquareQuadMeshFactory::buildUncomittedMesh()");
87 machRank_ = stk::parallel_machine_rank(parallelMach);
88 machSize_ = stk::parallel_machine_size(parallelMach);
100 const int maxFactor = 100;
103 int factors[maxFactor];
104 for (
int jj = 0; jj < maxFactor; jj++) factors[jj] = 0;
105 for (
int jj = 2; jj < maxFactor; jj++) {
108 int temp = ProcTemp/jj;
109 if (temp*jj == ProcTemp) {
119 for (
int jj = maxFactor-1; jj > 0; jj--) {
120 while (factors[jj] != 0) {
134 "Cannot build SquareQuadMeshFactory. The product of 'X Procs * Y Procs = " <<
xProcs_ <<
"*" << yProcs_ <<
" = " <<
xProcs_*yProcs_
135 <<
"' must equal the number of processors = " <<
machSize_
136 <<
"\n\n\t==> Run the simulation with an appropriate number of processors, i.e. #procs = " <<
xProcs_*yProcs_ <<
".\n");
149 PANZER_FUNC_TIME_MONITOR(
"panzer::SquareQuadMeshFactory::completeMeshConstruction()");
158 #ifndef ENABLE_UNIFORM
164 #ifndef ENABLE_UNIFORM
182 x0_ = paramList->
get<
double>(
"X0");
183 y0_ = paramList->
get<
double>(
"Y0");
185 xf_ = paramList->
get<
double>(
"Xf");
186 yf_ = paramList->
get<
double>(
"Yf");
207 if(defaultParams == Teuchos::null) {
210 defaultParams->
set<
double>(
"X0",0.0);
211 defaultParams->
set<
double>(
"Y0",0.0);
213 defaultParams->
set<
double>(
"Xf",1.0);
214 defaultParams->
set<
double>(
"Yf",1.0);
216 defaultParams->
set<
int>(
"X Blocks",1);
217 defaultParams->
set<
int>(
"Y Blocks",1);
219 defaultParams->
set<
int>(
"X Procs",-1);
220 defaultParams->
set<
int>(
"Y Procs",1);
222 defaultParams->
set<
int>(
"X Elements",5);
223 defaultParams->
set<
int>(
"Y Elements",5);
226 bcs.
set<
int>(
"Count",0);
229 return defaultParams;
243 typedef shards::Quadrilateral<4> QuadTopo;
244 const CellTopologyData * ctd = shards::getCellTopologyData<QuadTopo>();
245 const CellTopologyData * side_ctd = shards::CellTopology(ctd).getBaseCellTopologyData(1,0);
254 std::stringstream ebPostfix;
255 ebPostfix <<
"-" << bx <<
"_" << by;
265 #ifndef ENABLE_UNIFORM
272 std::stringstream ss;
273 ss <<
"vertical_" << bx-1;
277 std::stringstream ss;
278 ss <<
"horizontal_" << by-1;
292 for(
int xBlock=0;xBlock<
xBlocks_;xBlock++) {
293 for(
int yBlock=0;yBlock<
yBlocks_;yBlock++) {
306 int myXElems_start = sizeAndStartX.first;
307 int myXElems_end = myXElems_start+sizeAndStartX.second;
308 int myYElems_start = sizeAndStartY.first;
309 int myYElems_end = myYElems_start+sizeAndStartY.second;
313 double deltaX = (
xf_-
x0_)/
double(totalXElems);
314 double deltaY = (
yf_-
y0_)/
double(totalYElems);
316 std::vector<double> coord(2,0.0);
319 for(
int nx=myXElems_start;nx<myXElems_end+1;++nx) {
320 coord[0] = double(nx)*deltaX+
x0_;
321 for(
int ny=myYElems_start;ny<myYElems_end+1;++ny) {
322 coord[1] = double(ny)*deltaY+
y0_;
324 mesh.
addNode(ny*(totalXElems+1)+nx+1,coord);
328 std::stringstream blockName;
329 blockName <<
"eblock-" << xBlock <<
"_" << yBlock;
333 for(
int nx=myXElems_start;nx<myXElems_end;++nx) {
334 for(
int ny=myYElems_start;ny<myYElems_end;++ny) {
335 stk::mesh::EntityId gid = totalXElems*ny+nx+1;
336 std::vector<stk::mesh::EntityId> nodes(4);
337 nodes[0] = nx+1+ny*(totalXElems+1);
338 nodes[1] = nodes[0]+1;
339 nodes[2] = nodes[1]+(totalXElems+1);
340 nodes[3] = nodes[2]-1;
351 unsigned int minElements =
nXElems_/size;
352 unsigned int extra =
nXElems_ - minElements*size;
360 nume = minElements+1;
361 start = xProcLoc*(minElements+1);
365 start = extra*(minElements+1)+(xProcLoc-extra)*minElements;
368 return std::make_pair(start+
nXElems_*xBlock,nume);
374 unsigned int minElements =
nYElems_/size;
375 unsigned int extra =
nYElems_ - minElements*size;
383 nume = minElements+1;
384 start = yProcLoc*(minElements+1);
388 start = extra*(minElements+1)+(yProcLoc-extra)*minElements;
391 return std::make_pair(start+
nYElems_*yBlock,nume);
402 stk::mesh::Part * left = mesh.
getSideset(
"left");
403 stk::mesh::Part * right = mesh.
getSideset(
"right");
404 stk::mesh::Part * top = mesh.
getSideset(
"top");
405 stk::mesh::Part * bottom = mesh.
getSideset(
"bottom");
407 std::vector<stk::mesh::Part*> vertical;
408 std::vector<stk::mesh::Part*> horizontal;
410 std::stringstream ss;
411 ss <<
"vertical_" << bx-1;
412 vertical.push_back(mesh.
getSideset(ss.str()));
415 std::stringstream ss;
416 ss <<
"horizontal_" << by-1;
417 horizontal.push_back(mesh.
getSideset(ss.str()));
420 std::vector<stk::mesh::Entity> localElmts;
424 std::vector<stk::mesh::Entity>::const_iterator itr;
425 for(itr=localElmts.begin();itr!=localElmts.end();++itr) {
426 stk::mesh::Entity element = (*itr);
430 ny = (gid-1) / totalXElems;
431 nx = gid-ny*totalXElems-1;
436 if(nx+1==totalXElems) {
452 if(nx+1!=totalXElems && ((nx+1) %
nXElems_==0)) {
483 if(ny+1==totalYElems) {
501 if(ny+1!=totalYElems && ((ny+1) %
nYElems_==0)) {
520 stk::mesh::Part * lower_left = mesh.
getNodeset(
"lower_left");
521 stk::mesh::Part * origin = mesh.
getNodeset(
"origin");
530 stk::mesh::Entity node = bulkData->get_entity(mesh.
getNodeRank(),1);
549 return Teuchos::tuple(i,j);
void addNodeset(const std::string &name)
void initializeWithDefaults()
unsigned entityOwnerRank(stk::mesh::Entity entity) const
virtual Teuchos::RCP< STK_Interface > buildUncommitedMesh(stk::ParallelMachine parallelMach) const
void addEntityToNodeset(stk::mesh::Entity entity, stk::mesh::Part *nodeset)
Teuchos::Tuple< std::size_t, 2 > procTuple_
T & get(const std::string &name, T def_value)
stk::mesh::Part * getElementBlockPart(const std::string &name) const
get the block count
ParameterList & set(std::string const &name, T const &value, std::string const &docString="", RCP< const ParameterEntryValidator > const &validator=null)
std::pair< int, int > determineXElemSizeAndStart(int xBlock, unsigned int size, unsigned int rank) const
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
void addEntityToSideset(stk::mesh::Entity entity, stk::mesh::Part *sideset)
void buildElements(stk::ParallelMachine parallelMach, STK_Interface &mesh) const
stk::mesh::Part * getNodeset(const std::string &name) const
void initialize(stk::ParallelMachine parallelMach, bool setupIO=true)
stk::mesh::Part * getSideset(const std::string &name) const
void addNodeSets(STK_Interface &mesh) const
std::pair< int, int > determineYElemSizeAndStart(int yBlock, unsigned int size, unsigned int rank) const
void addElement(const Teuchos::RCP< ElementDescriptor > &ed, stk::mesh::Part *block)
Teuchos::RCP< STK_Interface > buildMesh(stk::ParallelMachine parallelMach) const
Build the mesh object.
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
SquareQuadMeshFactory(bool enableRebalance=false)
Constructor.
void rebalance(STK_Interface &mesh) const
Teuchos::RCP< stk::mesh::BulkData > getBulkData() const
void addSideSets(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)
void buildBlock(stk::ParallelMachine machRank, int xBlock, int yBlock, STK_Interface &mesh) const
Teuchos::Tuple< std::size_t, 2 > procRankToProcTuple(std::size_t procRank) const
what is the 2D tuple describe this processor distribution
void addNode(stk::mesh::EntityId gid, const std::vector< double > &coord)
void buildSubcells()
force the mesh to build subcells: edges and faces
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)
void buildMetaData(stk::ParallelMachine parallelMach, STK_Interface &mesh) const
std::vector< Teuchos::RCP< const PeriodicBC_MatcherBase > > periodicBCVec_
~SquareQuadMeshFactory()
Destructor.
void buildLocalElementIDs()
void setParameterList(const Teuchos::RCP< Teuchos::ParameterList > ¶mList)
From ParameterListAcceptor.
ParameterList & sublist(const std::string &name, bool mustAlreadyExist=false, const std::string &docString="")
void getMyElements(std::vector< stk::mesh::Entity > &elements) const
#define TEUCHOS_ASSERT(assertion_test)
stk::mesh::EntityRank getNodeRank() const
virtual void completeMeshConstruction(STK_Interface &mesh, stk::ParallelMachine parallelMach) const
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const
From ParameterListAcceptor.
void addElementBlock(const std::string &name, const CellTopologyData *ctData)
static void parsePeriodicBCList(const Teuchos::RCP< Teuchos::ParameterList > &pl, std::vector< Teuchos::RCP< const PeriodicBC_MatcherBase > > &periodicBC)