Panzer  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Panzer_STK_CubeTetMeshFactory.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 
47 using Teuchos::RCP;
48 using Teuchos::rcp;
49 
50 namespace panzer_stk {
51 
53 {
55 }
56 
59 {
60 }
61 
63 Teuchos::RCP<STK_Interface> CubeTetMeshFactory::buildMesh(stk::ParallelMachine parallelMach) const
64 {
65  PANZER_FUNC_TIME_MONITOR("panzer::CubeTetMeshFactory::buildMesh()");
66 
67  // build all meta data
68  RCP<STK_Interface> mesh = buildUncommitedMesh(parallelMach);
69 
70  // commit meta data
71  mesh->initialize(parallelMach);
72 
73  // build bulk data
74  completeMeshConstruction(*mesh,parallelMach);
75 
76  return mesh;
77 }
78 
80 {
81  PANZER_FUNC_TIME_MONITOR("panzer::CubeTetMeshFactory::buildUncomittedMesh()");
82 
83  RCP<STK_Interface> mesh = rcp(new STK_Interface(3));
84 
85  machRank_ = stk::parallel_machine_rank(parallelMach);
86  machSize_ = stk::parallel_machine_size(parallelMach);
87 
88  if (xProcs_ == -1 && yProcs_ == -1 && zProcs_ == -1) {
89  // copied from galeri
90  xProcs_ = yProcs_ = zProcs_ = Teuchos::as<int>(pow(Teuchos::as<double>(machSize_), 0.333334));
91 
92  if (xProcs_ * yProcs_ * zProcs_ != Teuchos::as<int>(machSize_)) {
93  // Simple method to find a set of processor assignments
94  xProcs_ = yProcs_ = zProcs_ = 1;
95 
96  // This means that this works correctly up to about maxFactor^3
97  // processors.
98  const int maxFactor = 50;
99 
100  int ProcTemp = machSize_;
101  int factors[maxFactor];
102  for (int jj = 0; jj < maxFactor; jj++) factors[jj] = 0;
103  for (int jj = 2; jj < maxFactor; jj++) {
104  bool flag = true;
105  while (flag) {
106  int temp = ProcTemp/jj;
107  if (temp*jj == ProcTemp) {
108  factors[jj]++;
109  ProcTemp = temp;
110 
111  } else {
112  flag = false;
113  }
114  }
115  }
116  xProcs_ = ProcTemp;
117  for (int jj = maxFactor-1; jj > 0; jj--) {
118  while (factors[jj] != 0) {
119  if ((xProcs_ <= yProcs_) && (xProcs_ <= zProcs_)) xProcs_ = xProcs_*jj;
120  else if ((yProcs_ <= xProcs_) && (yProcs_ <= zProcs_)) yProcs_ = yProcs_*jj;
121  else zProcs_ = zProcs_*jj;
122  factors[jj]--;
123  }
124  }
125  }
126 
127  } else if(xProcs_==-1) {
128  // default x only decomposition
129  xProcs_ = machSize_;
130  yProcs_ = 1;
131  zProcs_ = 1;
132  }
134  "Cannot build CubeTetMeshFactory, the product of \"X Procs\", \"Y Procs\", and \"Z Procs\""
135  " must equal the number of processors.");
137 
138  // build meta information: blocks and side set setups
139  buildMetaData(parallelMach,*mesh);
140 
141  mesh->addPeriodicBCs(periodicBCVec_);
142  mesh->setBoundingBoxSearchFlag(useBBoxSearch_);
143 
144  return mesh;
145 }
146 
147 void CubeTetMeshFactory::completeMeshConstruction(STK_Interface & mesh,stk::ParallelMachine parallelMach) const
148 {
149  PANZER_FUNC_TIME_MONITOR("panzer::CubeTetMeshFactory::completeMeshConstruction()");
150 
151  if(not mesh.isInitialized())
152  mesh.initialize(parallelMach);
153 
154  // add node and element information
155  buildElements(parallelMach,mesh);
156 
157  // finish up the edges and faces
158  mesh.buildSubcells();
159  mesh.buildLocalElementIDs();
160  if(createEdgeBlocks_) {
161  mesh.buildLocalEdgeIDs();
162  }
163  if(createFaceBlocks_) {
164  mesh.buildLocalFaceIDs();
165  }
166 
167  // now that edges are built, sidets can be added
168  addSideSets(mesh);
169  addNodeSets(mesh);
170 
171  mesh.beginModification();
172  if(createEdgeBlocks_) {
173  addEdgeBlocks(mesh);
174  }
175  if(createFaceBlocks_) {
176  addFaceBlocks(mesh);
177  }
178  mesh.endModification();
179 
180  // calls Stk_MeshFactory::rebalance
181  this->rebalance(mesh);
182 }
183 
186 {
188 
189  setMyParamList(paramList);
190 
191  x0_ = paramList->get<double>("X0");
192  y0_ = paramList->get<double>("Y0");
193  z0_ = paramList->get<double>("Z0");
194 
195  xf_ = paramList->get<double>("Xf");
196  yf_ = paramList->get<double>("Yf");
197  zf_ = paramList->get<double>("Zf");
198 
199  xBlocks_ = paramList->get<int>("X Blocks");
200  yBlocks_ = paramList->get<int>("Y Blocks");
201  zBlocks_ = paramList->get<int>("Z Blocks");
202 
203  xProcs_ = paramList->get<int>("X Procs");
204  yProcs_ = paramList->get<int>("Y Procs");
205  zProcs_ = paramList->get<int>("Z Procs");
206 
207  nXElems_ = paramList->get<int>("X Elements");
208  nYElems_ = paramList->get<int>("Y Elements");
209  nZElems_ = paramList->get<int>("Z Elements");
210 
211  createEdgeBlocks_ = paramList->get<bool>("Create Edge Blocks");
212  createFaceBlocks_ = paramList->get<bool>("Create Face Blocks");
213 
214  // read in periodic boundary conditions
215  parsePeriodicBCList(Teuchos::rcpFromRef(paramList->sublist("Periodic BCs")),periodicBCVec_,useBBoxSearch_);
216 }
217 
220 {
221  static RCP<Teuchos::ParameterList> defaultParams;
222 
223  // fill with default values
224  if(defaultParams == Teuchos::null) {
225  defaultParams = rcp(new Teuchos::ParameterList);
226 
227  defaultParams->set<double>("X0",0.0);
228  defaultParams->set<double>("Y0",0.0);
229  defaultParams->set<double>("Z0",0.0);
230 
231  defaultParams->set<double>("Xf",1.0);
232  defaultParams->set<double>("Yf",1.0);
233  defaultParams->set<double>("Zf",1.0);
234 
235  defaultParams->set<int>("X Blocks",1);
236  defaultParams->set<int>("Y Blocks",1);
237  defaultParams->set<int>("Z Blocks",1);
238 
239  defaultParams->set<int>("X Procs",-1);
240  defaultParams->set<int>("Y Procs",1);
241  defaultParams->set<int>("Z Procs",1);
242 
243  defaultParams->set<int>("X Elements",5);
244  defaultParams->set<int>("Y Elements",5);
245  defaultParams->set<int>("Z Elements",5);
246 
247  // default to false for backward compatibility
248  defaultParams->set<bool>("Create Edge Blocks",false,"Create edge blocks in the mesh");
249  defaultParams->set<bool>("Create Face Blocks",false,"Create face blocks in the mesh");
250 
251  Teuchos::ParameterList & bcs = defaultParams->sublist("Periodic BCs");
252  bcs.set<int>("Count",0); // no default periodic boundary conditions
253  }
254 
255  return defaultParams;
256 }
257 
259 {
260  // get valid parameters
262 
263  // set that parameter list
264  setParameterList(validParams);
265 
266  /* This is a tet mesh factory so all elements in all element blocks
267  * will be tet4. This means that all the edges will be line2 and
268  * all the faces will be tri3. The edge and face block names are
269  * hard coded to reflect this.
270  */
273 }
274 
275 void CubeTetMeshFactory::buildMetaData(stk::ParallelMachine /* parallelMach */, STK_Interface & mesh) const
276 {
277  typedef shards::Tetrahedron<4> TetTopo;
278  const CellTopologyData * ctd = shards::getCellTopologyData<TetTopo>();
279  const CellTopologyData * side_ctd = shards::CellTopology(ctd).getBaseCellTopologyData(2,0);
280  const CellTopologyData * edge_ctd = shards::CellTopology(ctd).getBaseCellTopologyData(1,0);
281  const CellTopologyData * face_ctd = shards::CellTopology(ctd).getBaseCellTopologyData(2,0);
282 
283  // build meta data
284  //mesh.setDimension(2);
285  for(int bx=0;bx<xBlocks_;bx++) {
286  for(int by=0;by<yBlocks_;by++) {
287  for(int bz=0;bz<zBlocks_;bz++) {
288 
289  std::stringstream ebPostfix;
290  ebPostfix << "-" << bx << "_" << by << "_" << bz;
291 
292  // add element blocks
293  mesh.addElementBlock("eblock"+ebPostfix.str(),ctd);
294  if(createEdgeBlocks_) {
295  mesh.addEdgeBlock("eblock"+ebPostfix.str(),
297  edge_ctd);
298  }
299  if(createFaceBlocks_) {
300  mesh.addFaceBlock("eblock"+ebPostfix.str(),
302  face_ctd);
303  }
304  }
305  }
306  }
307 
308  // add sidesets
309  mesh.addSideset("left",side_ctd);
310  mesh.addSideset("right",side_ctd);
311  mesh.addSideset("top",side_ctd);
312  mesh.addSideset("bottom",side_ctd);
313  mesh.addSideset("front",side_ctd);
314  mesh.addSideset("back",side_ctd);
315 
316  mesh.addNodeset("origin");
317 }
318 
319 void CubeTetMeshFactory::buildElements(stk::ParallelMachine parallelMach,STK_Interface & mesh) const
320 {
321  mesh.beginModification();
322  // build each block
323  for(int xBlock=0;xBlock<xBlocks_;xBlock++) {
324  for(int yBlock=0;yBlock<yBlocks_;yBlock++) {
325  for(int zBlock=0;zBlock<zBlocks_;zBlock++) {
326  buildBlock(parallelMach,xBlock,yBlock,zBlock,mesh);
327  }
328  }
329  }
330  mesh.endModification();
331 }
332 
333 void CubeTetMeshFactory::buildBlock(stk::ParallelMachine /* parallelMach */,int xBlock,int yBlock,int zBlock,STK_Interface & mesh) const
334 {
335  // grab this processors rank and machine size
336  std::pair<int,int> sizeAndStartX = determineXElemSizeAndStart(xBlock,xProcs_,machRank_);
337  std::pair<int,int> sizeAndStartY = determineYElemSizeAndStart(yBlock,yProcs_,machRank_);
338  std::pair<int,int> sizeAndStartZ = determineZElemSizeAndStart(zBlock,zProcs_,machRank_);
339 
340  int myXElems_start = sizeAndStartX.first;
341  int myXElems_end = myXElems_start+sizeAndStartX.second;
342  int myYElems_start = sizeAndStartY.first;
343  int myYElems_end = myYElems_start+sizeAndStartY.second;
344  int myZElems_start = sizeAndStartZ.first;
345  int myZElems_end = myZElems_start+sizeAndStartZ.second;
346 
347  int totalXElems = nXElems_*xBlocks_;
348  int totalYElems = nYElems_*yBlocks_;
349  int totalZElems = nZElems_*zBlocks_;
350 
351  double deltaX = (xf_-x0_)/double(totalXElems);
352  double deltaY = (yf_-y0_)/double(totalYElems);
353  double deltaZ = (zf_-z0_)/double(totalZElems);
354 
355  std::vector<double> coord(3,0.0);
356 
357  // build the nodes
358  for(int nx=myXElems_start;nx<myXElems_end+1;++nx) {
359  coord[0] = this->getMeshCoord(nx, deltaX, x0_);
360  for(int ny=myYElems_start;ny<myYElems_end+1;++ny) {
361  coord[1] = this->getMeshCoord(ny, deltaY, y0_);
362  for(int nz=myZElems_start;nz<myZElems_end+1;++nz) {
363  coord[2] = this->getMeshCoord(nz, deltaZ, z0_);
364 
365  mesh.addNode(nz*(totalYElems+1)*(totalXElems+1)+ny*(totalXElems+1)+nx+1,coord);
366  }
367  }
368  }
369 
370  std::stringstream blockName;
371  blockName << "eblock-" << xBlock << "_" << yBlock << "_" << zBlock;
372  stk::mesh::Part * block = mesh.getElementBlockPart(blockName.str());
373 
374  // build the elements
375  for(int nx=myXElems_start;nx<myXElems_end;++nx) {
376  for(int ny=myYElems_start;ny<myYElems_end;++ny) {
377  for(int nz=myZElems_start;nz<myZElems_end;++nz) {
378 
379  std::vector<stk::mesh::EntityId> nodes(8);
380  nodes[0] = nx+1+ny*(totalXElems+1) +nz*(totalYElems+1)*(totalXElems+1);
381  nodes[1] = nodes[0]+1;
382  nodes[2] = nodes[1]+(totalXElems+1);
383  nodes[3] = nodes[2]-1;
384  nodes[4] = nodes[0]+(totalYElems+1)*(totalXElems+1);
385  nodes[5] = nodes[1]+(totalYElems+1)*(totalXElems+1);
386  nodes[6] = nodes[2]+(totalYElems+1)*(totalXElems+1);
387  nodes[7] = nodes[3]+(totalYElems+1)*(totalXElems+1);
388 
389  buildTetsOnHex(Teuchos::tuple(totalXElems,totalYElems,totalZElems),
390  Teuchos::tuple(nx,ny,nz),
391  block,nodes,mesh);
392  }
393  }
394  }
395 }
396 
398  const Teuchos::Tuple<int,3> & element,
399  stk::mesh::Part * block,
400  const std::vector<stk::mesh::EntityId> & h_nodes,
401  STK_Interface & mesh) const
402 {
403  Teuchos::FancyOStream out(Teuchos::rcpFromRef(std::cout));
404  out.setShowProcRank(true);
405  out.setOutputToRootOnly(-1);
406 
407  int totalXElems = meshDesc[0]; int totalYElems = meshDesc[1]; int totalZElems = meshDesc[2];
408  int nx = element[0]; int ny = element[1]; int nz = element[2];
409 
410  stk::mesh::EntityId hex_id = totalXElems*totalYElems*nz+totalXElems*ny+nx+1;
411  stk::mesh::EntityId gid_0 = 12*(hex_id-1)+1;
412  std::vector<stk::mesh::EntityId> nodes(4);
413 
414  // add centroid node
415  stk::mesh::EntityId centroid = 0;
416  {
417  stk::mesh::EntityId largestNode = (totalXElems+1)*(totalYElems+1)*(totalZElems+1);
418  centroid = hex_id+largestNode;
419 
420  // compute average of coordinates
421  std::vector<double> coord(3,0.0);
422  for(std::size_t i=0;i<h_nodes.size();i++) {
423  const double * node_coord = mesh.getNodeCoordinates(h_nodes[i]);
424  coord[0] += node_coord[0];
425  coord[1] += node_coord[1];
426  coord[2] += node_coord[2];
427  }
428  coord[0] /= 8.0;
429  coord[1] /= 8.0;
430  coord[2] /= 8.0;
431 
432  mesh.addNode(centroid,coord);
433  }
434 
435  //
436  int idSet[][3] = { { 0, 1, 2}, // back
437  { 0, 2, 3},
438  { 0, 5, 1}, // bottom
439  { 0, 4, 5},
440  { 0, 7, 4}, // left
441  { 0, 3, 7},
442  { 6, 1, 5}, // right
443  { 6, 2, 1},
444  { 6, 3, 2}, // top
445  { 6, 7, 3},
446  { 6, 4, 7}, // front
447  { 6, 5, 4} };
448 
449  for(int i=0;i<12;i++) {
450  nodes[0] = h_nodes[idSet[i][0]];
451  nodes[1] = h_nodes[idSet[i][1]];
452  nodes[2] = h_nodes[idSet[i][2]];
453  nodes[3] = centroid;
454 
455  // add element to mesh
456  mesh.addElement(rcp(new ElementDescriptor(gid_0+i,nodes)),block);
457  }
458 }
459 
460 std::pair<int,int> CubeTetMeshFactory::determineXElemSizeAndStart(int xBlock,unsigned int size,unsigned int /* rank */) const
461 {
462  std::size_t xProcLoc = procTuple_[0];
463  unsigned int minElements = nXElems_/size;
464  unsigned int extra = nXElems_ - minElements*size;
465 
466  TEUCHOS_ASSERT(minElements>0);
467 
468  // first "extra" elements get an extra column of elements
469  // this determines the starting X index and number of elements
470  int nume=0, start=0;
471  if(xProcLoc<extra) {
472  nume = minElements+1;
473  start = xProcLoc*(minElements+1);
474  }
475  else {
476  nume = minElements;
477  start = extra*(minElements+1)+(xProcLoc-extra)*minElements;
478  }
479 
480  return std::make_pair(start+nXElems_*xBlock,nume);
481 }
482 
483 std::pair<int,int> CubeTetMeshFactory::determineYElemSizeAndStart(int yBlock,unsigned int size,unsigned int /* rank */) const
484 {
485  // int start = yBlock*nYElems_;
486  // return std::make_pair(start,nYElems_);
487 
488  std::size_t yProcLoc = procTuple_[1];
489  unsigned int minElements = nYElems_/size;
490  unsigned int extra = nYElems_ - minElements*size;
491 
492  TEUCHOS_ASSERT(minElements>0);
493 
494  // first "extra" elements get an extra column of elements
495  // this determines the starting X index and number of elements
496  int nume=0, start=0;
497  if(yProcLoc<extra) {
498  nume = minElements+1;
499  start = yProcLoc*(minElements+1);
500  }
501  else {
502  nume = minElements;
503  start = extra*(minElements+1)+(yProcLoc-extra)*minElements;
504  }
505 
506  return std::make_pair(start+nYElems_*yBlock,nume);
507 }
508 
509 std::pair<int,int> CubeTetMeshFactory::determineZElemSizeAndStart(int zBlock,unsigned int size,unsigned int /* rank */) const
510 {
511  // int start = zBlock*nZElems_;
512  // return std::make_pair(start,nZElems_);
513  std::size_t zProcLoc = procTuple_[2];
514  unsigned int minElements = nZElems_/size;
515  unsigned int extra = nZElems_ - minElements*size;
516 
517  TEUCHOS_ASSERT(minElements>0);
518 
519  // first "extra" elements get an extra column of elements
520  // this determines the starting X index and number of elements
521  int nume=0, start=0;
522  if(zProcLoc<extra) {
523  nume = minElements+1;
524  start = zProcLoc*(minElements+1);
525  }
526  else {
527  nume = minElements;
528  start = extra*(minElements+1)+(zProcLoc-extra)*minElements;
529  }
530 
531  return std::make_pair(start+nZElems_*zBlock,nume);
532 }
533 
535 {
536  mesh.beginModification();
537  const stk::mesh::EntityRank side_rank = mesh.getSideRank();
538 
539  std::size_t totalXElems = nXElems_*xBlocks_;
540  std::size_t totalYElems = nYElems_*yBlocks_;
541  std::size_t totalZElems = nZElems_*zBlocks_;
542 
543  // get all part vectors
544  stk::mesh::Part * left = mesh.getSideset("left");
545  stk::mesh::Part * right = mesh.getSideset("right");
546  stk::mesh::Part * top = mesh.getSideset("top");
547  stk::mesh::Part * bottom = mesh.getSideset("bottom");
548  stk::mesh::Part * front = mesh.getSideset("front");
549  stk::mesh::Part * back = mesh.getSideset("back");
550 
551  std::vector<stk::mesh::Entity> localElmts;
552  mesh.getMyElements(localElmts);
553 
554  // gid = totalXElems*totalYElems*nz+totalXElems*ny+nx+1
555 
556  // loop over elements adding sides to sidesets
557  std::vector<stk::mesh::Entity>::const_iterator itr;
558  for(itr=localElmts.begin();itr!=localElmts.end();++itr) {
559  stk::mesh::Entity element = (*itr);
560  stk::mesh::EntityId gid = mesh.elementGlobalId(element);
561 
562  // get hex global id
563  stk::mesh::EntityId h_gid = (gid-1)/12+1;
564  stk::mesh::EntityId t_offset = gid - (12*(h_gid-1)+1);
565 
566  std::size_t nx,ny,nz;
567  nz = (h_gid-1) / (totalXElems*totalYElems);
568  h_gid = (h_gid-1)-nz*(totalXElems*totalYElems);
569  ny = h_gid / totalXElems;
570  nx = h_gid-ny*totalXElems;
571 
572  if(nz==0 && (t_offset==0 || t_offset==1)) {
573  stk::mesh::Entity side = mesh.findConnectivityById(element, side_rank, 3);
574 
575  // on the back
576  if(mesh.entityOwnerRank(side)==machRank_)
577  mesh.addEntityToSideset(side,back);
578  }
579  if(nz+1==totalZElems && (t_offset==10 || t_offset==11)) {
580  stk::mesh::Entity side = mesh.findConnectivityById(element, side_rank, 3);
581 
582  // on the front
583  if(mesh.entityOwnerRank(side)==machRank_)
584  mesh.addEntityToSideset(side,front);
585  }
586 
587  if(ny==0 && (t_offset==2 || t_offset==3)) {
588  stk::mesh::Entity side = mesh.findConnectivityById(element, side_rank, 3);
589 
590  // on the bottom
591  if(mesh.entityOwnerRank(side)==machRank_)
592  mesh.addEntityToSideset(side,bottom);
593  }
594  if(ny+1==totalYElems && (t_offset==8 || t_offset==9)) {
595  stk::mesh::Entity side = mesh.findConnectivityById(element, side_rank, 3);
596 
597  // on the top
598  if(mesh.entityOwnerRank(side)==machRank_)
599  mesh.addEntityToSideset(side,top);
600  }
601 
602  if(nx==0 && (t_offset==4 || t_offset==5)) {
603  stk::mesh::Entity side = mesh.findConnectivityById(element, side_rank, 3);
604 
605  // on the left
606  if(mesh.entityOwnerRank(side)==machRank_)
607  mesh.addEntityToSideset(side,left);
608  }
609  if(nx+1==totalXElems && (t_offset==6 || t_offset==7)) {
610  stk::mesh::Entity side = mesh.findConnectivityById(element, side_rank, 3);
611 
612  // on the right
613  if(mesh.entityOwnerRank(side)==machRank_)
614  mesh.addEntityToSideset(side,right);
615  }
616  }
617 
618  mesh.endModification();
619 }
620 
622 {
623  mesh.beginModification();
624 
625  // get all part vectors
626  stk::mesh::Part * origin = mesh.getNodeset("origin");
627 
629  if(machRank_==0)
630  {
631  // add zero node to origin node set
632  stk::mesh::Entity node = bulkData->get_entity(mesh.getNodeRank(),1);
633  mesh.addEntityToNodeset(node,origin);
634  }
635 
636  mesh.endModification();
637 }
638 
639 // Pre-Condition: call beginModification() before entry
640 // Post-Condition: call endModification() after exit
642 {
645 
646  stk::mesh::Part * edge_block = mesh.getEdgeBlock(edgeBlockName_);
647 
648  stk::mesh::Selector owned_block = metaData->locally_owned_part();
649 
650  std::vector<stk::mesh::Entity> edges;
651  bulkData->get_entities(mesh.getEdgeRank(), owned_block, edges);
652  mesh.addEntitiesToEdgeBlock(edges, edge_block);
653 }
654 
655 // Pre-Condition: call beginModification() before entry
656 // Post-Condition: call endModification() after exit
658 {
661 
662  stk::mesh::Part * face_block = mesh.getFaceBlock(faceBlockName_);
663 
664  stk::mesh::Selector owned_block = metaData->locally_owned_part();
665 
666  std::vector<stk::mesh::Entity> faces;
667  bulkData->get_entities(mesh.getFaceRank(), owned_block, faces);
668  mesh.addEntitiesToFaceBlock(faces, face_block);
669 }
670 
673 {
674  std::size_t i=0,j=0,k=0;
675 
676  k = procRank/(xProcs_*yProcs_); procRank = procRank % (xProcs_*yProcs_);
677  j = procRank/xProcs_; procRank = procRank % xProcs_;
678  i = procRank;
679 
680  return Teuchos::tuple(i,j,k);
681 }
682 
683 } // end panzer_stk
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const
From ParameterListAcceptor.
std::pair< int, int > determineYElemSizeAndStart(int yBlock, unsigned int size, unsigned int rank) const
void addNodeset(const std::string &name)
unsigned entityOwnerRank(stk::mesh::Entity entity) const
void addEntityToNodeset(stk::mesh::Entity entity, stk::mesh::Part *nodeset)
static void parsePeriodicBCList(const Teuchos::RCP< Teuchos::ParameterList > &pl, std::vector< Teuchos::RCP< const PeriodicBC_MatcherBase > > &periodicBC, bool &useBBoxSearch)
basic_FancyOStream & setShowProcRank(const bool showProcRank)
T & get(const std::string &name, T def_value)
Teuchos::RCP< STK_Interface > buildMesh(stk::ParallelMachine parallelMach) const
Build the mesh object.
stk::mesh::Part * getElementBlockPart(const std::string &name) const
get the block part
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 addFaceBlocks(STK_Interface &mesh) const
std::pair< int, int > determineZElemSizeAndStart(int zBlock, unsigned int size, unsigned int rank) 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)
const double * getNodeCoordinates(stk::mesh::EntityId nodeId) const
stk::mesh::Part * getSideset(const std::string &name) const
virtual Teuchos::RCP< STK_Interface > buildUncommitedMesh(stk::ParallelMachine parallelMach) const
void buildElements(stk::ParallelMachine parallelMach, STK_Interface &mesh) const
void addSideSets(STK_Interface &mesh) 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
void buildBlock(stk::ParallelMachine machRank, int xBlock, int yBlock, int zBlock, STK_Interface &mesh) const
Teuchos::RCP< stk::mesh::BulkData > getBulkData() const
stk::mesh::EntityId elementGlobalId(std::size_t lid) const
void setMyParamList(const RCP< ParameterList > &paramList)
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 > &paramList)
From ParameterListAcceptor.
std::pair< int, int > determineXElemSizeAndStart(int xBlock, unsigned int size, unsigned int rank) const
stk::mesh::EntityRank getFaceRank() const
void addEdgeBlocks(STK_Interface &mesh) 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::RCP< stk::mesh::MetaData > getMetaData() const
bool isInitialized() const
Has initialize been called on this mesh object?
Teuchos::Tuple< std::size_t, 3 > procRankToProcTuple(std::size_t procRank) const
what is the 3D tuple describe this processor distribution
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
void buildMetaData(stk::ParallelMachine parallelMach, STK_Interface &mesh) const
stk::mesh::Part * getEdgeBlock(const std::string &name) const
get the block part
void buildTetsOnHex(const Teuchos::Tuple< int, 3 > &meshDesc, const Teuchos::Tuple< int, 3 > &element, stk::mesh::Part *block, const std::vector< stk::mesh::EntityId > &h_nodes, STK_Interface &mesh) const
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)
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="")
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::Tuple< std::size_t, 3 > procTuple_
void addElementBlock(const std::string &name, const CellTopologyData *ctData)
void addNodeSets(STK_Interface &mesh) const