Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Remove hardcoded inverse orientation #16806

Merged
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Jump to
Jump to file
Failed to load files.
Diff view
Diff view
8 changes: 6 additions & 2 deletions include/deal.II/grid/reference_cell.h
Original file line number Diff line number Diff line change
Expand Up @@ -1877,7 +1877,7 @@ ReferenceCell::face_to_cell_vertices(
if (get_dimension() == 1)
Assert(combined_face_orientation ==
ReferenceCell::default_combined_face_orientation(),
ExcMessage("In 1D, all cells must have the default orientation."));
ExcMessage("In 1D, all faces must have the default orientation."));
else
AssertIndexRange(combined_face_orientation, n_face_orientations(face));

Expand Down Expand Up @@ -1988,9 +1988,13 @@ ReferenceCell::standard_to_real_face_vertex(
switch (this->kind)
{
case ReferenceCells::Vertex:
case ReferenceCells::Line:
DEAL_II_NOT_IMPLEMENTED();
break;
case ReferenceCells::Line:
Assert(face_orientation == default_combined_face_orientation(),
ExcMessage(
"In 1D, all faces must have the default orientation."));
return vertex;
case ReferenceCells::Triangle:
case ReferenceCells::Quadrilateral:
return line_vertex_permutations[face_orientation][vertex];
Expand Down
54 changes: 8 additions & 46 deletions source/grid/grid_tools.cc
Original file line number Diff line number Diff line change
Expand Up @@ -4791,32 +4791,6 @@ namespace GridTools
// 1) determine for each vertex a vertex it coincides with and
// put it into a map
{
static const int lookup_table_2d[2][2] =
// flip:
{
{0, 1}, // false
{1, 0} // true
};

static const int lookup_table_3d[2][2][2][4] =
// orientation flip rotation
{{{
{0, 2, 1, 3}, // false false false
{2, 3, 0, 1} // false false true
},
{
{3, 1, 2, 0}, // false true false
{1, 0, 3, 2} // false true true
}},
{{
{0, 1, 2, 3}, // true false false
{1, 3, 0, 2} // true false true
},
{
{3, 2, 1, 0}, // true true false
{2, 0, 3, 1} // true true true
}}};

// loop over all periodic face pairs
for (const auto &pair : tria.get_periodic_face_map())
{
Expand All @@ -4826,34 +4800,22 @@ namespace GridTools
const auto face_a = pair.first.first->face(pair.first.second);
const auto face_b =
pair.second.first.first->face(pair.second.first.second);
const auto reference_cell = pair.first.first->reference_cell();
const auto face_reference_cell = face_a->reference_cell();
const unsigned char combined_orientation = pair.second.second;
const unsigned char inverse_combined_orientation =
face_reference_cell.get_inverse_combined_orientation(
combined_orientation);

AssertDimension(face_a->n_vertices(), face_b->n_vertices());

// loop over all vertices on face
for (unsigned int i = 0; i < face_a->n_vertices(); ++i)
{
const auto [face_orientation, face_rotation, face_flip] =
::dealii::internal::split_face_orientation(
combined_orientation);

// find the right local vertex index for the second face
unsigned int j = 0;
switch (dim)
{
case 1:
j = i;
break;
case 2:
j = lookup_table_2d[face_flip][i];
break;
case 3:
j = lookup_table_3d[face_orientation][face_flip]
[face_rotation][i];
break;
default:
AssertThrow(false, ExcNotImplemented());
}
const unsigned int j =
reference_cell.standard_to_real_face_vertex(
i, pair.first.second, inverse_combined_orientation);

// get vertex indices and store in map
const auto vertex_a = face_a->vertex_index(i);
Expand Down
53 changes: 10 additions & 43 deletions source/grid/tria.cc
Original file line number Diff line number Diff line change
Expand Up @@ -1657,6 +1657,10 @@ namespace
const FaceIterator face_1 = cell_1->face(n_face_1);
const FaceIterator face_2 = cell_2->face(n_face_2);

const unsigned char inverse_orientation =
face_1->reference_cell().get_inverse_combined_orientation(orientation);

#ifdef DEBUG
const auto [face_orientation, face_rotation, face_flip] =
internal::split_face_orientation(orientation);

Expand All @@ -1670,6 +1674,7 @@ namespace
ExcMessage("The supplied orientation "
"(face_orientation, face_flip, face_rotation) "
"is invalid for 2d"));
#endif

Assert(face_1 != face_2, ExcMessage("face_1 and face_2 are equal!"));

Expand Down Expand Up @@ -1727,36 +1732,6 @@ namespace
}
else // dim == 2 || dim == 3
{
// A lookup table on how to go through the child cells depending on the
// orientation:
// see Documentation of GeometryInfo for details

static const int lookup_table_2d[2][2] =
// orientation:
{
{1, 0}, // false
{0, 1} // true
};

static const int lookup_table_3d[2][2][2][4] =
// orientation flip rotation
{{{
{0, 2, 1, 3}, // false false false
{2, 3, 0, 1} // false false true
},
{
{3, 1, 2, 0}, // false true false
{1, 0, 3, 2} // false true true
}},
{{
{0, 1, 2, 3}, // true false false
{1, 3, 0, 2} // true false true
},
{
{3, 2, 1, 0}, // true true false
{2, 0, 3, 1} // true true true
}}};

if (cell_1->has_children())
{
if (cell_2->has_children())
Expand All @@ -1771,24 +1746,16 @@ namespace
GeometryInfo<dim>::max_children_per_face,
ExcNotImplemented());

const auto reference_cell = cell_1->reference_cell();

for (unsigned int i = 0;
i < GeometryInfo<dim>::max_children_per_face;
++i)
{
// Lookup the index for the second face
unsigned int j = 0;
switch (dim)
{
case 2:
j = lookup_table_2d[face_orientation][i];
break;
case 3:
j = lookup_table_3d[face_orientation][face_flip]
[face_rotation][i];
break;
default:
AssertThrow(false, ExcNotImplemented());
}
const unsigned int j =
reference_cell.standard_to_real_face_vertex(
i, n_face_1, inverse_orientation);

// find subcell ids that belong to the subface indices
unsigned int child_cell_1 =
Expand Down