45 #include "Teuchos_RCP.hpp" 46 #include "Teuchos_Comm.hpp" 47 #include "Teuchos_Assert.hpp" 49 #include "Panzer_Workset_Builder.hpp" 52 #include "Phalanx_KokkosDeviceTypes.hpp" 54 #include <unordered_set> 55 #include <unordered_map> 60 namespace partitioning_utilities
63 template<
typename LO,
typename GO>
66 const std::vector<LO> & owned_parent_cells,
83 const int num_owned_cells = owned_parent_cells.size();
84 TEUCHOS_TEST_FOR_EXCEPT_MSG(num_owned_cells == 0,
"panzer::partitioning_utilities::setupSubLocalMeshInfo : Input parent subcells must exist (owned_parent_cells)");
87 TEUCHOS_TEST_FOR_EXCEPT_MSG(num_parent_owned_cells <= 0,
"panzer::partitioning_utilities::setupSubLocalMeshInfo : Input parent info must contain owned cells");
94 TEUCHOS_ASSERT(static_cast<int>(parent_info.
cell_to_faces.extent(0)) == num_parent_total_cells);
95 const int num_faces_per_cell = parent_info.
cell_to_faces.extent(1);
99 std::vector<LO> ghstd_parent_cells;
100 std::vector<LO> virtual_parent_cells;
105 std::unordered_set<LO> owned_parent_cells_set(owned_parent_cells.begin(), owned_parent_cells.end());
112 const int virtual_parent_cell_offset = num_parent_owned_cells + num_parent_ghstd_cells;
114 std::unordered_set<LO> ghstd_parent_cells_set;
115 std::unordered_set<LO> virtual_parent_cells_set;
116 for(
int i=0;i<num_owned_cells;++i){
117 const LO parent_cell_index = owned_parent_cells[i];
118 for(
int local_face_index=0;local_face_index<num_faces_per_cell;++local_face_index){
119 const LO parent_face = parent_info.
cell_to_faces(parent_cell_index, local_face_index);
127 const LO neighbor_parent_side = (parent_info.
face_to_cells(parent_face,0) == parent_cell_index) ? 1 : 0;
130 const LO neighbor_parent_cell = parent_info.
face_to_cells(parent_face, neighbor_parent_side);
133 TEUCHOS_ASSERT(neighbor_parent_cell >= 0);
136 if(neighbor_parent_cell >= virtual_parent_cell_offset){
137 virtual_parent_cells_set.insert(neighbor_parent_cell);
138 }
else if(neighbor_parent_cell >= num_parent_owned_cells){
141 ghstd_parent_cells_set.insert(neighbor_parent_cell);
145 if(owned_parent_cells_set.find(neighbor_parent_cell) == owned_parent_cells_set.end()){
147 ghstd_parent_cells_set.insert(neighbor_parent_cell);
156 ghstd_parent_cells.assign(ghstd_parent_cells_set.begin(), ghstd_parent_cells_set.end());
157 virtual_parent_cells.assign(virtual_parent_cells_set.begin(), virtual_parent_cells_set.end());
161 const int num_ghstd_cells = ghstd_parent_cells.size();
162 const int num_virtual_cells = virtual_parent_cells.size();
163 const int num_real_cells = num_owned_cells + num_ghstd_cells;
164 const int num_total_cells = num_real_cells + num_virtual_cells;
166 std::vector<LO> all_parent_cells;
167 all_parent_cells.insert(all_parent_cells.end(), owned_parent_cells.begin(), owned_parent_cells.end());
168 all_parent_cells.insert(all_parent_cells.end(), ghstd_parent_cells.begin(), ghstd_parent_cells.end());
169 all_parent_cells.insert(all_parent_cells.end(), virtual_parent_cells.begin(), virtual_parent_cells.end());
178 TEUCHOS_ASSERT(static_cast<int>(parent_info.
cell_vertices.extent(0)) == num_parent_total_cells);
179 TEUCHOS_ASSERT(static_cast<int>(parent_info.
local_cells.extent(0)) == num_parent_total_cells);
180 TEUCHOS_ASSERT(static_cast<int>(parent_info.
global_cells.extent(0)) == num_parent_total_cells);
182 const int num_vertices_per_cell = parent_info.
cell_vertices.extent(1);
186 sub_info.
global_cells = Kokkos::View<GO*>(
"global_cells", num_total_cells);
187 sub_info.
local_cells = Kokkos::View<LO*>(
"local_cells", num_total_cells);
188 sub_info.
cell_vertices = Kokkos::View<double***, PHX::Device>(
"cell_vertices", num_total_cells, num_vertices_per_cell, num_dims);
189 for(
int cell=0;cell<num_total_cells;++cell){
190 const LO parent_cell = all_parent_cells[cell];
193 for(
int vertex=0;vertex<num_vertices_per_cell;++vertex){
194 for(
int dim=0;dim<num_dims;++dim){
207 face_t(LO c0, LO c1, LO sc0, LO sc1)
222 std::vector<face_t> faces;
226 std::unordered_map<LO,std::unordered_map<LO, std::pair<LO,LO> > > faces_set;
227 for(
int owned_cell=0;owned_cell<num_owned_cells;++owned_cell){
228 const LO owned_parent_cell = owned_parent_cells[owned_cell];
229 for(
int local_face=0;local_face<num_faces_per_cell;++local_face){
230 const LO parent_face = parent_info.
cell_to_faces(owned_parent_cell,local_face);
237 const LO neighbor_side = (parent_info.
face_to_cells(parent_face,0) == owned_parent_cell) ? 1 : 0;
239 const LO neighbor_parent_cell = parent_info.
face_to_cells(parent_face, neighbor_side);
240 const LO neighbor_subcell_index = parent_info.
face_to_lidx(parent_face, neighbor_side);
243 auto itr = std::find(all_parent_cells.begin(), all_parent_cells.end(), neighbor_parent_cell);
244 TEUCHOS_TEST_FOR_EXCEPT_MSG(itr == all_parent_cells.end(),
"panzer_stk::setupSubLocalMeshInfo : Neighbor cell was not found in owned, ghosted, or virtual cells");
246 const LO neighbor_cell = std::distance(all_parent_cells.begin(), itr);
248 LO cell_0, cell_1, subcell_index_0, subcell_index_1;
249 if(owned_cell < neighbor_cell){
251 subcell_index_0 = local_face;
252 cell_1 = neighbor_cell;
253 subcell_index_1 = neighbor_subcell_index;
256 subcell_index_1 = local_face;
257 cell_0 = neighbor_cell;
258 subcell_index_0 = neighbor_subcell_index;
262 faces_set[cell_0][subcell_index_0] = std::pair<LO,LO>(cell_1, subcell_index_1);
266 for(
const auto & cell_pair : faces_set){
267 const LO cell_0 = cell_pair.first;
268 for(
const auto & subcell_pair : cell_pair.second){
269 const LO subcell_index_0 = subcell_pair.first;
270 const LO cell_1 = subcell_pair.second.first;
271 const LO subcell_index_1 = subcell_pair.second.second;
272 faces.push_back(face_t(cell_0,cell_1,subcell_index_0,subcell_index_1));
277 const int num_faces = faces.size();
279 sub_info.
face_to_cells = Kokkos::View<LO*[2]>(
"face_to_cells", num_faces);
280 sub_info.
face_to_lidx = Kokkos::View<LO*[2]>(
"face_to_lidx", num_faces);
281 sub_info.
cell_to_faces = Kokkos::View<LO**>(
"cell_to_faces", num_total_cells, num_faces_per_cell);
286 for(
int face_index=0;face_index<num_faces;++face_index){
287 const face_t & face = faces[face_index];
292 sub_info.
cell_to_faces(face.cell_0,face.subcell_index_0) = face_index;
293 sub_info.
cell_to_faces(face.cell_1,face.subcell_index_1) = face_index;
295 sub_info.
face_to_lidx(face_index,0) = face.subcell_index_0;
296 sub_info.
face_to_lidx(face_index,1) = face.subcell_index_1;
308 template<
typename LO,
typename GO>
311 const int splitting_size,
316 TEUCHOS_ASSERT(splitting_size != 0);
324 if(splitting_size < 0){
327 std::vector<LO> partition_cells;
330 partition_cells.push_back(i);
334 if(partition_cells.size() > 0){
341 std::vector<LO> partition_cells;
342 partition_cells.reserve(splitting_size);
347 LO partition_start_index = 0;
348 for(LO partition=0;partition<num_partitions;++partition){
351 const LO partition_end_index = std::min(partition_start_index + splitting_size, num_owned_cells);
353 partition_cells.resize(partition_end_index - partition_start_index,-1);
354 for(
int i=partition_start_index;i<partition_end_index;++i){
355 partition_cells[i] = i;
359 if(partition_cells.size() > 0){
364 partition_start_index = partition_end_index;
371 template<
typename LO,
typename GO>
378 TEUCHOS_ASSERT(
false);
383 template <
typename LO,
typename GO>
397 const std::string & element_block_name = description.
getElementBlock();
405 const auto & sideset_map = mesh_info.
sidesets.at(element_block_name);
407 const std::string & sideset_name = description.
getSideset();
410 if(sideset_map.find(sideset_name) == sideset_map.end())
416 panzer::partitioning_utilities::splitMeshInfo<LO,GO>(sideset_info, description.
getWorksetSize(), partitions);
418 for(
auto & partition : partitions){
420 partition.element_block_name = element_block_name;
421 partition.cell_topology = sideset_info.cell_topology;
447 for(
auto & partition : partitions){
448 partition.element_block_name = element_block_name;
461 const std::vector<int> & owned_parent_cells,
464 #ifndef PANZER_ORDINAL64_IS_INT 468 const std::vector<int> & owned_parent_cells,
476 std::vector<panzer::LocalMeshPartition<int,int> > & partitions);
478 #ifndef PANZER_ORDINAL64_IS_INT 483 std::vector<panzer::LocalMeshPartition<int,panzer::Ordinal64> > & partitions);
Backwards compatibility mode that ignores the worksetSize in the WorksetDescriptor.
std::map< std::string, LocalMeshBlockInfo< LO, GO > > element_blocks
Kokkos::View< GO * > global_cells
Teuchos::RCP< const shards::CellTopology > cell_topology
const std::string & getSideset(const int block=0) const
Get sideset name.
void setupSubLocalMeshInfo(const panzer::LocalMeshInfoBase< LO, GO > &parent_info, const std::vector< LO > &owned_parent_cells, panzer::LocalMeshInfoBase< LO, GO > &sub_info)
void splitMeshInfo(const panzer::LocalMeshInfoBase< LO, GO > &mesh_info, const int splitting_size, std::vector< panzer::LocalMeshPartition< LO, GO > > &partitions)
std::map< std::string, std::map< std::string, LocalMeshSidesetInfo< LO, GO > > > sidesets
Kokkos::View< LO * > local_cells
bool requiresPartitioning() const
Do we need to partition the local mesh prior to generating worksets.
int getWorksetSize() const
Get the requested workset size (default -2 (workset size is set elsewhere), -1 (largest possible work...
bool useSideset() const
This descriptor is for a side set.
Kokkos::View< LO *[2]> face_to_cells
Kokkos::View< double ***, PHX::Device > cell_vertices
Kokkos::View< LO ** > cell_to_faces
Kokkos::View< LO *[2]> face_to_lidx
void partitionMeshInfo(const panzer::LocalMeshInfoBase< LO, GO > &, const size_t, std::vector< panzer::LocalMeshPartition< LO, GO >> &)
void generateLocalMeshPartitions(const panzer::LocalMeshInfo< LO, GO > &mesh_info, const panzer::WorksetDescriptor &description, std::vector< panzer::LocalMeshPartition< LO, GO > > &partitions)
const std::string & getElementBlock(const int block=0) const
Get element block name.
Workset size is set to the total number of local elements in the MPI process.