60 #include "Phalanx_KokkosDeviceTypes.hpp" 62 #include "Teuchos_Assert.hpp" 64 #include "Tpetra_Import.hpp" 65 #include "Tpetra_CrsMatrix.hpp" 66 #include "Tpetra_RowMatrixTransposer.hpp" 71 #include <unordered_set> 83 template <
typename LO,
typename GO>
86 Kokkos::View<GO*> & globals)
89 std::vector<shards::CellTopology> elementBlockTopologies;
91 const shards::CellTopology & topology = elementBlockTopologies[0];
94 for(
const auto & other_topology : elementBlockTopologies){
99 if(topology.getDimension() == 1){
101 }
else if(topology.getDimension() == 2){
103 }
else if(topology.getDimension() == 3){
111 std::vector<std::string> block_ids;
114 std::size_t totalSize = 0;
115 for (std::size_t which_blk=0;which_blk<block_ids.size();which_blk++) {
117 const std::vector<LO> & localIDs = conn.
getElementBlock(block_ids[which_blk]);
118 totalSize += localIDs.size();
120 globals = Kokkos::View<GO*>(
"global_cells",totalSize);
122 for (std::size_t
id=0;
id<totalSize; ++id) {
128 globals(
id) = connectivity[0];
137 template <
typename LO,
typename GO>
142 std::vector<shards::CellTopology> elementBlockTopologies;
144 const shards::CellTopology & topology = elementBlockTopologies[0];
147 for(
const auto & other_topology : elementBlockTopologies){
155 std::vector<std::string> block_ids;
159 std::size_t totalCells=0, maxNodes=0;
160 for (std::size_t which_blk=0;which_blk<block_ids.size();which_blk++) {
162 const std::vector<LO> & localIDs = conn.
getElementBlock(block_ids[which_blk]);
165 totalCells += localIDs.size();
166 maxNodes = maxNodes<Teuchos::as<std::size_t>(thisSize) ? Teuchos::as<std::size_t>(thisSize) : maxNodes;
168 globals = Kokkos::View<GO**>(
"cell_to_node",totalCells,maxNodes);
171 for (std::size_t
id=0;
id<totalCells; ++id) {
175 for(
int n=0;n<nodeCnt;n++)
176 globals(
id,n) = connectivity[n];
182 template <
typename LO,
typename GO>
185 Kokkos::View<const GO**> cells_to_nodes)
202 typedef Tpetra::Map<LO,GO> map_type;
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));
211 Kokkos::View<GO*> node_ids(
"global_nodes",global_nodes.size());
213 for(
auto itr=global_nodes.begin();itr!=global_nodes.end();++itr,++i)
219 return rcp(
new map_type(-1,node_ids,0,comm));
225 template <
typename LO,
typename GO>
228 Kokkos::View<const GO*> owned_cells,
229 Kokkos::View<const GO**> owned_cells_to_nodes)
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;
238 TEUCHOS_ASSERT(owned_cells.extent(0)==owned_cells_to_nodes.extent(0));
243 auto node_map = buildNodeMap<LO,GO>(comm,owned_cells_to_nodes);
246 auto unique_node_map = Tpetra::createOneToOne(node_map);
249 RCP<map_type> cell_map =
rcp(
new map_type(-1,owned_cells,0,comm));
253 RCP<crs_type> cell_to_node;
257 cell_to_node =
rcp(
new crs_type(cell_map,0));
258 cell_to_node->resumeFill();
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);
269 for(
unsigned int j=0;j<num_nodes_per_cell;j++) {
272 local_node_indexes[j] = Teuchos::as<LO>(j);
273 global_node_indexes[j] = owned_cells_to_nodes(i,j);
277 cell_to_node->insertGlobalValues(global_cell_index,global_node_indexes,local_node_indexes);
279 cell_to_node->fillComplete(unique_node_map,cell_map);
284 RCP<crs_type> node_to_cell;
288 Tpetra::RowMatrixTransposer<LO,LO,GO> transposer(cell_to_node);
291 auto trans = transposer.createTranspose();
295 RCP<import_type>
import =
rcp(
new import_type(unique_node_map,node_map));
299 node_to_cell =
rcp(
new crs_type(node_map,0));
302 node_to_cell->doImport(*trans,*
import,Tpetra::INSERT);
306 node_to_cell->fillComplete(cell_map,unique_node_map);
315 template <
typename GO>
316 Kokkos::View<const GO*>
318 Kokkos::View<const GO*> cells,
319 Kokkos::View<const GO**> cells_to_nodes)
321 typedef Tpetra::CrsMatrix<int,int,GO> crs_type;
353 std::unordered_set<GO> unique_cells;
357 for(
size_t i=0;i<cells.extent(0);i++) {
358 unique_cells.insert(cells(i));
362 std::set<GO> ghstd_cells_set;
366 auto node_map = node_to_cell->getMap()->getMyGlobalIndices();
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));
376 node_to_cell->getGlobalRowCopy(global_node_index,indices,values,numEntries);
378 for(
auto index : indices) {
381 if(unique_cells.find(index)==unique_cells.end()) {
382 ghstd_cells_set.insert(index);
383 unique_cells.insert(index);
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;
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)
412 typedef Tpetra::MultiVector<double,int,GO> mvec_type;
413 typedef typename mvec_type::dual_view_type dual_view_type;
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);
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));
427 auto owned_vertices_view = owned_vertices_mv->template getLocalView<dual_view_type>();
428 for(
size_t i=0;i<owned_cell_cnt;i++) {
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);
437 ghstd_vertices_mv->doImport(*owned_vertices_mv,importer,Tpetra::INSERT);
440 Kokkos::DynRankView<double,PHX::Device> ghstd_vertices(
"ghstd_vertices",ghstd_cell_cnt,vertices_per_cell,space_dim);
442 auto ghstd_vertices_view = ghstd_vertices_mv->template getLocalView<dual_view_type>();
443 for(
size_t i=0;i<ghstd_cell_cnt;i++) {
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);
451 return ghstd_vertices;
454 template<
typename LO,
typename GO>
459 const std::string & element_block_name,
470 const shards::CellTopology & topology = *(mesh.
getCellTopology(element_block_name));
472 if(topology.getDimension() == 1){
474 }
else if(topology.getDimension() == 2){
476 }
else if(topology.getDimension() == 3){
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;
489 owned_block_cells.push_back(parent_owned_cell);
498 panzer::partitioning_utilities::setupSubLocalMeshInfo<LO,GO>(mesh_info, owned_block_cells, block_info);
503 template<
typename LO,
typename GO>
508 const std::string & element_block_name,
509 const std::string & sideset_name,
516 const size_t face_subcell_dimension = mesh.
getCellTopology(element_block_name)->getDimension()-1;
518 Kokkos::DynRankView<double,PHX::Device> data_allocation;
521 std::vector<stk::mesh::Entity> side_entities;
527 mesh.
getMySides(sideset_name, element_block_name, side_entities);
529 }
catch(STK_Interface::SidesetException & e) {
530 std::stringstream ss;
531 std::vector<std::string> sideset_names;
535 ss << e.what() <<
"\nChoose existing sideset:\n";
536 for(
const auto & optional_sideset_name : sideset_names){
537 ss <<
"\t\"" << optional_sideset_name <<
"\"\n";
542 }
catch(STK_Interface::ElementBlockException & e) {
543 std::stringstream ss;
544 std::vector<std::string> element_block_names;
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";
555 }
catch(std::logic_error & e) {
556 std::stringstream ss;
557 ss << e.what() <<
"\nUnrecognized logic error.\n";
566 std::map<std::pair<size_t,size_t>, std::vector<size_t> > local_cell_indexes_by_subcell;
575 std::vector<stk::mesh::Entity> elements;
576 std::vector<size_t> subcell_indexes;
577 std::vector<size_t> subcell_dimensions;
579 const LO num_elements = subcell_dimensions.size();
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]);
588 local_cell_indexes_by_subcell[std::pair<size_t,size_t>(subcell_dimension, subcell_index)].push_back(element_local_index);
596 std::unordered_set<LO> owned_parent_cells_set;
597 std::map<LO,std::vector<LO> > owned_parent_cell_index_map;
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){
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);
617 for(
const auto & pair : owned_parent_cell_index_map){
618 owned_parent_cells_set.insert(pair.first);
623 const LO num_owned_cells = owned_parent_cell_index_map.size();
632 face_t(LO c0, LO c1, LO sc0, LO sc1)
647 std::unordered_set<LO> ghstd_parent_cells_set, virtual_parent_cells_set;
648 std::vector<face_t> faces;
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;
655 for(
const LO & subcell_index : subcell_indexes){
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;
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);
665 faces.push_back(face_t(local_cell, other_side_cell, subcell_index, other_side_subcell_index));
667 if(other_side_cell >= parent_virtual_cell_offset){
668 virtual_parent_cells_set.insert(other_side_cell);
670 ghstd_parent_cells_set.insert(other_side_cell);
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());
686 const LO num_vertices_per_cell = mesh_info.
cell_vertices.extent(1);
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);
694 for(LO i=0; i<num_total_cells; ++i){
695 const LO parent_cell = all_cells[i];
698 for(LO j=0; j<num_vertices_per_cell; ++j){
699 for(LO k=0; k<num_dims; ++k){
707 const LO num_faces = faces.size();
708 const LO num_faces_per_cell = mesh_info.
cell_to_faces.extent(1);
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);
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;
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));
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;
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;
741 template <
typename LO,
typename GO>
750 typedef Tpetra::Map<LO,GO> map_type;
751 typedef Tpetra::Import<LO,GO> import_type;
774 Kokkos::View<GO**> owned_cell_to_nodes;
775 buildCellToNodes(conn, owned_cell_to_nodes);
779 Kokkos::View<GO*> owned_cells;
780 buildCellGlobalIDs(conn, owned_cells);
784 Kokkos::View<const GO*> ghstd_cells = buildGhostedCellOneRing<GO>(comm,owned_cells,owned_cell_to_nodes);
793 RCP<import_type> cellimport_own2ghst =
rcp(
new import_type(owned_cell_map,ghstd_cell_map));
797 std::vector<std::size_t> localCells(owned_cells.extent(0),-1);
798 for(
size_t i=0;i<localCells.size();i++){
803 std::vector<std::string> element_block_names;
806 const shards::CellTopology & cell_topology = *(mesh.
getCellTopology(element_block_names[0]));
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);
813 Kokkos::DynRankView<double,PHX::Device> owned_vertices(
"owned_vertices",localCells.size(),vertices_per_cell,space_dim);
817 Kokkos::DynRankView<double,PHX::Device> ghstd_vertices = buildGhostedVertices(*cellimport_own2ghst,owned_vertices);
822 std::unordered_map<GO,int> global_to_local;
823 global_to_local[-1] = -1;
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));
831 faceToElement->initialize(conn);
832 auto elems_by_face = faceToElement->getFaceToElementsMap();
833 auto face_to_lidx = faceToElement->getFaceToCellLocalIdxMap();
835 const int num_owned_cells =owned_cells.extent(0);
836 const int num_ghstd_cells =ghstd_cells.extent(0);
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);
856 const LO num_virtual_cells = all_boundary_faces.size();
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);
865 Kokkos::View<GO*> virtual_cells = Kokkos::View<GO*>(
"virtual_cells",num_virtual_cells);
867 const int num_ranks = comm->getSize();
868 const int rank = comm->getRank();
870 std::vector<GO> owned_cell_distribution(num_ranks,0);
872 std::vector<GO> my_owned_cell_distribution(num_ranks,0);
873 my_owned_cell_distribution[rank] = num_owned_cells;
875 Teuchos::reduceAll(*comm,Teuchos::REDUCE_SUM, num_ranks, my_owned_cell_distribution.data(),owned_cell_distribution.data());
878 std::vector<GO> virtual_cell_distribution(num_ranks,0);
880 std::vector<GO> my_virtual_cell_distribution(num_ranks,0);
881 my_virtual_cell_distribution[rank] = num_virtual_cells;
883 Teuchos::reduceAll(*comm,Teuchos::REDUCE_SUM, num_ranks, my_virtual_cell_distribution.data(),virtual_cell_distribution.data());
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];
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];
896 for(
int i=0;i<num_virtual_cells;++i){
897 virtual_cells(i) = global_virtual_start_idx + GO(i);
903 Kokkos::View<LO*[2]> face_to_cells = Kokkos::View<LO*[2]>(
"face_to_cells",num_total_faces);
906 Kokkos::View<LO*[2]> face_to_localidx = Kokkos::View<LO*[2]>(
"face_to_localidx",num_total_faces);
909 Kokkos::View<LO**> cell_to_face = Kokkos::View<LO**>(
"cell_to_face",num_total_cells,faces_per_cell);
912 Kokkos::deep_copy(cell_to_face,-1);
916 int virtual_cell_index = num_real_cells;
917 for(
size_t f=0;f<elems_by_face.extent(0);f++) {
919 const GO global_c0 = elems_by_face(f,0);
920 const GO global_c1 = elems_by_face(f,1);
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());
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);
936 c0 = virtual_cell_index++;
943 cell_to_face(c0,lidx0) = f;
949 c1 = virtual_cell_index++;
956 cell_to_face(c1,lidx1) = f;
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;
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;
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);
993 for(
int i=0;i<num_owned_cells;++i){
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;
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;
1008 mesh_info.
cell_vertices = Kokkos::View<double***, PHX::Device>(
"cell_vertices",num_total_cells,vertices_per_cell,space_dim);
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){
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);
1031 for(
int i=0;i<num_virtual_cells;++i){
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);
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);
1042 for(
int j=0;j<vertices_per_cell;++j){
1043 for(
int k=0;k<space_dim;++k){
1055 std::vector<std::string> sideset_names;
1058 for(
const std::string & element_block_name : element_block_names){
1060 setupLocalMeshBlockInfo(mesh, conn, mesh_info, element_block_name, block_info);
1063 for(
const std::string & sideset_name : sideset_names){
1065 setupLocalMeshSidesetInfo(mesh, conn, mesh_info, element_block_name, sideset_name, sideset_info);
1080 #ifndef PANZER_ORDINAL64_IS_INT 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)
std::string element_block_name
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
std::string element_block_name
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)