Panzer  Version of the Day
Panzer_STK_SquareQuadMeshFactory.cpp
Go to the documentation of this file.
1 // @HEADER
2 // ***********************************************************************
3 //
4 // Panzer: A partial differential equation assembly
5 // engine for strongly coupled complex multiphysics systems
6 // Copyright (2011) Sandia Corporation
7 //
8 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
9 // the U.S. Government retains certain rights in this software.
10 //
11 // Redistribution and use in source and binary forms, with or without
12 // modification, are permitted provided that the following conditions are
13 // met:
14 //
15 // 1. Redistributions of source code must retain the above copyright
16 // notice, this list of conditions and the following disclaimer.
17 //
18 // 2. Redistributions in binary form must reproduce the above copyright
19 // notice, this list of conditions and the following disclaimer in the
20 // documentation and/or other materials provided with the distribution.
21 //
22 // 3. Neither the name of the Corporation nor the names of the
23 // contributors may be used to endorse or promote products derived from
24 // this software without specific prior written permission.
25 //
26 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37 //
38 // Questions? Contact Roger P. Pawlowski (rppawlo@sandia.gov) and
39 // Eric C. Cyr (eccyr@sandia.gov)
40 // ***********************************************************************
41 // @HEADER
42 
44 #include <Teuchos_TimeMonitor.hpp>
45 #include <PanzerAdaptersSTK_config.hpp>
46 #include "Teuchos_StandardParameterEntryValidators.hpp" // for plist validation
47 
48 // #define ENABLE_UNIFORM
49 
50 using Teuchos::RCP;
51 using Teuchos::rcp;
52 
53 namespace panzer_stk {
54 
56 {
58 }
59 
62 {
63 }
64 
66 Teuchos::RCP<STK_Interface> SquareQuadMeshFactory::buildMesh(stk::ParallelMachine parallelMach) const
67 {
68  PANZER_FUNC_TIME_MONITOR("panzer::SquareQuadMeshFactory::buildMesh()");
69 
70  // build all meta data
71  RCP<STK_Interface> mesh = buildUncommitedMesh(parallelMach);
72 
73  // commit meta data
74  mesh->initialize(parallelMach);
75 
76  // build bulk data
77  completeMeshConstruction(*mesh,parallelMach);
78 
79  return mesh;
80 }
81 
82 Teuchos::RCP<STK_Interface> SquareQuadMeshFactory::buildUncommitedMesh(stk::ParallelMachine parallelMach) const
83 {
84  PANZER_FUNC_TIME_MONITOR("panzer::SquareQuadMeshFactory::buildUncomittedMesh()");
85 
86  RCP<STK_Interface> mesh = rcp(new STK_Interface(2));
87 
88  machRank_ = stk::parallel_machine_rank(parallelMach);
89  machSize_ = stk::parallel_machine_size(parallelMach);
90 
91  if (xProcs_ == -1 && yProcs_ == -1) {
92  // copied from galeri
93  xProcs_ = yProcs_ = Teuchos::as<int>(pow(Teuchos::as<double>(machSize_), 0.5));
94 
95  if (xProcs_ * yProcs_ != Teuchos::as<int>(machSize_)) {
96  // Simple method to find a set of processor assignments
97  xProcs_ = yProcs_ = 1;
98 
99  // This means that this works correctly up to about maxFactor^2
100  // processors.
101  const int maxFactor = 100;
102 
103  int ProcTemp = machSize_;
104  int factors[maxFactor];
105  for (int jj = 0; jj < maxFactor; jj++) factors[jj] = 0;
106  for (int jj = 2; jj < maxFactor; jj++) {
107  bool flag = true;
108  while (flag) {
109  int temp = ProcTemp/jj;
110  if (temp*jj == ProcTemp) {
111  factors[jj]++;
112  ProcTemp = temp;
113 
114  } else {
115  flag = false;
116  }
117  }
118  }
119  xProcs_ = ProcTemp;
120  for (int jj = maxFactor-1; jj > 0; jj--) {
121  while (factors[jj] != 0) {
122  if (xProcs_ <= yProcs_) xProcs_ = xProcs_*jj;
123  else yProcs_ = yProcs_*jj;
124  factors[jj]--;
125  }
126  }
127  }
128 
129  } else if(xProcs_==-1) {
130  // default x only decomposition
131  xProcs_ = machSize_;
132  yProcs_ = 1;
133  }
134  TEUCHOS_TEST_FOR_EXCEPTION(int(machSize_) != xProcs_ * yProcs_, std::logic_error,
135  "Cannot build SquareQuadMeshFactory. The product of 'X Procs * Y Procs = " << xProcs_ << "*" << yProcs_ << " = " << xProcs_*yProcs_
136  << "' must equal the number of processors = " << machSize_
137  << "\n\n\t==> Run the simulation with an appropriate number of processors, i.e. #procs = " << xProcs_*yProcs_ << ".\n");
139 
140  // build meta information: blocks and side set setups
141  buildMetaData(parallelMach,*mesh);
142 
143  mesh->addPeriodicBCs(periodicBCVec_);
144 
145  return mesh;
146 }
147 
148 void SquareQuadMeshFactory::completeMeshConstruction(STK_Interface & mesh,stk::ParallelMachine parallelMach) const
149 {
150  PANZER_FUNC_TIME_MONITOR("panzer::SquareQuadMeshFactory::completeMeshConstruction()");
151 
152  if(not mesh.isInitialized())
153  mesh.initialize(parallelMach);
154 
155  // add node and element information
156  buildElements(parallelMach,mesh);
157 
158  // finish up the edges
159 #ifndef ENABLE_UNIFORM
160  mesh.buildSubcells();
161 #endif
162  mesh.buildLocalElementIDs();
163  if(createEdgeBlocks_) {
164  mesh.buildLocalEdgeIDs();
165  }
166 
167  // now that edges are built, sidsets can be added
168 #ifndef ENABLE_UNIFORM
169  addSideSets(mesh);
170 #endif
171 
172  // add nodesets
173  addNodeSets(mesh);
174 
175  if(createEdgeBlocks_) {
176  addEdgeBlocks(mesh);
177  }
178 
179  // calls Stk_MeshFactory::rebalance
180  this->rebalance(mesh);
181 }
182 
184 void SquareQuadMeshFactory::setParameterList(const Teuchos::RCP<Teuchos::ParameterList> & paramList)
185 {
186  paramList->validateParametersAndSetDefaults(*getValidParameters(),0);
187 
188  setMyParamList(paramList);
189 
190  x0_ = paramList->get<double>("X0");
191  y0_ = paramList->get<double>("Y0");
192 
193  xf_ = paramList->get<double>("Xf");
194  yf_ = paramList->get<double>("Yf");
195 
196  xBlocks_ = paramList->get<int>("X Blocks");
197  yBlocks_ = paramList->get<int>("Y Blocks");
198 
199  nXElems_ = paramList->get<int>("X Elements");
200  nYElems_ = paramList->get<int>("Y Elements");
201 
202  xProcs_ = paramList->get<int>("X Procs");
203  yProcs_ = paramList->get<int>("Y Procs");
204 
205  offsetGIDs_ = (paramList->get<std::string>("Offset mesh GIDs above 32-bit int limit") == "ON") ? true : false;
206 
207  createEdgeBlocks_ = paramList->get<bool>("Create Edge Blocks");
208 
209  // read in periodic boundary conditions
210  parsePeriodicBCList(Teuchos::rcpFromRef(paramList->sublist("Periodic BCs")),periodicBCVec_);
211 }
212 
214 Teuchos::RCP<const Teuchos::ParameterList> SquareQuadMeshFactory::getValidParameters() const
215 {
216  static RCP<Teuchos::ParameterList> defaultParams;
217 
218  // fill with default values
219  if(defaultParams == Teuchos::null) {
220  defaultParams = rcp(new Teuchos::ParameterList);
221 
222  defaultParams->set<double>("X0",0.0);
223  defaultParams->set<double>("Y0",0.0);
224 
225  defaultParams->set<double>("Xf",1.0);
226  defaultParams->set<double>("Yf",1.0);
227 
228  defaultParams->set<int>("X Blocks",1);
229  defaultParams->set<int>("Y Blocks",1);
230 
231  defaultParams->set<int>("X Procs",-1);
232  defaultParams->set<int>("Y Procs",1);
233 
234  defaultParams->set<int>("X Elements",5);
235  defaultParams->set<int>("Y Elements",5);
236 
237  // default to false for backward compatibility
238  defaultParams->set<bool>("Create Edge Blocks",false,"Create edge blocks in the mesh");
239 
240  Teuchos::setStringToIntegralParameter<int>(
241  "Offset mesh GIDs above 32-bit int limit",
242  "OFF",
243  "If 64-bit GIDs are supported, the mesh element and node global indices will start at a value greater than 32-bit limit.",
244  Teuchos::tuple<std::string>("OFF", "ON"),
245  defaultParams.get());
246 
247  Teuchos::ParameterList & bcs = defaultParams->sublist("Periodic BCs");
248  bcs.set<int>("Count",0); // no default periodic boundary conditions
249  }
250 
251  return defaultParams;
252 }
253 
255 {
256  // get valid parameters
257  RCP<Teuchos::ParameterList> validParams = rcp(new Teuchos::ParameterList(*getValidParameters()));
258 
259  // set that parameter list
260  setParameterList(validParams);
261 }
262 
263 void SquareQuadMeshFactory::buildMetaData(stk::ParallelMachine /* parallelMach */, STK_Interface & mesh) const
264 {
265  typedef shards::Quadrilateral<4> QuadTopo;
266  const CellTopologyData * ctd = shards::getCellTopologyData<QuadTopo>();
267  const CellTopologyData * side_ctd = shards::CellTopology(ctd).getBaseCellTopologyData(1,0);
268 
269  // build meta data
270  //mesh.setDimension(2);
271  for(int bx=0;bx<xBlocks_;bx++) {
272  for(int by=0;by<yBlocks_;by++) {
273 
274  // add this element block
275  {
276  std::stringstream ebPostfix;
277  ebPostfix << "-" << bx << "_" << by;
278 
279  // add element blocks
280  mesh.addElementBlock("eblock"+ebPostfix.str(),ctd);
281  }
282 
283  }
284  }
285 
286  // add sidesets
287 #ifndef ENABLE_UNIFORM
288  mesh.addSideset("left",side_ctd);
289  mesh.addSideset("right",side_ctd);
290  mesh.addSideset("top",side_ctd);
291  mesh.addSideset("bottom",side_ctd);
292 
293  for(int bx=1;bx<xBlocks_;bx++) {
294  std::stringstream ss;
295  ss << "vertical_" << bx-1;
296  mesh.addSideset(ss.str(),side_ctd);
297  }
298  for(int by=1;by<yBlocks_;by++) {
299  std::stringstream ss;
300  ss << "horizontal_" << by-1;
301  mesh.addSideset(ss.str(),side_ctd);
302  }
303 #endif
304 
305  // add nodesets
306  mesh.addNodeset("lower_left");
307  mesh.addNodeset("origin");
308 
309  if(createEdgeBlocks_) {
310  const CellTopologyData * edge_ctd = shards::CellTopology(ctd).getBaseCellTopologyData(1,0);
312  }
313 }
314 
315 void SquareQuadMeshFactory::buildElements(stk::ParallelMachine parallelMach,STK_Interface & mesh) const
316 {
317  mesh.beginModification();
318  // build each block
319  for(int xBlock=0;xBlock<xBlocks_;xBlock++) {
320  for(int yBlock=0;yBlock<yBlocks_;yBlock++) {
321  buildBlock(parallelMach,xBlock,yBlock,mesh);
322  }
323  }
324  mesh.endModification();
325 }
326 
327 void SquareQuadMeshFactory::buildBlock(stk::ParallelMachine /* parallelMach */,int xBlock,int yBlock,STK_Interface & mesh) const
328 {
329  // grab this processors rank and machine size
330  std::pair<int,int> sizeAndStartX = determineXElemSizeAndStart(xBlock,xProcs_,machRank_);
331  std::pair<int,int> sizeAndStartY = determineYElemSizeAndStart(yBlock,yProcs_,machRank_);
332 
333  int myXElems_start = sizeAndStartX.first;
334  int myXElems_end = myXElems_start+sizeAndStartX.second;
335  int myYElems_start = sizeAndStartY.first;
336  int myYElems_end = myYElems_start+sizeAndStartY.second;
337  int totalXElems = nXElems_*xBlocks_;
338  int totalYElems = nYElems_*yBlocks_;
339 
340  double deltaX = (xf_-x0_)/double(totalXElems);
341  double deltaY = (yf_-y0_)/double(totalYElems);
342 
343  std::vector<double> coord(2,0.0);
344 
345  offset_ = 0;
346  if (offsetGIDs_) {
347  if (std::numeric_limits<panzer::GlobalOrdinal>::max() > std::numeric_limits<unsigned int>::max())
348  offset_ = panzer::GlobalOrdinal(std::numeric_limits<unsigned int>::max()) + 1;
349  }
350 
351  // build the nodes
352  for(int nx=myXElems_start;nx<myXElems_end+1;++nx) {
353  coord[0] = this->getMeshCoord(nx, deltaX, x0_);
354  for(int ny=myYElems_start;ny<myYElems_end+1;++ny) {
355  coord[1] = this->getMeshCoord(ny, deltaY, y0_);
356 
357  mesh.addNode(ny*(totalXElems+1)+nx+1+offset_,coord);
358  }
359  }
360 
361  std::stringstream blockName;
362  blockName << "eblock-" << xBlock << "_" << yBlock;
363  stk::mesh::Part * block = mesh.getElementBlockPart(blockName.str());
364 
365  // build the elements
366  for(int nx=myXElems_start;nx<myXElems_end;++nx) {
367  for(int ny=myYElems_start;ny<myYElems_end;++ny) {
368  stk::mesh::EntityId gid = totalXElems*ny+nx+1 + offset_;
369  std::vector<stk::mesh::EntityId> nodes(4);
370  nodes[0] = nx+1+ny*(totalXElems+1);
371  nodes[1] = nodes[0]+1;
372  nodes[2] = nodes[1]+(totalXElems+1);
373  nodes[3] = nodes[2]-1;
374 
375  for (int i=0; i < 4; ++i)
376  nodes[i] += offset_;
377 
378  RCP<ElementDescriptor> ed = rcp(new ElementDescriptor(gid,nodes));
379  mesh.addElement(ed,block);
380  }
381  }
382 }
383 
384 std::pair<int,int> SquareQuadMeshFactory::determineXElemSizeAndStart(int xBlock,unsigned int size,unsigned int /* rank */) const
385 {
386  std::size_t xProcLoc = procTuple_[0];
387  unsigned int minElements = nXElems_/size;
388  unsigned int extra = nXElems_ - minElements*size;
389 
390  TEUCHOS_ASSERT(minElements>0);
391 
392  // first "extra" elements get an extra column of elements
393  // this determines the starting X index and number of elements
394  int nume=0, start=0;
395  if(xProcLoc<extra) {
396  nume = minElements+1;
397  start = xProcLoc*(minElements+1);
398  }
399  else {
400  nume = minElements;
401  start = extra*(minElements+1)+(xProcLoc-extra)*minElements;
402  }
403 
404  return std::make_pair(start+nXElems_*xBlock,nume);
405 }
406 
407 std::pair<int,int> SquareQuadMeshFactory::determineYElemSizeAndStart(int yBlock,unsigned int size,unsigned int /* rank */) const
408 {
409  std::size_t yProcLoc = procTuple_[1];
410  unsigned int minElements = nYElems_/size;
411  unsigned int extra = nYElems_ - minElements*size;
412 
413  TEUCHOS_ASSERT(minElements>0);
414 
415  // first "extra" elements get an extra column of elements
416  // this determines the starting X index and number of elements
417  int nume=0, start=0;
418  if(yProcLoc<extra) {
419  nume = minElements+1;
420  start = yProcLoc*(minElements+1);
421  }
422  else {
423  nume = minElements;
424  start = extra*(minElements+1)+(yProcLoc-extra)*minElements;
425  }
426 
427  return std::make_pair(start+nYElems_*yBlock,nume);
428 }
429 
431 {
432  mesh.beginModification();
433 
434  std::size_t totalXElems = nXElems_*xBlocks_;
435  std::size_t totalYElems = nYElems_*yBlocks_;
436 
437  // get all part vectors
438  stk::mesh::Part * left = mesh.getSideset("left");
439  stk::mesh::Part * right = mesh.getSideset("right");
440  stk::mesh::Part * top = mesh.getSideset("top");
441  stk::mesh::Part * bottom = mesh.getSideset("bottom");
442 
443  std::vector<stk::mesh::Part*> vertical;
444  std::vector<stk::mesh::Part*> horizontal;
445  for(int bx=1;bx<xBlocks_;bx++) {
446  std::stringstream ss;
447  ss << "vertical_" << bx-1;
448  vertical.push_back(mesh.getSideset(ss.str()));
449  }
450  for(int by=1;by<yBlocks_;by++) {
451  std::stringstream ss;
452  ss << "horizontal_" << by-1;
453  horizontal.push_back(mesh.getSideset(ss.str()));
454  }
455 
456  std::vector<stk::mesh::Entity> localElmts;
457  mesh.getMyElements(localElmts);
458 
459  // loop over elements adding edges to sidesets
460  std::vector<stk::mesh::Entity>::const_iterator itr;
461  for(itr=localElmts.begin();itr!=localElmts.end();++itr) {
462  stk::mesh::Entity element = (*itr);
463  stk::mesh::EntityId gid = mesh.elementGlobalId(element);
464 
465  // reverse the offset for local gid numbering scheme
466  gid -= offset_;
467 
468  std::size_t nx,ny;
469  ny = (gid-1) / totalXElems;
470  nx = gid-ny*totalXElems-1;
471 
472  // vertical boundaries
474 
475  if(nx+1==totalXElems) {
476  stk::mesh::Entity edge = mesh.findConnectivityById(element, stk::topology::EDGE_RANK, 1);
477 
478  // on the right
479  if(mesh.entityOwnerRank(edge)==machRank_)
480  mesh.addEntityToSideset(edge,right);
481  }
482 
483  if(nx==0) {
484  stk::mesh::Entity edge = mesh.findConnectivityById(element, stk::topology::EDGE_RANK, 3);
485 
486  // on the left
487  if(mesh.entityOwnerRank(edge)==machRank_)
488  mesh.addEntityToSideset(edge,left);
489  }
490 
491  if(nx+1!=totalXElems && ((nx+1) % nXElems_==0)) {
492  stk::mesh::Entity edge = mesh.findConnectivityById(element, stk::topology::EDGE_RANK, 1);
493 
494  // on the right
495  if(mesh.entityOwnerRank(edge)==machRank_) {
496  int index = (nx+1)/nXElems_-1;
497  mesh.addEntityToSideset(edge,vertical[index]);
498  }
499  }
500 
501  if(nx!=0 && (nx % nXElems_==0)) {
502  stk::mesh::Entity edge = mesh.findConnectivityById(element, stk::topology::EDGE_RANK, 3);
503 
504  // on the left
505  if(mesh.entityOwnerRank(edge)==machRank_) {
506  int index = nx/nXElems_-1;
507  mesh.addEntityToSideset(edge,vertical[index]);
508  }
509  }
510 
511  // horizontal boundaries
513 
514  if(ny==0) {
515  stk::mesh::Entity edge = mesh.findConnectivityById(element, stk::topology::EDGE_RANK, 0);
516 
517  // on the bottom
518  if(mesh.entityOwnerRank(edge)==machRank_)
519  mesh.addEntityToSideset(edge,bottom);
520  }
521 
522  if(ny+1==totalYElems) {
523  stk::mesh::Entity edge = mesh.findConnectivityById(element, stk::topology::EDGE_RANK, 2);
524 
525  // on the top
526  if(mesh.entityOwnerRank(edge)==machRank_)
527  mesh.addEntityToSideset(edge,top);
528  }
529 
530  if(ny!=0 && (ny % nYElems_==0)) {
531  stk::mesh::Entity edge = mesh.findConnectivityById(element, stk::topology::EDGE_RANK, 0);
532 
533  // on the bottom
534  if(mesh.entityOwnerRank(edge)==machRank_) {
535  int index = ny/nYElems_-1;
536  mesh.addEntityToSideset(edge,horizontal[index]);
537  }
538  }
539 
540  if(ny+1!=totalYElems && ((ny+1) % nYElems_==0)) {
541  stk::mesh::Entity edge = mesh.findConnectivityById(element, stk::topology::EDGE_RANK, 2);
542 
543  // on the top
544  if(mesh.entityOwnerRank(edge)==machRank_) {
545  int index = (ny+1)/nYElems_-1;
546  mesh.addEntityToSideset(edge,horizontal[index]);
547  }
548  }
549  }
550 
551  mesh.endModification();
552 }
553 
555 {
556  mesh.beginModification();
557 
558  // get all part vectors
559  stk::mesh::Part * lower_left = mesh.getNodeset("lower_left");
560  stk::mesh::Part * origin = mesh.getNodeset("origin");
561 
562  // std::vector<stk::mesh::Entity> localElmts;
563  // mesh.getMyElements(localElmts);
564 
565  Teuchos::RCP<stk::mesh::BulkData> bulkData = mesh.getBulkData();
566  if(machRank_==0)
567  {
568  // add zero node to lower_left node set
569  stk::mesh::Entity node = bulkData->get_entity(mesh.getNodeRank(),1 + offset_);
570  mesh.addEntityToNodeset(node,lower_left);
571 
572  // add zero node to origin node set
573  mesh.addEntityToNodeset(node,origin);
574  }
575 
576  mesh.endModification();
577 }
578 
580 {
581  mesh.beginModification();
582 
583  stk::mesh::Part * edge_block = mesh.getEdgeBlock(panzer_stk::STK_Interface::edgeBlockString);
584 
585  Teuchos::RCP<stk::mesh::BulkData> bulkData = mesh.getBulkData();
586  Teuchos::RCP<stk::mesh::MetaData> metaData = mesh.getMetaData();
587 
588  std::vector<stk::mesh::Entity> edges;
589  bulkData->get_entities(mesh.getEdgeRank(),metaData->locally_owned_part(),edges);
590  mesh.addEntitiesToEdgeBlock(edges, edge_block);
591 
592  mesh.endModification();
593 }
594 
596 Teuchos::Tuple<std::size_t,2> SquareQuadMeshFactory::procRankToProcTuple(std::size_t procRank) const
597 {
598  std::size_t i=0,j=0;
599 
600  j = procRank/xProcs_;
601  procRank = procRank % xProcs_;
602  i = procRank;
603 
604  return Teuchos::tuple(i,j);
605 }
606 
607 } // end panzer_stk
Teuchos::RCP< stk::mesh::MetaData > getMetaData() const
stk::mesh::Part * getElementBlockPart(const std::string &name) const
get the block part
void buildMetaData(stk::ParallelMachine parallelMach, STK_Interface &mesh) const
void addNodeset(const std::string &name)
void buildBlock(stk::ParallelMachine machRank, int xBlock, int yBlock, STK_Interface &mesh) const
void addEntityToNodeset(stk::mesh::Entity entity, stk::mesh::Part *nodeset)
void getMyElements(std::vector< stk::mesh::Entity > &elements) const
void addEntityToSideset(stk::mesh::Entity entity, stk::mesh::Part *sideset)
stk::mesh::EntityRank getNodeRank() const
Teuchos::RCP< stk::mesh::BulkData > getBulkData() const
virtual Teuchos::RCP< STK_Interface > buildUncommitedMesh(stk::ParallelMachine parallelMach) const
stk::mesh::Entity findConnectivityById(stk::mesh::Entity src, stk::mesh::EntityRank tgt_rank, unsigned rel_id) const
bool offsetGIDs_
If true, offset mesh GIDs to exercise 32-bit limits.
stk::mesh::Part * getEdgeBlock(const std::string &name) const
get the block part
stk::mesh::EntityRank getEdgeRank() const
stk::mesh::Part * getSideset(const std::string &name) const
void addElement(const Teuchos::RCP< ElementDescriptor > &ed, stk::mesh::Part *block)
stk::mesh::EntityId elementGlobalId(std::size_t lid) const
void initialize(stk::ParallelMachine parallelMach, bool setupIO=true, const bool buildRefinementSupport=false)
SquareQuadMeshFactory(bool enableRebalance=false)
Constructor.
virtual void completeMeshConstruction(STK_Interface &mesh, stk::ParallelMachine parallelMach) const
double getMeshCoord(const int nx, const double deltaX, const double x0) const
bool isInitialized() const
Has initialize been called on this mesh object?
static const std::string edgeBlockString
std::pair< int, int > determineYElemSizeAndStart(int yBlock, unsigned int size, unsigned int rank) 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
void buildElements(stk::ParallelMachine parallelMach, STK_Interface &mesh) const
void addSideset(const std::string &name, const CellTopologyData *ctData)
void addEdgeBlock(const std::string &name, const CellTopologyData *ctData)
std::pair< int, int > determineXElemSizeAndStart(int xBlock, unsigned int size, unsigned int rank) const
std::vector< Teuchos::RCP< const PeriodicBC_MatcherBase > > periodicBCVec_
Teuchos::Tuple< std::size_t, 2 > procRankToProcTuple(std::size_t procRank) const
what is the 2D tuple describe this processor distribution
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const
From ParameterListAcceptor.
void setParameterList(const Teuchos::RCP< Teuchos::ParameterList > &paramList)
From ParameterListAcceptor.
Teuchos::RCP< STK_Interface > buildMesh(stk::ParallelMachine parallelMach) const
Build the mesh object.
unsigned entityOwnerRank(stk::mesh::Entity entity) const
void rebalance(STK_Interface &mesh) const
void addElementBlock(const std::string &name, const CellTopologyData *ctData)
stk::mesh::Part * getNodeset(const std::string &name) const
static void parsePeriodicBCList(const Teuchos::RCP< Teuchos::ParameterList > &pl, std::vector< Teuchos::RCP< const PeriodicBC_MatcherBase > > &periodicBC)