Panzer  Version of the Day
Panzer_STK_LocalMeshUtilities.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 "Panzer_STK_Interface.hpp"
47 
48 #include "Panzer_HashUtils.hpp"
49 #include "Panzer_LocalMeshInfo.hpp"
51 #include "Panzer_FaceToElement.hpp"
52 
53 #include "Panzer_ConnManager.hpp"
54 #include "Panzer_FieldPattern.hpp"
59 
60 #include "Phalanx_KokkosDeviceTypes.hpp"
61 
62 #include "Teuchos_Assert.hpp"
63 
64 #include "Tpetra_Import.hpp"
65 #include "Tpetra_CrsMatrix.hpp"
66 #include "Tpetra_RowMatrixTransposer.hpp"
67 
68 #include <string>
69 #include <map>
70 #include <vector>
71 #include <unordered_set>
72 
73 namespace panzer_stk
74 {
75 
76 // No external access
77 namespace
78 {
79 
83 template <typename LO,typename GO>
84 void
85 buildCellGlobalIDs(panzer::ConnManager<LO,GO> & conn,
86  Kokkos::View<GO*> & globals)
87 {
88  // extract topologies, and build global connectivity...currently assuming only one topology
89  std::vector<shards::CellTopology> elementBlockTopologies;
90  conn.getElementBlockTopologies(elementBlockTopologies);
91  const shards::CellTopology & topology = elementBlockTopologies[0];
92 
93  // FIXME: We assume that all element blocks have the same topology.
94  for(const auto & other_topology : elementBlockTopologies){
95  TEUCHOS_ASSERT(other_topology.getKey() == topology.getKey());
96  }
97 
99  if(topology.getDimension() == 1){
100  cell_pattern = Teuchos::rcp(new panzer::EdgeFieldPattern(topology));
101  } else if(topology.getDimension() == 2){
102  cell_pattern = Teuchos::rcp(new panzer::FaceFieldPattern(topology));
103  } else if(topology.getDimension() == 3){
104  cell_pattern = Teuchos::rcp(new panzer::ElemFieldPattern(topology));
105  }
106 
107 // panzer::EdgeFieldPattern cell_pattern(elementBlockTopologies[0]);
108  conn.buildConnectivity(*cell_pattern);
109 
110  // calculate total number of local cells
111  std::vector<std::string> block_ids;
112  conn.getElementBlockIds(block_ids);
113 
114  std::size_t totalSize = 0;
115  for (std::size_t which_blk=0;which_blk<block_ids.size();which_blk++) {
116  // get the elem to face mapping
117  const std::vector<LO> & localIDs = conn.getElementBlock(block_ids[which_blk]);
118  totalSize += localIDs.size();
119  }
120  globals = Kokkos::View<GO*>("global_cells",totalSize);
121 
122  for (std::size_t id=0;id<totalSize; ++id) {
123  // sanity check
124  int n_conn = conn.getConnectivitySize(id);
125  TEUCHOS_ASSERT(n_conn==1);
126 
127  const GO * connectivity = conn.getConnectivity(id);
128  globals(id) = connectivity[0];
129  }
130 
131 // print_view_1D("buildCellGlobalIDs : globals",globals);
132 }
133 
137 template <typename LO,typename GO>
138 void
139 buildCellToNodes(panzer::ConnManager<LO,GO> & conn, Kokkos::View<GO**> & globals)
140 {
141  // extract topologies, and build global connectivity...currently assuming only one topology
142  std::vector<shards::CellTopology> elementBlockTopologies;
143  conn.getElementBlockTopologies(elementBlockTopologies);
144  const shards::CellTopology & topology = elementBlockTopologies[0];
145 
146  // FIXME: We assume that all element blocks have the same topology.
147  for(const auto & other_topology : elementBlockTopologies){
148  TEUCHOS_ASSERT(other_topology.getKey() == topology.getKey());
149  }
150 
151  panzer::NodalFieldPattern pattern(topology);
152  conn.buildConnectivity(pattern);
153 
154  // calculate total number of local cells
155  std::vector<std::string> block_ids;
156  conn.getElementBlockIds(block_ids);
157 
158  // compute total cells and maximum nodes
159  std::size_t totalCells=0, maxNodes=0;
160  for (std::size_t which_blk=0;which_blk<block_ids.size();which_blk++) {
161  // get the elem to face mapping
162  const std::vector<LO> & localIDs = conn.getElementBlock(block_ids[which_blk]);
163  LO thisSize = conn.getConnectivitySize(localIDs[0]);
164 
165  totalCells += localIDs.size();
166  maxNodes = maxNodes<Teuchos::as<std::size_t>(thisSize) ? Teuchos::as<std::size_t>(thisSize) : maxNodes;
167  }
168  globals = Kokkos::View<GO**>("cell_to_node",totalCells,maxNodes);
169 
170  // build connectivity array
171  for (std::size_t id=0;id<totalCells; ++id) {
172  const GO * connectivity = conn.getConnectivity(id);
173  LO nodeCnt = conn.getConnectivitySize(id);
174 
175  for(int n=0;n<nodeCnt;n++)
176  globals(id,n) = connectivity[n];
177  }
178 
179 // print_view("buildCellToNodes : globals",globals);
180 }
181 
182 template <typename LO,typename GO>
184 buildNodeMap(const Teuchos::RCP<const Teuchos::Comm<int> > & comm,
185  Kokkos::View<const GO**> cells_to_nodes)
186 {
187  using Teuchos::RCP;
188  using Teuchos::rcp;
189 
190  /*
191 
192  This function converts
193 
194  cells_to_nodes(local cell, local node) = global node index
195 
196  to a map describing which global nodes are found on this process
197 
198  Note that we have to ensure that that the global nodes found on this process are unique.
199 
200  */
201 
202  typedef Tpetra::Map<LO,GO> map_type;
203 
204  // get locally unique global ids
205  std::set<GO> global_nodes;
206  for(unsigned int i=0;i<cells_to_nodes.extent(0);i++)
207  for(unsigned int j=0;j<cells_to_nodes.extent(1);j++)
208  global_nodes.insert(cells_to_nodes(i,j));
209 
210  // build local vector contribution
211  Kokkos::View<GO*> node_ids("global_nodes",global_nodes.size());
212  int i = 0;
213  for(auto itr=global_nodes.begin();itr!=global_nodes.end();++itr,++i)
214  node_ids(i) = *itr;
215 
216 // print_view("buildNodeMap : cells_to_nodes",cells_to_nodes);
217 // print_view_1D("buildNodeMap : node_ids",node_ids);
218 
219  return rcp(new map_type(-1,node_ids,0,comm));
220 }
221 
225 template <typename LO,typename GO>
227 buildNodeToCellMatrix(const Teuchos::RCP<const Teuchos::Comm<int> > & comm,
228  Kokkos::View<const GO*> owned_cells,
229  Kokkos::View<const GO**> owned_cells_to_nodes)
230 {
231  using Teuchos::RCP;
232  using Teuchos::rcp;
233 
234  typedef Tpetra::Map<LO,GO> map_type;
235  typedef Tpetra::CrsMatrix<LO,LO,GO> crs_type;
236  typedef Tpetra::Import<LO,GO> import_type;
237 
238  TEUCHOS_ASSERT(owned_cells.extent(0)==owned_cells_to_nodes.extent(0));
239 
240  // build a unque node map to use with fill complete
241 
242  // This map identifies all nodes linked to cells on this process
243  auto node_map = buildNodeMap<LO,GO>(comm,owned_cells_to_nodes);
244 
245  // This map identifies nodes owned by this process
246  auto unique_node_map = Tpetra::createOneToOne(node_map);
247 
248  // This map identifies the cells owned by this process
249  RCP<map_type> cell_map = rcp(new map_type(-1,owned_cells,0,comm));
250 
251  // Create a CRS matrix that stores a pointless value for every global node that belongs to a global cell
252  // This is essentially another way to store cells_to_nodes
253  RCP<crs_type> cell_to_node;
254  {
255 
256  // The matrix is indexed by (global cell, global node) = local node
257  cell_to_node = rcp(new crs_type(cell_map,0));
258  cell_to_node->resumeFill();
259 
260  // fill in the cell to node matrix
261  const unsigned int num_local_cells = owned_cells_to_nodes.extent(0);
262  const unsigned int num_nodes_per_cell = owned_cells_to_nodes.extent(1);
263  std::vector<LO> local_node_indexes(num_nodes_per_cell);
264  std::vector<GO> global_node_indexes(num_nodes_per_cell);
265  for(unsigned int i=0;i<num_local_cells;i++) {
266  const GO global_cell_index = owned_cells(i);
267  // std::vector<LO> vals(cells_to_nodes.extent(1));
268  // std::vector<GO> cols(cells_to_nodes.extent(1));
269  for(unsigned int j=0;j<num_nodes_per_cell;j++) {
270  // vals[j] = Teuchos::as<LO>(j);
271  // cols[j] = cells_to_nodes(i,j);
272  local_node_indexes[j] = Teuchos::as<LO>(j);
273  global_node_indexes[j] = owned_cells_to_nodes(i,j);
274  }
275 
276  // cell_to_node_mat->insertGlobalValues(cells(i),cols,vals);
277  cell_to_node->insertGlobalValues(global_cell_index,global_node_indexes,local_node_indexes);
278  }
279  cell_to_node->fillComplete(unique_node_map,cell_map);
280 
281  }
282 
283  // Transpose the cell_to_node array to create the node_to_cell array
284  RCP<crs_type> node_to_cell;
285  {
286  // Create an object designed to transpose the (global cell, global node) matrix to give
287  // a (global node, global cell) matrix
288  Tpetra::RowMatrixTransposer<LO,LO,GO> transposer(cell_to_node);
289 
290  // Create the transpose crs matrix
291  auto trans = transposer.createTranspose();
292 
293  // We want to import the portion of the transposed matrix relating to all nodes on this process
294  // The importer must import nodes required by this process (node_map) from the unique nodes (nodes living on a process)
295  RCP<import_type> import = rcp(new import_type(unique_node_map,node_map));
296 
297  // Create the crs matrix to store (ghost node, global cell) array
298  // This CRS matrix will have overlapping rows for shared global nodes
299  node_to_cell = rcp(new crs_type(node_map,0));
300 
301  // Transfer data from the transpose array (associated with unique_node_map) to node_to_cell (associated with node_map)
302  node_to_cell->doImport(*trans,*import,Tpetra::INSERT);
303 
304  // Set the fill - basicially locks down the structured of the CRS matrix - required before doing some operations
305  //node_to_cell->fillComplete();
306  node_to_cell->fillComplete(cell_map,unique_node_map);
307  }
308 
309  // Return the node to cell array
310  return node_to_cell;
311 }
312 
315 template <typename GO>
316 Kokkos::View<const GO*>
317 buildGhostedCellOneRing(const Teuchos::RCP<const Teuchos::Comm<int> > & comm,
318  Kokkos::View<const GO*> cells,
319  Kokkos::View<const GO**> cells_to_nodes)
320 {
321  typedef Tpetra::CrsMatrix<int,int,GO> crs_type;
322 
323  // cells : (local cell index) -> global cell index
324  // cells_to_nodes : (local cell index, local node_index) -> global node index
325 
326  /*
327 
328  This function creates a list of global indexes relating to ghost cells
329 
330  It is a little misleading how it does this, but the idea is to use the indexing of a CRS matrix to identify
331  what global cells are connected to which global nodes. The values in the CRS matrix are meaningless, however,
332  the fact that they are filled allows us to ping what index combinations exist.
333 
334  To do this we are going to use cell 'nodes' which could also be cell 'vertices'. It is unclear.
335 
336  First we construct an array that stores that global node indexes make up the connectivity for a given global cell index (order doesn't matter)
337 
338  cell_to_node : cell_to_node[global cell index][global node index] = some value (not important, just has to exist)
339 
340  We are then going to transpose that array
341 
342  node_to_cell : node_to_cell[global node index][global cell index] = some value (not important, just has to exist)
343 
344  The coloring of the node_to_cell array tells us what global cells are connected to a given global node.
345 
346 
347  */
348 
349  // the node to cell matrix: Row = Global Node Id, Cell = Global Cell Id, Value = Cell Local Node Id
350  Teuchos::RCP<crs_type> node_to_cell = buildNodeToCellMatrix<int,GO>(comm,cells,cells_to_nodes);
351 
352  // the set of cells already known
353  std::unordered_set<GO> unique_cells;
354 
355  // mark all owned cells as already known, e.g. and not in the list of
356  // ghstd cells to be constructed
357  for(size_t i=0;i<cells.extent(0);i++) {
358  unique_cells.insert(cells(i));
359  }
360 
361  // The set of ghost cells that share a global node with an owned cell
362  std::set<GO> ghstd_cells_set;
363 
364  // Get a list of global node indexes associated with the cells owned by this process
365 // auto node_map = node_to_cell->getRangeMap()->getMyGlobalIndices();
366  auto node_map = node_to_cell->getMap()->getMyGlobalIndices();
367 
368  // Iterate through the global node indexes associated with this process
369  for(size_t i=0;i<node_map.extent(0);i++) {
370  const GO global_node_index = node_map(i);
371  size_t numEntries = node_to_cell->getNumEntriesInGlobalRow(node_map(i));
372  Teuchos::Array<GO> indices(numEntries);
373  Teuchos::Array<int> values(numEntries);
374 
375  // Copy the row for a global node index into a local vector
376  node_to_cell->getGlobalRowCopy(global_node_index,indices,values,numEntries);
377 
378  for(auto index : indices) {
379  // if this is a new index (not owned, not previously found ghstd index
380  // add it to the list of ghstd cells
381  if(unique_cells.find(index)==unique_cells.end()) {
382  ghstd_cells_set.insert(index);
383  unique_cells.insert(index); // make sure you don't find it again
384  }
385  }
386  }
387 
388  // build an array containing only the ghstd cells
389  int indx = 0;
390  Kokkos::View<GO*> ghstd_cells("ghstd_cells",ghstd_cells_set.size());
391  for(auto global_cell_index : ghstd_cells_set) {
392  ghstd_cells(indx) = global_cell_index;
393  indx++;
394  }
395 
396 // print_view_1D("ghstd_cells",ghstd_cells);
397 
398  return ghstd_cells;
399 }
400 
404 template <typename GO>
405 Kokkos::DynRankView<double,PHX::Device>
406 buildGhostedVertices(const Tpetra::Import<int,GO> & importer,
407  Kokkos::DynRankView<const double,PHX::Device> owned_vertices)
408 {
409  using Teuchos::RCP;
410  using Teuchos::rcp;
411 
412  typedef Tpetra::MultiVector<double,int,GO> mvec_type;
413  typedef typename mvec_type::dual_view_type dual_view_type;
414 
415  size_t owned_cell_cnt = importer.getSourceMap()->getNodeNumElements();
416  size_t ghstd_cell_cnt = importer.getTargetMap()->getNodeNumElements();
417  int vertices_per_cell = owned_vertices.extent(1);
418  int space_dim = owned_vertices.extent(2);
419 
420  TEUCHOS_ASSERT(owned_vertices.extent(0)==owned_cell_cnt);
421 
422  // build vertex multivector
423  RCP<mvec_type> owned_vertices_mv = rcp(new mvec_type(importer.getSourceMap(),vertices_per_cell*space_dim));
424  RCP<mvec_type> ghstd_vertices_mv = rcp(new mvec_type(importer.getTargetMap(),vertices_per_cell*space_dim));
425 
426  {
427  auto owned_vertices_view = owned_vertices_mv->template getLocalView<dual_view_type>();
428  for(size_t i=0;i<owned_cell_cnt;i++) {
429  int l = 0;
430  for(int j=0;j<vertices_per_cell;j++)
431  for(int k=0;k<space_dim;k++,l++)
432  owned_vertices_view(i,l) = owned_vertices(i,j,k);
433  }
434  }
435 
436  // communicate ghstd vertices
437  ghstd_vertices_mv->doImport(*owned_vertices_mv,importer,Tpetra::INSERT);
438 
439  // copy multivector into ghstd vertex structure
440  Kokkos::DynRankView<double,PHX::Device> ghstd_vertices("ghstd_vertices",ghstd_cell_cnt,vertices_per_cell,space_dim);
441  {
442  auto ghstd_vertices_view = ghstd_vertices_mv->template getLocalView<dual_view_type>();
443  for(size_t i=0;i<ghstd_cell_cnt;i++) {
444  int l = 0;
445  for(int j=0;j<vertices_per_cell;j++)
446  for(int k=0;k<space_dim;k++,l++)
447  ghstd_vertices(i,j,k) = ghstd_vertices_view(i,l);
448  }
449  }
450 
451  return ghstd_vertices;
452 } // end build ghstd vertices
453 
454 template<typename LO, typename GO>
455 void
456 setupLocalMeshBlockInfo(const panzer_stk::STK_Interface & mesh,
458  const panzer::LocalMeshInfo<LO,GO> & mesh_info,
459  const std::string & element_block_name,
461 {
462 
463  // This function identifies all cells in mesh_info that belong to element_block_name
464  // and creates a block_info from it.
465 
466  const int num_parent_owned_cells = mesh_info.num_owned_cells;
467 
468  // Make sure connectivity is setup for interfaces between cells
469  {
470  const shards::CellTopology & topology = *(mesh.getCellTopology(element_block_name));
472  if(topology.getDimension() == 1){
473  cell_pattern = Teuchos::rcp(new panzer::EdgeFieldPattern(topology));
474  } else if(topology.getDimension() == 2){
475  cell_pattern = Teuchos::rcp(new panzer::FaceFieldPattern(topology));
476  } else if(topology.getDimension() == 3){
477  cell_pattern = Teuchos::rcp(new panzer::ElemFieldPattern(topology));
478  }
479 
480  conn.buildConnectivity(*cell_pattern);
481  }
482 
483  std::vector<LO> owned_block_cells;
484  for(int parent_owned_cell=0;parent_owned_cell<num_parent_owned_cells;++parent_owned_cell){
485  const LO local_cell = mesh_info.local_cells(parent_owned_cell);
486  const bool is_in_block = conn.getBlockId(local_cell) == element_block_name;
487 
488  if(is_in_block){
489  owned_block_cells.push_back(parent_owned_cell);
490  }
491 
492  }
493 
494  block_info.num_owned_cells = owned_block_cells.size();
495  block_info.element_block_name = element_block_name;
496  block_info.cell_topology = mesh.getCellTopology(element_block_name);
497 
498  panzer::partitioning_utilities::setupSubLocalMeshInfo<LO,GO>(mesh_info, owned_block_cells, block_info);
499 
500 }
501 
502 
503 template<typename LO, typename GO>
504 void
505 setupLocalMeshSidesetInfo(const panzer_stk::STK_Interface & mesh,
506  panzer::ConnManager<LO,GO>& /* conn */,
507  const panzer::LocalMeshInfo<LO,GO> & mesh_info,
508  const std::string & element_block_name,
509  const std::string & sideset_name,
511 {
512 
513  // This function identifies all cells in mesh_info that belong to element_block_name
514  // and creates a block_info from it.
515 
516  const size_t face_subcell_dimension = mesh.getCellTopology(element_block_name)->getDimension()-1;
517 
518  Kokkos::DynRankView<double,PHX::Device> data_allocation;
519 
520  // This is a list of all entities that make up the 'side'
521  std::vector<stk::mesh::Entity> side_entities;
522 
523  // Grab the side entities associated with this sideset on the mesh
524  // Note: Throws exception if element block or sideset doesn't exist
525  try{
526 
527  mesh.getMySides(sideset_name, element_block_name, side_entities);
528 
529  } catch(STK_Interface::SidesetException & e) {
530  std::stringstream ss;
531  std::vector<std::string> sideset_names;
532  mesh.getSidesetNames(sideset_names);
533 
534  // build an error message
535  ss << e.what() << "\nChoose existing sideset:\n";
536  for(const auto & optional_sideset_name : sideset_names){
537  ss << "\t\"" << optional_sideset_name << "\"\n";
538  }
539 
540  TEUCHOS_TEST_FOR_EXCEPTION_PURE_MSG(true,std::logic_error,ss.str());
541 
542  } catch(STK_Interface::ElementBlockException & e) {
543  std::stringstream ss;
544  std::vector<std::string> element_block_names;
545  mesh.getElementBlockNames(element_block_names);
546 
547  // build an error message
548  ss << e.what() << "\nChoose existing element block:\n";
549  for(const auto & optional_element_block_name : element_block_names){
550  ss << "\t\"" << optional_element_block_name << "\"\n";
551  }
552 
553  TEUCHOS_TEST_FOR_EXCEPTION_PURE_MSG(true,std::logic_error,ss.str());
554 
555  } catch(std::logic_error & e) {
556  std::stringstream ss;
557  ss << e.what() << "\nUnrecognized logic error.\n";
558 
559  TEUCHOS_TEST_FOR_EXCEPTION_PURE_MSG(true,std::logic_error,ss.str());
560 
561  }
562 
563  // We now have a list of sideset entities, lets unwrap them and create the sideset_info!
564 
565  // (subcell_dimension, subcell_index) -> list of cells
566  std::map<std::pair<size_t,size_t>, std::vector<size_t> > local_cell_indexes_by_subcell;
567 
568  // List of local subcell indexes indexed by element:
569  // For example: a Tet (element) would have
570  // - 4 triangular faces (subcell_index 0-3, subcell_dimension=2)
571  // - 6 edges (subcell_index 0-5, subcell_dimension=1)
572  // - 4 vertices (subcell_index 0-3, subcell_dimension=0)
573  // Another example: a Line (element) would have
574  // - 2 vertices (subcell_index 0-1, subcell_dimension=0)
575  std::vector<stk::mesh::Entity> elements;
576  std::vector<size_t> subcell_indexes;
577  std::vector<size_t> subcell_dimensions;
578  panzer_stk::workset_utils::getSideElementCascade(mesh, element_block_name, side_entities, subcell_dimensions, subcell_indexes, elements);
579  const LO num_elements = subcell_dimensions.size();
580 
581  // build local cell_ids, mapped by local side id
582  for(LO element_index=0; element_index<num_elements; ++element_index) {
583  const size_t subcell_dimension = subcell_dimensions[element_index];
584  const size_t subcell_index = subcell_indexes[element_index];
585  const size_t element_local_index = mesh.elementLocalId(elements[element_index]);
586 
587  // Add subcell to map
588  local_cell_indexes_by_subcell[std::pair<size_t,size_t>(subcell_dimension, subcell_index)].push_back(element_local_index);
589  }
590 
591  // We now have a list of all local cells 'touching' the side set, as well as the faces, edges, and vertices that define the 'touch'
592 
593  // For our purposes we only really want the cells and faces associated with the side
594 
595  // Each cell can have multiple faces on the sideset
596  std::unordered_set<LO> owned_parent_cells_set;
597  std::map<LO,std::vector<LO> > owned_parent_cell_index_map;
598  {
599  // We use a set to avoid duplicates
600  for(const auto & cell_subcell_pair : local_cell_indexes_by_subcell){
601  const std::pair<size_t,size_t> & subcell_definition = cell_subcell_pair.first;
602  const LO subcell_index = LO(subcell_definition.second);
603  if(subcell_definition.first == face_subcell_dimension){
604  const std::vector<size_t> & local_cell_indexes_for_subcell = cell_subcell_pair.second;
605  for(const size_t & local_cell_index : local_cell_indexes_for_subcell){
606 
607  // TODO: Super slow - fix this!!
608  for(LO i=0;i<LO(mesh_info.local_cells.extent(0)); ++i){
609  if(mesh_info.local_cells(i) == LO(local_cell_index)){
610  owned_parent_cell_index_map[i].push_back(subcell_index);
611  }
612  }
613  }
614  }
615  }
616 
617  for(const auto & pair : owned_parent_cell_index_map){
618  owned_parent_cells_set.insert(pair.first);
619  }
620 
621  }
622 
623  const LO num_owned_cells = owned_parent_cell_index_map.size();
624 
625  sideset_info.element_block_name = element_block_name;
626  sideset_info.sideset_name = sideset_name;
627  sideset_info.cell_topology = mesh.getCellTopology(element_block_name);
628 
629  sideset_info.num_owned_cells = num_owned_cells;
630 
631  struct face_t{
632  face_t(LO c0, LO c1, LO sc0, LO sc1)
633  {
634  cell_0=c0;
635  cell_1=c1;
636  subcell_index_0=sc0;
637  subcell_index_1=sc1;
638  }
639  LO cell_0;
640  LO cell_1;
641  LO subcell_index_0;
642  LO subcell_index_1;
643  };
644 
645 
646  // Figure out how many cells on the other side of the sideset are ghost or virtual
647  std::unordered_set<LO> ghstd_parent_cells_set, virtual_parent_cells_set;
648  std::vector<face_t> faces;
649  {
650  LO parent_virtual_cell_offset = mesh_info.num_owned_cells + mesh_info.num_ghstd_cells;
651  for(const auto & local_cell_index_pair : owned_parent_cell_index_map){
652  const LO local_cell = local_cell_index_pair.first;
653  const std::vector<LO> & subcell_indexes = local_cell_index_pair.second;
654 
655  for(const LO & subcell_index : subcell_indexes){
656 
657  const LO face = mesh_info.cell_to_faces(local_cell, subcell_index);
658  const LO face_other_side = (mesh_info.face_to_cells(face,0) == local_cell) ? 1 : 0;
659 
660  TEUCHOS_ASSERT(subcell_index == mesh_info.face_to_lidx(face, 1-face_other_side));
661 
662  const LO other_side_cell = mesh_info.face_to_cells(face, face_other_side);
663  const LO other_side_subcell_index = mesh_info.face_to_lidx(face, face_other_side);
664 
665  faces.push_back(face_t(local_cell, other_side_cell, subcell_index, other_side_subcell_index));
666 
667  if(other_side_cell >= parent_virtual_cell_offset){
668  virtual_parent_cells_set.insert(other_side_cell);
669  } else {
670  ghstd_parent_cells_set.insert(other_side_cell);
671  }
672  }
673  }
674  }
675 
676  std::vector<LO> all_cells;
677  all_cells.insert(all_cells.end(),owned_parent_cells_set.begin(),owned_parent_cells_set.end());
678  all_cells.insert(all_cells.end(),ghstd_parent_cells_set.begin(),ghstd_parent_cells_set.end());
679  all_cells.insert(all_cells.end(),virtual_parent_cells_set.begin(),virtual_parent_cells_set.end());
680 
681  sideset_info.num_ghstd_cells = ghstd_parent_cells_set.size();
682  sideset_info.num_virtual_cells = virtual_parent_cells_set.size();
683 
684  const LO num_real_cells = sideset_info.num_owned_cells + sideset_info.num_ghstd_cells;
685  const LO num_total_cells = num_real_cells + sideset_info.num_virtual_cells;
686  const LO num_vertices_per_cell = mesh_info.cell_vertices.extent(1);
687  const LO num_dims = mesh_info.cell_vertices.extent(2);
688 
689  sideset_info.global_cells = Kokkos::View<GO*>("global_cells", num_total_cells);
690  sideset_info.local_cells = Kokkos::View<LO*>("local_cells", num_total_cells);
691  sideset_info.cell_vertices = Kokkos::View<double***,PHX::Device>("cell_vertices", num_total_cells, num_vertices_per_cell, num_dims);
692  Kokkos::deep_copy(sideset_info.cell_vertices,0.);
693 
694  for(LO i=0; i<num_total_cells; ++i){
695  const LO parent_cell = all_cells[i];
696  sideset_info.local_cells(i) = mesh_info.local_cells(parent_cell);
697  sideset_info.global_cells(i) = mesh_info.global_cells(parent_cell);
698  for(LO j=0; j<num_vertices_per_cell; ++j){
699  for(LO k=0; k<num_dims; ++k){
700  sideset_info.cell_vertices(i,j,k) = mesh_info.cell_vertices(parent_cell,j,k);
701  }
702  }
703  }
704 
705  // Now we have to set the connectivity for the faces.
706 
707  const LO num_faces = faces.size();
708  const LO num_faces_per_cell = mesh_info.cell_to_faces.extent(1);
709 
710  sideset_info.face_to_cells = Kokkos::View<LO*[2]>("face_to_cells", num_faces);
711  sideset_info.face_to_lidx = Kokkos::View<LO*[2]>("face_to_lidx", num_faces);
712  sideset_info.cell_to_faces = Kokkos::View<LO**>("cell_to_faces", num_total_cells, num_faces_per_cell);
713 
714  // Default the system with invalid cell index - this will be most of the entries
715  Kokkos::deep_copy(sideset_info.cell_to_faces, -1);
716 
717  for(int face_index=0;face_index<num_faces;++face_index){
718  const face_t & face = faces[face_index];
719  const LO & cell_0 = face.cell_0;
720  const LO & cell_1 = face.cell_1;
721 
722  // TODO: Super expensive... need to find a better way
723  const LO sideset_cell_0 = std::distance(all_cells.begin(), std::find(all_cells.begin(), all_cells.end(), cell_0));
724  const LO sideset_cell_1 = std::distance(all_cells.begin(), std::find(all_cells.begin()+num_owned_cells, all_cells.end(), cell_1));
725 
726  sideset_info.face_to_cells(face_index,0) = sideset_cell_0;
727  sideset_info.face_to_cells(face_index,1) = sideset_cell_1;
728 
729  sideset_info.face_to_lidx(face_index,0) = face.subcell_index_0;
730  sideset_info.face_to_lidx(face_index,1) = face.subcell_index_1;
731 
732  sideset_info.cell_to_faces(sideset_cell_0,face.subcell_index_0) = face_index;
733  sideset_info.cell_to_faces(sideset_cell_1,face.subcell_index_1) = face_index;
734 
735  }
736 
737 }
738 
739 }
740 
741 template <typename LO, typename GO>
742 void
744  panzer::LocalMeshInfo<LO,GO> & mesh_info)
745 {
746  using Teuchos::RCP;
747  using Teuchos::rcp;
748 
749  //typedef Tpetra::CrsMatrix<int,LO,GO> crs_type;
750  typedef Tpetra::Map<LO,GO> map_type;
751  typedef Tpetra::Import<LO,GO> import_type;
752  //typedef Tpetra::MultiVector<double,LO,GO> mvec_type;
753  //typedef Tpetra::MultiVector<GO,LO,GO> ordmvec_type;
754 
755  // Make sure the STK interface is valid
757 
758  // This is required by some of the STK stuff
759  TEUCHOS_ASSERT(typeid(LO) == typeid(int));
760 
762 
763  TEUCHOS_FUNC_TIME_MONITOR("panzer_stk::generateLocalMeshInfo");
764 
765  // This horrible line of code is required since the connection manager only takes rcps of a mesh
766  RCP<const panzer_stk::STK_Interface> mesh_rcp = Teuchos::rcpFromRef(mesh);
767  // We're allowed to do this since the connection manager only exists in this scope... even though it is also an RCP...
768 
769  // extract topology handle
771  panzer::ConnManager<LO,GO> & conn = *conn_rcp;
772 
773  // build cell to node map
774  Kokkos::View<GO**> owned_cell_to_nodes;
775  buildCellToNodes(conn, owned_cell_to_nodes);
776 
777  // build the local to global cell ID map
779  Kokkos::View<GO*> owned_cells;
780  buildCellGlobalIDs(conn, owned_cells);
781 
782  // get neighboring cells
784  Kokkos::View<const GO*> ghstd_cells = buildGhostedCellOneRing<GO>(comm,owned_cells,owned_cell_to_nodes);
785 
786  // build cell maps
788 
789  RCP<map_type> owned_cell_map = rcp(new map_type(-1,owned_cells, 0,comm));
790  RCP<map_type> ghstd_cell_map = rcp(new map_type(-1,ghstd_cells,0,comm));
791 
792  // build importer: cell importer, owned to ghstd
793  RCP<import_type> cellimport_own2ghst = rcp(new import_type(owned_cell_map,ghstd_cell_map));
794 
795  // read all the vertices associated with these elements, get ghstd contributions
797  std::vector<std::size_t> localCells(owned_cells.extent(0),-1);
798  for(size_t i=0;i<localCells.size();i++){
799  localCells[i] = i;
800  }
801 
802  // TODO: This all needs to be rewritten for when element blocks have different cell topologies
803  std::vector<std::string> element_block_names;
804  mesh.getElementBlockNames(element_block_names);
805 
806  const shards::CellTopology & cell_topology = *(mesh.getCellTopology(element_block_names[0]));
807 
808  const int space_dim = cell_topology.getDimension();
809  const int vertices_per_cell = cell_topology.getVertexCount();
810  const int faces_per_cell = cell_topology.getSubcellCount(space_dim-1);
811 
812  // Kokkos::View<double***> owned_vertices("owned_vertices",localCells.size(),vertices_per_cell,space_dim);
813  Kokkos::DynRankView<double,PHX::Device> owned_vertices("owned_vertices",localCells.size(),vertices_per_cell,space_dim);
814  mesh.getElementVerticesNoResize(localCells,owned_vertices);
815 
816  // this builds a ghstd vertex array
817  Kokkos::DynRankView<double,PHX::Device> ghstd_vertices = buildGhostedVertices(*cellimport_own2ghst,owned_vertices);
818 
819  // build edge to cell neighbor mapping
821 
822  std::unordered_map<GO,int> global_to_local;
823  global_to_local[-1] = -1; // this is the "no neighbor" flag
824  for(size_t i=0;i<owned_cells.extent(0);i++)
825  global_to_local[owned_cells(i)] = i;
826  for(size_t i=0;i<ghstd_cells.extent(0);i++)
827  global_to_local[ghstd_cells(i)] = i+Teuchos::as<int>(owned_cells.extent(0));
828 
829  // this class comes from Mini-PIC and Matt B
831  faceToElement->initialize(conn);
832  auto elems_by_face = faceToElement->getFaceToElementsMap();
833  auto face_to_lidx = faceToElement->getFaceToCellLocalIdxMap();
834 
835  const int num_owned_cells =owned_cells.extent(0);
836  const int num_ghstd_cells =ghstd_cells.extent(0);
837 
838 // print_view("elems_by_face",elems_by_face);
839 
840  // We also need to consider faces that connect to cells that do not exist, but are needed for boundary conditions
841  // We dub them virtual cell since there should be no geometry associated with them, or topology really
842  // They exist only for datastorage so that they are consistent with 'real' cells from an algorithm perspective
843 
844  // Each virtual face (face linked to a '-1' cell) requires a virtual cell (i.e. turn the '-1' into a virtual cell)
845  // Virtual cells are those that do not exist but are connected to an owned cell
846  // Note - in the future, ghosted cells will also need to connect to virtual cells at boundary conditions, but for the moment we will ignore this.
847 
848  // Iterate over all faces and identify the faces connected to a potential virtual cell
849  std::vector<int> all_boundary_faces;
850  const int num_faces = elems_by_face.extent(0);
851  for(int face=0;face<num_faces;++face){
852  if(elems_by_face(face,0) < 0 or elems_by_face(face,1) < 0){
853  all_boundary_faces.push_back(face);
854  }
855  }
856  const LO num_virtual_cells = all_boundary_faces.size();
857 
858  // total cells and faces include owned, ghosted, and virtual
859  const LO num_real_cells = num_owned_cells + num_ghstd_cells;
860  const LO num_total_cells = num_real_cells + num_virtual_cells;
861  const LO num_total_faces = elems_by_face.extent(0);
862 
863  // Create some global indexes associated with the virtual cells
864  // Note: We are assuming that virtual cells belong to ranks and are not 'shared' - this will change later on
865  Kokkos::View<GO*> virtual_cells = Kokkos::View<GO*>("virtual_cells",num_virtual_cells);
866  {
867  const int num_ranks = comm->getSize();
868  const int rank = comm->getRank();
869 
870  std::vector<GO> owned_cell_distribution(num_ranks,0);
871  {
872  std::vector<GO> my_owned_cell_distribution(num_ranks,0);
873  my_owned_cell_distribution[rank] = num_owned_cells;
874 
875  Teuchos::reduceAll(*comm,Teuchos::REDUCE_SUM, num_ranks, my_owned_cell_distribution.data(),owned_cell_distribution.data());
876  }
877 
878  std::vector<GO> virtual_cell_distribution(num_ranks,0);
879  {
880  std::vector<GO> my_virtual_cell_distribution(num_ranks,0);
881  my_virtual_cell_distribution[rank] = num_virtual_cells;
882 
883  Teuchos::reduceAll(*comm,Teuchos::REDUCE_SUM, num_ranks, my_virtual_cell_distribution.data(),virtual_cell_distribution.data());
884  }
885 
886  GO num_global_real_cells=0;
887  for(int i=0;i<num_ranks;++i){
888  num_global_real_cells+=owned_cell_distribution[i];
889  }
890 
891  GO global_virtual_start_idx = num_global_real_cells;
892  for(int i=0;i<rank;++i){
893  global_virtual_start_idx += virtual_cell_distribution[i];
894  }
895 
896  for(int i=0;i<num_virtual_cells;++i){
897  virtual_cells(i) = global_virtual_start_idx + GO(i);
898  }
899 
900  }
901 
902  // Lookup cells connected to a face
903  Kokkos::View<LO*[2]> face_to_cells = Kokkos::View<LO*[2]>("face_to_cells",num_total_faces);
904 
905  // Lookup local face indexes given cell and left/right state (0/1)
906  Kokkos::View<LO*[2]> face_to_localidx = Kokkos::View<LO*[2]>("face_to_localidx",num_total_faces);
907 
908  // Lookup face index given a cell and local face index
909  Kokkos::View<LO**> cell_to_face = Kokkos::View<LO**>("cell_to_face",num_total_cells,faces_per_cell);
910 
911  // initialize with negative one cells that are not associated with a face
912  Kokkos::deep_copy(cell_to_face,-1);
913 
914  // Transfer information from 'faceToElement' datasets to local arrays
915  {
916  int virtual_cell_index = num_real_cells;
917  for(size_t f=0;f<elems_by_face.extent(0);f++) {
918 
919  const GO global_c0 = elems_by_face(f,0);
920  const GO global_c1 = elems_by_face(f,1);
921 
922  // make sure that no bonus cells get in here
923  TEUCHOS_ASSERT(global_to_local.find(global_c0)!=global_to_local.end());
924  TEUCHOS_ASSERT(global_to_local.find(global_c1)!=global_to_local.end());
925 
926  auto c0 = global_to_local[global_c0];
927  auto lidx0 = face_to_lidx(f,0);
928  auto c1 = global_to_local[global_c1];
929  auto lidx1 = face_to_lidx(f,1);
930 
931  // Test for virtual cells
932 
933  // Left cell
934  if(c0 < 0){
935  // Virtual cell - create it!
936  c0 = virtual_cell_index++;
937 
938  // We need the subcell_index to line up between real and virtual cell
939  // This way the face has the same geometry... though the face normal
940  // will point in the wrong direction
941  lidx0 = lidx1;
942  }
943  cell_to_face(c0,lidx0) = f;
944 
945 
946  // Right cell
947  if(c1 < 0){
948  // Virtual cell - create it!
949  c1 = virtual_cell_index++;
950 
951  // We need the subcell_index to line up between real and virtual cell
952  // This way the face has the same geometry... though the face normal
953  // will point in the wrong direction
954  lidx1 = lidx0;
955  }
956  cell_to_face(c1,lidx1) = f;
957 
958  // Faces point from low cell index to high cell index
959  if(c0<c1){
960  face_to_cells(f,0) = c0;
961  face_to_localidx(f,0) = lidx0;
962  face_to_cells(f,1) = c1;
963  face_to_localidx(f,1) = lidx1;
964  } else {
965  face_to_cells(f,0) = c1;
966  face_to_localidx(f,0) = lidx1;
967  face_to_cells(f,1) = c0;
968  face_to_localidx(f,1) = lidx0;
969  }
970 
971  // We should avoid having two virtual cells linked together.
972  TEUCHOS_ASSERT(c0<num_real_cells or c1<num_real_cells)
973 
974  }
975  }
976 
977  // at this point all the data structures have been built, so now we can "do" DG.
978  // There are two of everything, an "owned" data structured corresponding to "owned"
979  // cells. And a "ghstd" data structure corresponding to ghosted cells
981 
982  mesh_info.cell_to_faces = cell_to_face;
983  mesh_info.face_to_cells = face_to_cells; // faces
984  mesh_info.face_to_lidx = face_to_localidx;
985 
986  mesh_info.num_owned_cells = owned_cells.extent(0);
987  mesh_info.num_ghstd_cells = ghstd_cells.extent(0);
988  mesh_info.num_virtual_cells = virtual_cells.extent(0);
989 
990  mesh_info.global_cells = Kokkos::View<GO*>("global_cell_indices",num_total_cells);
991  mesh_info.local_cells = Kokkos::View<LO*>("local_cell_indices",num_total_cells);
992 
993  for(int i=0;i<num_owned_cells;++i){
994  mesh_info.global_cells(i) = owned_cells(i);
995  mesh_info.local_cells(i) = i;
996  }
997 
998  for(int i=0;i<num_ghstd_cells;++i){
999  mesh_info.global_cells(i+num_owned_cells) = ghstd_cells(i);
1000  mesh_info.local_cells(i+num_owned_cells) = i+num_owned_cells;
1001  }
1002 
1003  for(int i=0;i<num_virtual_cells;++i){
1004  mesh_info.global_cells(i+num_real_cells) = virtual_cells(i);
1005  mesh_info.local_cells(i+num_real_cells) = i+num_real_cells;
1006  }
1007 
1008  mesh_info.cell_vertices = Kokkos::View<double***, PHX::Device>("cell_vertices",num_total_cells,vertices_per_cell,space_dim);
1009 
1010  // Initialize coordinates to zero
1011  Kokkos::deep_copy(mesh_info.cell_vertices, 0.);
1012 
1013  for(int i=0;i<num_owned_cells;++i){
1014  for(int j=0;j<vertices_per_cell;++j){
1015  for(int k=0;k<space_dim;++k){
1016  mesh_info.cell_vertices(i,j,k) = owned_vertices(i,j,k);
1017  }
1018  }
1019  }
1020 
1021  for(int i=0;i<num_ghstd_cells;++i){
1022  for(int j=0;j<vertices_per_cell;++j){
1023  for(int k=0;k<space_dim;++k){
1024  mesh_info.cell_vertices(i+num_owned_cells,j,k) = ghstd_vertices(i,j,k);
1025  }
1026  }
1027  }
1028 
1029  // This will backfire at some point, but we're going to make the virtual cell have the same geometry as the cell it interfaces with
1030  // This way we can define a virtual cell geometry without extruding the face outside of the domain
1031  for(int i=0;i<num_virtual_cells;++i){
1032 
1033  const LO virtual_cell = i+num_real_cells;
1034  bool exists = false;
1035  for(int local_face=0; local_face<faces_per_cell; ++local_face){
1036  const LO face = cell_to_face(virtual_cell, local_face);
1037  if(face >= 0){
1038  exists = true;
1039  const LO other_side = (face_to_cells(face, 0) == virtual_cell) ? 1 : 0;
1040  const LO real_cell = face_to_cells(face,other_side);
1041  TEUCHOS_ASSERT(real_cell < num_real_cells);
1042  for(int j=0;j<vertices_per_cell;++j){
1043  for(int k=0;k<space_dim;++k){
1044  mesh_info.cell_vertices(virtual_cell,j,k) = mesh_info.cell_vertices(real_cell,j,k);
1045  }
1046  }
1047  break;
1048  }
1049  }
1050 
1051  TEUCHOS_TEST_FOR_EXCEPT_MSG(!exists, "panzer_stk::generateLocalMeshInfo : Virtual cell is not linked to real cell");
1052  }
1053 
1054  // Setup element blocks and sidesets
1055  std::vector<std::string> sideset_names;
1056  mesh.getSidesetNames(sideset_names);
1057 
1058  for(const std::string & element_block_name : element_block_names){
1059  panzer::LocalMeshBlockInfo<LO,GO> & block_info = mesh_info.element_blocks[element_block_name];
1060  setupLocalMeshBlockInfo(mesh, conn, mesh_info, element_block_name, block_info);
1061 
1062  // Setup sidesets
1063  for(const std::string & sideset_name : sideset_names){
1064  panzer::LocalMeshSidesetInfo<LO,GO> & sideset_info = mesh_info.sidesets[element_block_name][sideset_name];
1065  setupLocalMeshSidesetInfo(mesh, conn, mesh_info, element_block_name, sideset_name, sideset_info);
1066  }
1067 
1068  }
1069 
1070 }
1071 
1072 }
1073 
1074 // Explicit instantiation
1075 template
1076 void
1077 panzer_stk::generateLocalMeshInfo<int,int>(const panzer_stk::STK_Interface & mesh,
1078  panzer::LocalMeshInfo<int,int> & mesh_info);
1079 
1080 #ifndef PANZER_ORDINAL64_IS_INT
1081 template
1082 void
1083 panzer_stk::generateLocalMeshInfo<int,panzer::Ordinal64>(const panzer_stk::STK_Interface & mesh,
1085 #endif
Teuchos::RCP< const shards::CellTopology > cell_topology
#define TEUCHOS_FUNC_TIME_MONITOR(FUNCNAME)
void getMySides(const std::string &sideName, std::vector< stk::mesh::Entity > &sides) const
virtual void getElementBlockTopologies(std::vector< shards::CellTopology > &elementBlockTopologies) const =0
virtual const GlobalOrdinal * getConnectivity(LocalOrdinal localElmtId) const =0
std::size_t elementLocalId(stk::mesh::Entity elmt) const
virtual const std::vector< LocalOrdinal > & getElementBlock(const std::string &blockID) const =0
std::map< std::string, LocalMeshBlockInfo< LO, GO > > element_blocks
Kokkos::View< GO * > global_cells
Teuchos::RCP< const shards::CellTopology > cell_topology
virtual void buildConnectivity(const FieldPattern &fp)=0
virtual LocalOrdinal getConnectivitySize(LocalOrdinal localElmtId) const =0
virtual void getElementBlockIds(std::vector< std::string > &elementBlockIds) const =0
void getElementBlockNames(std::vector< std::string > &names) const
virtual std::string getBlockId(LocalOrdinal localElmtId) const =0
std::map< std::string, std::map< std::string, LocalMeshSidesetInfo< LO, GO > > > sidesets
Kokkos::View< LO * > local_cells
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
#define TEUCHOS_TEST_FOR_EXCEPT_MSG(throw_exception_test, msg)
bool isInitialized() const
Has initialize been called on this mesh object?
TEUCHOS_DEPRECATED void reduceAll(const Comm< Ordinal > &comm, const EReductionType reductType, const Packet &send, Packet *globalReduct)
void getElementVerticesNoResize(const std::vector< std::size_t > &localIds, ArrayT &vertices) const
Kokkos::View< LO *[2]> face_to_cells
Kokkos::View< double ***, PHX::Device > cell_vertices
Teuchos::RCP< const Teuchos::Comm< int > > getComm() const
get the comm associated with this mesh
Kokkos::View< LO ** > cell_to_faces
Teuchos::RCP< const shards::CellTopology > getCellTopology(const std::string &eBlock) const
void getSidesetNames(std::vector< std::string > &name) const
Kokkos::View< LO *[2]> face_to_lidx
void generateLocalMeshInfo(const panzer_stk::STK_Interface &mesh, panzer::LocalMeshInfo< LO, GO > &mesh_info)
#define TEUCHOS_ASSERT(assertion_test)
void getSideElementCascade(const panzer_stk::STK_Interface &mesh, const std::string &blockId, const std::vector< stk::mesh::Entity > &sides, std::vector< std::size_t > &localSubcellDim, std::vector< std::size_t > &localSubcellIds, std::vector< stk::mesh::Entity > &elements)
#define TEUCHOS_TEST_FOR_EXCEPTION_PURE_MSG(throw_exception_test, Exception, msg)