62 #include "Phalanx_KokkosDeviceTypes.hpp"
64 #include "Teuchos_Assert.hpp"
67 #include "Tpetra_Import.hpp"
72 #include <unordered_set>
84 Kokkos::DynRankView<double,PHX::Device>
85 buildGhostedNodes(
const Tpetra::Import<int,panzer::GlobalOrdinal,panzer::TpetraNodeType> & importer,
86 Kokkos::DynRankView<const double,PHX::Device> owned_nodes)
91 typedef Tpetra::MultiVector<double,int,panzer::GlobalOrdinal,panzer::TpetraNodeType> mvec_type;
92 typedef typename mvec_type::dual_view_type dual_view_type;
94 size_t owned_cell_cnt = importer.getSourceMap()->getLocalNumElements();
95 size_t ghstd_cell_cnt = importer.getTargetMap()->getLocalNumElements();
96 int nodes_per_cell = owned_nodes.extent(1);
97 int space_dim = owned_nodes.extent(2);
102 RCP<mvec_type> owned_nodes_mv =
rcp(
new mvec_type(importer.getSourceMap(),nodes_per_cell*space_dim));
103 RCP<mvec_type> ghstd_nodes_mv =
rcp(
new mvec_type(importer.getTargetMap(),nodes_per_cell*space_dim));
106 auto owned_nodes_view = owned_nodes_mv->getLocalViewDevice(Tpetra::Access::OverwriteAll);
107 Kokkos::parallel_for(owned_cell_cnt, KOKKOS_LAMBDA (
size_t i) {
109 for(
int j=0;j<nodes_per_cell;j++)
110 for(
int k=0;k<space_dim;k++,l++)
111 owned_nodes_view(i,l) = owned_nodes(i,j,k);
116 ghstd_nodes_mv->doImport(*owned_nodes_mv,importer,Tpetra::INSERT);
119 Kokkos::DynRankView<double,PHX::Device> ghstd_nodes(
"ghstd_nodes",ghstd_cell_cnt,nodes_per_cell,space_dim);
121 auto ghstd_nodes_view = ghstd_nodes_mv->getLocalViewDevice(Tpetra::Access::ReadOnly);
122 Kokkos::parallel_for(ghstd_cell_cnt, KOKKOS_LAMBDA (
size_t i) {
124 for(
int j=0;j<nodes_per_cell;j++)
125 for(
int k=0;k<space_dim;k++,l++)
126 ghstd_nodes(i,j,k) = ghstd_nodes_view(i,l);
137 const std::string & element_block_name,
148 const shards::CellTopology & topology = *(mesh.
getCellTopology(element_block_name));
150 if(topology.getDimension() == 1){
152 }
else if(topology.getDimension() == 2){
154 }
else if(topology.getDimension() == 3){
159 PANZER_FUNC_TIME_MONITOR(
"Build connectivity");
164 std::vector<panzer::LocalOrdinal> owned_block_cells;
165 auto local_cells_h = Kokkos::create_mirror_view(mesh_info.
local_cells);
166 Kokkos::deep_copy(local_cells_h, mesh_info.
local_cells);
167 for(
int parent_owned_cell=0;parent_owned_cell<num_parent_owned_cells;++parent_owned_cell){
168 const panzer::LocalOrdinal local_cell = local_cells_h(parent_owned_cell);
169 const bool is_in_block = conn.
getBlockId(local_cell) == element_block_name;
172 owned_block_cells.push_back(parent_owned_cell);
177 if ( owned_block_cells.size() == 0 )
183 PANZER_FUNC_TIME_MONITOR(
"panzer::partitioning_utilities::setupSubLocalMeshInfo");
193 const std::string & element_block_name,
194 const std::string & sideset_name,
199 using LocalOrdinal = panzer::LocalOrdinal;
200 using ParentOrdinal = int;
201 using SubcellOrdinal = short;
207 std::vector<stk::mesh::Entity> side_entities;
213 mesh.
getAllSides(sideset_name, element_block_name, side_entities);
215 }
catch(STK_Interface::SidesetException & e) {
216 std::stringstream ss;
217 std::vector<std::string> sideset_names;
221 ss << e.what() <<
"\nChoose existing sideset:\n";
222 for(
const auto & optional_sideset_name : sideset_names){
223 ss <<
"\t\"" << optional_sideset_name <<
"\"\n";
228 }
catch(STK_Interface::ElementBlockException & e) {
229 std::stringstream ss;
230 std::vector<std::string> element_block_names;
234 ss << e.what() <<
"\nChoose existing element block:\n";
235 for(
const auto & optional_element_block_name : element_block_names){
236 ss <<
"\t\"" << optional_element_block_name <<
"\"\n";
241 }
catch(std::logic_error & e) {
242 std::stringstream ss;
243 ss << e.what() <<
"\nUnrecognized logic error.\n";
250 std::map<ParentOrdinal,std::vector<SubcellOrdinal> > owned_parent_cell_to_subcell_indexes;
254 const size_t face_subcell_dimension =
static_cast<size_t>(mesh.
getCellTopology(element_block_name)->getDimension()-1);
264 std::vector<stk::mesh::Entity> elements;
265 std::vector<size_t> subcell_indexes;
266 std::vector<size_t> subcell_dimensions;
268 const size_t num_elements = subcell_dimensions.size();
271 std::unordered_map<LocalOrdinal,ParentOrdinal> local_to_parent_lookup;
272 auto local_cells_h = Kokkos::create_mirror_view(mesh_info.
local_cells);
273 Kokkos::deep_copy(local_cells_h, mesh_info.
local_cells);
274 for(ParentOrdinal parent_index=0; parent_index<static_cast<ParentOrdinal>(mesh_info.
local_cells.extent(0)); ++parent_index)
275 local_to_parent_lookup[local_cells_h(parent_index)] = parent_index;
279 for(
size_t element_index=0; element_index<num_elements; ++element_index) {
280 const size_t subcell_dimension = subcell_dimensions[element_index];
283 if(subcell_dimension == face_subcell_dimension){
284 const SubcellOrdinal subcell_index =
static_cast<SubcellOrdinal
>(subcell_indexes[element_index]);
285 const LocalOrdinal element_local_index =
static_cast<LocalOrdinal
>(mesh.
elementLocalId(elements[element_index]));
288 const auto itr = local_to_parent_lookup.find(element_local_index);
290 const ParentOrdinal element_parent_index = itr->second;
292 owned_parent_cell_to_subcell_indexes[element_parent_index].push_back(subcell_index);
299 const panzer::LocalOrdinal num_owned_cells = owned_parent_cell_to_subcell_indexes.size();
308 face_t(
const ParentOrdinal c0,
309 const ParentOrdinal c1,
310 const SubcellOrdinal sc0,
311 const SubcellOrdinal sc1)
318 ParentOrdinal cell_0;
319 ParentOrdinal cell_1;
320 SubcellOrdinal subcell_index_0;
321 SubcellOrdinal subcell_index_1;
326 std::unordered_set<panzer::LocalOrdinal> owned_parent_cells_set, ghstd_parent_cells_set, virtual_parent_cells_set;
327 std::vector<face_t> faces;
329 auto cell_to_faces_h = Kokkos::create_mirror_view(mesh_info.
cell_to_faces);
330 auto face_to_cells_h = Kokkos::create_mirror_view(mesh_info.
face_to_cells);
331 auto face_to_lidx_h = Kokkos::create_mirror_view(mesh_info.
face_to_lidx);
334 Kokkos::deep_copy(face_to_lidx_h, mesh_info.
face_to_lidx);
336 for(
const auto & local_cell_index_pair : owned_parent_cell_to_subcell_indexes){
338 const ParentOrdinal parent_cell = local_cell_index_pair.first;
339 const auto & subcell_indexes = local_cell_index_pair.second;
341 owned_parent_cells_set.insert(parent_cell);
343 for(
const SubcellOrdinal & subcell_index : subcell_indexes){
345 const LocalOrdinal face = cell_to_faces_h(parent_cell, subcell_index);
346 const LocalOrdinal face_other_side = (face_to_cells_h(face,0) == parent_cell) ? 1 : 0;
348 TEUCHOS_ASSERT(subcell_index == face_to_lidx_h(face, 1-face_other_side));
350 const ParentOrdinal other_side_cell = face_to_cells_h(face, face_other_side);
351 const SubcellOrdinal other_side_subcell_index = face_to_lidx_h(face, face_other_side);
353 faces.push_back(face_t(parent_cell, other_side_cell, subcell_index, other_side_subcell_index));
355 if(other_side_cell >= parent_virtual_cell_offset){
356 virtual_parent_cells_set.insert(other_side_cell);
358 ghstd_parent_cells_set.insert(other_side_cell);
364 std::vector<ParentOrdinal> all_cells;
365 all_cells.insert(all_cells.end(),owned_parent_cells_set.begin(),owned_parent_cells_set.end());
366 all_cells.insert(all_cells.end(),ghstd_parent_cells_set.begin(),ghstd_parent_cells_set.end());
367 all_cells.insert(all_cells.end(),virtual_parent_cells_set.begin(),virtual_parent_cells_set.end());
373 const LocalOrdinal num_total_cells = num_real_cells + sideset_info.
num_virtual_cells;
374 const LocalOrdinal num_nodes_per_cell = mesh_info.
cell_nodes.extent(1);
375 const LocalOrdinal num_dims = mesh_info.
cell_nodes.extent(2);
379 sideset_info.
global_cells = PHX::View<panzer::GlobalOrdinal*>(
"global_cells", num_total_cells);
380 sideset_info.
local_cells = PHX::View<LocalOrdinal*>(
"local_cells", num_total_cells);
381 sideset_info.
cell_nodes = PHX::View<double***>(
"cell_nodes", num_total_cells, num_nodes_per_cell, num_dims);
382 Kokkos::deep_copy(sideset_info.
cell_nodes,0.);
384 typename PHX::View<ParentOrdinal*>::HostMirror all_cells_h(
"all_cells_h", num_total_cells);
385 PHX::View<ParentOrdinal*> all_cells_d(
"all_cells_d", num_total_cells);
386 for(LocalOrdinal i=0; i<num_total_cells; ++i)
387 all_cells_h(i) = all_cells[i];
388 Kokkos::deep_copy(all_cells_d, all_cells_h);
389 Kokkos::parallel_for(num_total_cells, KOKKOS_LAMBDA (LocalOrdinal i) {
390 const ParentOrdinal parent_cell = all_cells_d(i);
393 for(LocalOrdinal j=0; j<num_nodes_per_cell; ++j)
394 for(LocalOrdinal k=0; k<num_dims; ++k)
401 const LocalOrdinal num_faces = faces.size();
402 const LocalOrdinal num_faces_per_cell = mesh_info.
cell_to_faces.extent(1);
404 sideset_info.
face_to_cells = PHX::View<LocalOrdinal*[2]>(
"face_to_cells", num_faces);
405 sideset_info.
face_to_lidx = PHX::View<LocalOrdinal*[2]>(
"face_to_lidx", num_faces);
406 sideset_info.
cell_to_faces = PHX::View<LocalOrdinal**>(
"cell_to_faces", num_total_cells, num_faces_per_cell);
407 auto cell_to_faces_h = Kokkos::create_mirror_view(sideset_info.
cell_to_faces);
408 auto face_to_cells_h = Kokkos::create_mirror_view(sideset_info.
face_to_cells);
409 auto face_to_lidx_h = Kokkos::create_mirror_view(sideset_info.
face_to_lidx);
412 Kokkos::deep_copy(cell_to_faces_h, -1);
415 std::unordered_map<ParentOrdinal,ParentOrdinal> parent_to_sideset_lookup;
416 for(
unsigned int i=0; i<all_cells.size(); ++i)
417 parent_to_sideset_lookup[all_cells[i]] = i;
419 for(
int face_index=0;face_index<num_faces;++face_index){
420 const face_t & face = faces[face_index];
421 const ParentOrdinal & parent_cell_0 = face.cell_0;
422 const ParentOrdinal & parent_cell_1 = face.cell_1;
425 const auto itr_0 = parent_to_sideset_lookup.find(parent_cell_0);
427 const ParentOrdinal sideset_cell_0 = itr_0->second;
429 const auto itr_1 = parent_to_sideset_lookup.find(parent_cell_1);
431 const ParentOrdinal sideset_cell_1 = itr_1->second;
436 face_to_cells_h(face_index,0) = sideset_cell_0;
437 face_to_cells_h(face_index,1) = sideset_cell_1;
439 face_to_lidx_h(face_index,0) = face.subcell_index_0;
440 face_to_lidx_h(face_index,1) = face.subcell_index_1;
442 cell_to_faces_h(sideset_cell_0,face.subcell_index_0) = face_index;
443 cell_to_faces_h(sideset_cell_1,face.subcell_index_1) = face_index;
446 Kokkos::deep_copy(sideset_info.
cell_to_faces, cell_to_faces_h);
447 Kokkos::deep_copy(sideset_info.
face_to_cells, face_to_cells_h);
448 Kokkos::deep_copy(sideset_info.
face_to_lidx, face_to_lidx_h );
463 typedef Tpetra::Map<panzer::LocalOrdinal,panzer::GlobalOrdinal,panzer::TpetraNodeType> map_type;
464 typedef Tpetra::Import<panzer::LocalOrdinal,panzer::GlobalOrdinal,panzer::TpetraNodeType> import_type;
469 auto & mesh_info = *mesh_info_rcp;
487 PHX::View<panzer::GlobalOrdinal*> owned_cells, ghost_cells, virtual_cells;
497 RCP<import_type> cellimport_own2ghst =
rcp(
new import_type(owned_cell_map,ghstd_cell_map));
503 std::vector<std::string> element_block_names;
506 const shards::CellTopology & cell_topology = *(mesh.
getCellTopology(element_block_names[0]));
508 const int space_dim = cell_topology.getDimension();
509 const int nodes_per_cell = cell_topology.getNodeCount();
510 const int faces_per_cell = cell_topology.getSubcellCount(space_dim-1);
512 Kokkos::DynRankView<double,PHX::Device> owned_nodes(
"owned_nodes",owned_cells.extent(0),nodes_per_cell,space_dim);
515 for(
size_t i=0;i<localCells.size();i++)
521 Kokkos::DynRankView<double,PHX::Device> ghstd_nodes = buildGhostedNodes(*cellimport_own2ghst,owned_nodes);
526 auto owned_cells_h = Kokkos::create_mirror_view(owned_cells);
527 auto ghost_cells_h = Kokkos::create_mirror_view(ghost_cells);
528 Kokkos::deep_copy(owned_cells_h, owned_cells);
529 Kokkos::deep_copy(ghost_cells_h, ghost_cells);
530 std::unordered_map<panzer::GlobalOrdinal,int> global_to_local;
531 global_to_local[-1] = -1;
532 for(
size_t i=0;i<owned_cells.extent(0);i++)
533 global_to_local[owned_cells_h(i)] = i;
534 for(
size_t i=0;i<ghost_cells.extent(0);i++)
535 global_to_local[ghost_cells_h(i)] = i+Teuchos::as<int>(owned_cells.extent(0));
539 faceToElement->initialize(conn, comm);
540 auto elems_by_face = faceToElement->getFaceToElementsMap();
541 auto face_to_lidx = faceToElement->getFaceToCellLocalIdxMap();
553 const panzer::LocalOrdinal num_owned_cells = owned_cells.extent(0);
554 const panzer::LocalOrdinal num_ghstd_cells = ghost_cells.extent(0);
555 const panzer::LocalOrdinal num_virtual_cells = virtual_cells.extent(0);
558 const panzer::LocalOrdinal num_real_cells = num_owned_cells + num_ghstd_cells;
559 const panzer::LocalOrdinal num_total_cells = num_real_cells + num_virtual_cells;
560 const panzer::LocalOrdinal num_total_faces = elems_by_face.extent(0);
563 PHX::View<panzer::LocalOrdinal*[2]> face_to_cells(
"face_to_cells",num_total_faces);
566 PHX::View<panzer::LocalOrdinal*[2]> face_to_localidx(
"face_to_localidx",num_total_faces);
569 PHX::View<panzer::LocalOrdinal**> cell_to_face(
"cell_to_face",num_total_cells,faces_per_cell);
573 PANZER_FUNC_TIME_MONITOR_DIFF(
"Transer faceToElement to local",TransferFaceToElementLocal);
575 int virtual_cell_index = num_real_cells;
576 auto elems_by_face_h = Kokkos::create_mirror_view(elems_by_face);
577 auto face_to_lidx_h = Kokkos::create_mirror_view(face_to_lidx);
578 auto face_to_cells_h = Kokkos::create_mirror_view(face_to_cells);
579 auto face_to_localidx_h = Kokkos::create_mirror_view(face_to_localidx);
580 auto cell_to_face_h = Kokkos::create_mirror_view(cell_to_face);
581 Kokkos::deep_copy(elems_by_face_h, elems_by_face);
582 Kokkos::deep_copy(face_to_lidx_h, face_to_lidx);
584 Kokkos::deep_copy(cell_to_face_h, -1);
585 for(
size_t f=0;f<elems_by_face.extent(0);f++) {
587 const panzer::GlobalOrdinal global_c0 = elems_by_face_h(f,0);
588 const panzer::GlobalOrdinal global_c1 = elems_by_face_h(f,1);
591 TEUCHOS_ASSERT(global_to_local.find(global_c0)!=global_to_local.end());
592 TEUCHOS_ASSERT(global_to_local.find(global_c1)!=global_to_local.end());
594 auto c0 = global_to_local[global_c0];
595 auto lidx0 = face_to_lidx_h(f,0);
596 auto c1 = global_to_local[global_c1];
597 auto lidx1 = face_to_lidx_h(f,1);
604 c0 = virtual_cell_index++;
611 cell_to_face_h(c0,lidx0) = f;
617 c1 = virtual_cell_index++;
624 cell_to_face_h(c1,lidx1) = f;
628 face_to_cells_h(f,0) = c0;
629 face_to_localidx_h(f,0) = lidx0;
630 face_to_cells_h(f,1) = c1;
631 face_to_localidx_h(f,1) = lidx1;
633 face_to_cells_h(f,0) = c1;
634 face_to_localidx_h(f,0) = lidx1;
635 face_to_cells_h(f,1) = c0;
636 face_to_localidx_h(f,1) = lidx0;
643 Kokkos::deep_copy(face_to_cells, face_to_cells_h);
644 Kokkos::deep_copy(face_to_localidx, face_to_localidx_h);
645 Kokkos::deep_copy(cell_to_face, cell_to_face_h);
652 PANZER_FUNC_TIME_MONITOR_DIFF(
"Assign Indices",AssignIndices);
664 mesh_info.
global_cells = PHX::View<panzer::GlobalOrdinal*>(
"global_cell_indices",num_total_cells);
665 mesh_info.
local_cells = PHX::View<panzer::LocalOrdinal*>(
"local_cell_indices",num_total_cells);
667 Kokkos::parallel_for(num_owned_cells,KOKKOS_LAMBDA (
int i) {
672 Kokkos::parallel_for(num_ghstd_cells,KOKKOS_LAMBDA (
int i) {
673 mesh_info.
global_cells(i+num_owned_cells) = ghost_cells(i);
674 mesh_info.
local_cells(i+num_owned_cells) = i+num_owned_cells;
677 Kokkos::parallel_for(num_virtual_cells,KOKKOS_LAMBDA (
int i) {
678 mesh_info.
global_cells(i+num_real_cells) = virtual_cells(i);
679 mesh_info.
local_cells(i+num_real_cells) = i+num_real_cells;
682 mesh_info.
cell_nodes = PHX::View<double***>(
"cell_nodes",num_total_cells,nodes_per_cell,space_dim);
687 Kokkos::parallel_for(num_owned_cells,KOKKOS_LAMBDA (
int i) {
688 for(
int j=0;j<nodes_per_cell;++j){
689 for(
int k=0;k<space_dim;++k){
690 mesh_info.
cell_nodes(i,j,k) = owned_nodes(i,j,k);
695 Kokkos::parallel_for(num_ghstd_cells,KOKKOS_LAMBDA (
int i) {
696 for(
int j=0;j<nodes_per_cell;++j){
697 for(
int k=0;k<space_dim;++k){
698 mesh_info.
cell_nodes(i+num_owned_cells,j,k) = ghstd_nodes(i,j,k);
707 PANZER_FUNC_TIME_MONITOR_DIFF(
"Assign geometry traits",AssignGeometryTraits);
708 Kokkos::parallel_for(num_virtual_cells,KOKKOS_LAMBDA (
int i) {
709 const panzer::LocalOrdinal virtual_cell = i+num_real_cells;
710 for(
int local_face=0; local_face<faces_per_cell; ++local_face){
711 const panzer::LocalOrdinal face = cell_to_face(virtual_cell, local_face);
713 const panzer::LocalOrdinal other_side = (face_to_cells(face, 0) == virtual_cell) ? 1 : 0;
714 const panzer::LocalOrdinal real_cell = face_to_cells(face,other_side);
715 for(
int j=0;j<nodes_per_cell;++j){
716 for(
int k=0;k<space_dim;++k){
729 std::vector<std::string> sideset_names;
732 for(
const std::string & element_block_name : element_block_names){
733 PANZER_FUNC_TIME_MONITOR_DIFF(
"Set up setupLocalMeshBlockInfo",SetupLocalMeshBlockInfo);
735 setupLocalMeshBlockInfo(mesh, conn, mesh_info, element_block_name, block_info);
741 for(
const std::string & sideset_name : sideset_names){
742 PANZER_FUNC_TIME_MONITOR_DIFF(
"Setup LocalMeshSidesetInfo",SetupLocalMeshSidesetInfo);
744 setupLocalMeshSidesetInfo(mesh, conn, mesh_info, element_block_name, sideset_name, sideset_info);
752 return mesh_info_rcp;
void getSidesetNames(std::vector< std::string > &name) const
void getElementNodesNoResize(const std::vector< std::size_t > &localIds, ArrayT &nodes) const
std::map< std::string, std::map< std::string, LocalMeshSidesetInfo > > sidesets
std::string element_block_name
void getElementBlockNames(std::vector< std::string > &names) const
Teuchos::RCP< panzer::LocalMeshInfo > generateLocalMeshInfo(const panzer_stk::STK_Interface &mesh)
Create a structure containing information about the local portion of a given element block...
Teuchos::RCP< const shards::CellTopology > cell_topology
PHX::View< panzer::LocalOrdinal * > local_cells
void fillLocalCellIDs(const Teuchos::RCP< const Teuchos::Comm< int >> &comm, panzer::ConnManager &conn, PHX::View< panzer::GlobalOrdinal * > &owned_cells, PHX::View< panzer::GlobalOrdinal * > &ghost_cells, PHX::View< panzer::GlobalOrdinal * > &virtual_cells)
Get the owned, ghost and virtual global cell ids for this process.
PHX::View< panzer::LocalOrdinal *[2]> face_to_cells
panzer::LocalOrdinal num_owned_cells
std::size_t elementLocalId(stk::mesh::Entity elmt) const
panzer::LocalOrdinal num_virtual_cells
Teuchos::RCP< const shards::CellTopology > cell_topology
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
PHX::View< panzer::LocalOrdinal *[2]> face_to_lidx
Pure virtual base class for supplying mesh connectivity information to the DOF Manager.
panzer::LocalOrdinal num_ghstd_cells
std::map< std::string, LocalMeshBlockInfo > element_blocks
bool isInitialized() const
Has initialize been called on this mesh object?
void getAllSides(const std::string &sideName, std::vector< stk::mesh::Entity > &sides) const
std::string element_block_name
Teuchos::RCP< const Teuchos::Comm< int > > getComm() const
get the comm associated with this mesh
PHX::View< double *** > cell_nodes
#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)
void setupSubLocalMeshInfo(const panzer::LocalMeshInfoBase &parent_info, const std::vector< panzer::LocalOrdinal > &owned_parent_cells, panzer::LocalMeshInfoBase &sub_info)
PHX::View< panzer::GlobalOrdinal * > global_cells
virtual void buildConnectivity(const FieldPattern &fp)=0
Teuchos::RCP< const shards::CellTopology > getCellTopology(const std::string &eBlock) const
#define TEUCHOS_FUNC_TIME_MONITOR_DIFF(FUNCNAME, DIFF)
virtual std::string getBlockId(LocalOrdinal localElmtId) const =0
#define TEUCHOS_TEST_FOR_EXCEPTION_PURE_MSG(throw_exception_test, Exception, msg)
PHX::View< panzer::LocalOrdinal ** > cell_to_faces