From 0d0b243b8658ea0d663440b50612229796abb756 Mon Sep 17 00:00:00 2001 From: David Wells Date: Wed, 20 Dec 2023 09:49:41 -0500 Subject: [PATCH] Use face_orientation to encode 2D orientation information. More progress towards the goal of consistently encoding orientations in the library. Ultimately, 0 should be the default orientation in all dimensions. With this change only orientations 0 and 1 are valid in 2D (though 1 is the default), whereas before it was 1 and 3. --- .../incompatibilities/20231220DavidWells | 4 +++ include/deal.II/base/geometry_info.h | 20 ++++++------ include/deal.II/grid/reference_cell.h | 8 +++++ .../deal.II/grid/tria_accessor.templates.h | 6 ++-- source/distributed/tria.cc | 6 ++-- source/dofs/dof_tools_constraints.cc | 2 +- source/fe/fe_q_base.cc | 6 ++-- source/grid/grid_tools_dof_handlers.cc | 14 ++------ source/grid/tria.cc | 10 +++--- tests/base/geometry_info_2.output | 16 +++++----- tests/base/geometry_info_8.output | 32 +++++++++---------- tests/dofs/dof_tools_21_b.output | 4 +-- tests/dofs/dof_tools_21_b_x.cc | 4 +-- tests/dofs/dof_tools_21_b_x_q3.cc | 4 +-- tests/dofs/dof_tools_21_b_y.cc | 4 +-- tests/dofs/dof_tools_21_c.output | 4 +-- tests/fe/face_to_cell_q1_2d.cc | 16 ++++++---- tests/fe/face_to_cell_q1_2d.output | 24 +++++++------- tests/fe/face_to_cell_q2_2d.cc | 16 ++++++---- tests/fe/face_to_cell_q2_2d.output | 24 +++++++------- tests/fe/face_to_cell_q2xq2_2d.cc | 16 ++++++---- tests/fe/face_to_cell_q2xq2_2d.output | 24 +++++++------- tests/fe/face_to_cell_q3_2d.cc | 16 ++++++---- tests/fe/face_to_cell_q3_2d.output | 24 +++++++------- tests/fe/face_to_cell_q3xq4_2d.cc | 16 ++++++---- tests/fe/face_to_cell_q3xq4_2d.output | 24 +++++++------- tests/fe/face_to_cell_q4_2d.cc | 16 ++++++---- tests/fe/face_to_cell_q4_2d.output | 24 +++++++------- tests/grid/grid_tools_06.output | 2 +- 29 files changed, 206 insertions(+), 180 deletions(-) create mode 100644 doc/news/changes/incompatibilities/20231220DavidWells diff --git a/doc/news/changes/incompatibilities/20231220DavidWells b/doc/news/changes/incompatibilities/20231220DavidWells new file mode 100644 index 0000000000..2ad2e248e7 --- /dev/null +++ b/doc/news/changes/incompatibilities/20231220DavidWells @@ -0,0 +1,4 @@ +Changed: Meshes in 2D now use the face_orientation boolean (instead of the +face_flip boolean) to represent their orientation. +
+(David Wells, 2023/12/20) diff --git a/include/deal.II/base/geometry_info.h b/include/deal.II/base/geometry_info.h index 9b7f66d907..715adfaa00 100644 --- a/include/deal.II/base/geometry_info.h +++ b/include/deal.II/base/geometry_info.h @@ -4410,8 +4410,8 @@ inline unsigned int GeometryInfo<2>::child_cell_on_face(const RefinementCase<2> &ref_case, const unsigned int face, const unsigned int subface, - const bool /*face_orientation*/, - const bool face_flip, + const bool face_orientation, + const bool /*face_flip*/, const bool /*face_rotation*/, const RefinementCase<1> &) { @@ -4430,19 +4430,19 @@ GeometryInfo<2>::child_cell_on_face(const RefinementCase<2> &ref_case, [/* number of different ways to refine a cell */ 4] [/* faces_per_cell */ 4][/* max_children_per_face */ 2] = { { - // Normal orientation (face_flip = false) - {{0, 0}, {1, 1}, {0, 1}, {0, 1}}, // cut_x - {{0, 1}, {0, 1}, {0, 0}, {1, 1}}, // cut_y - {{0, 2}, {1, 3}, {0, 1}, {2, 3}} // cut_xy, i.e., isotropic - }, - { - // Flipped orientation (face_flip = true) + // Flipped orientation (face_orientation = false) {{0, 0}, {1, 1}, {1, 0}, {1, 0}}, // cut_x {{1, 0}, {1, 0}, {0, 0}, {1, 1}}, // cut_y {{2, 0}, {3, 1}, {1, 0}, {3, 2}} // cut_xy, i.e., isotropic + }, + { + // Normal orientation (face_orientation = true) + {{0, 0}, {1, 1}, {0, 1}, {0, 1}}, // cut_x + {{0, 1}, {0, 1}, {0, 0}, {1, 1}}, // cut_y + {{0, 2}, {1, 3}, {0, 1}, {2, 3}} // cut_xy, i.e., isotropic }}; - return subcells[face_flip][ref_case - 1][face][subface]; + return subcells[face_orientation][ref_case - 1][face][subface]; } diff --git a/include/deal.II/grid/reference_cell.h b/include/deal.II/grid/reference_cell.h index f7ab3ec133..dbde9f4121 100644 --- a/include/deal.II/grid/reference_cell.h +++ b/include/deal.II/grid/reference_cell.h @@ -1865,6 +1865,14 @@ ReferenceCell::face_to_cell_vertices( { AssertIndexRange(face, n_faces()); AssertIndexRange(vertex, face_reference_cell(face).n_vertices()); + // TODO: once the default orientation is switched to 0 then we can remove this + // special case for 1D. + if (get_dimension() == 1) + Assert(combined_face_orientation == + ReferenceCell::default_combined_face_orientation(), + ExcMessage("In 1D, all cells must have the default orientation.")); + else + AssertIndexRange(combined_face_orientation, n_face_orientations(face)); switch (this->kind) { diff --git a/include/deal.II/grid/tria_accessor.templates.h b/include/deal.II/grid/tria_accessor.templates.h index cd8f4d0ac1..ddac77f8c7 100644 --- a/include/deal.II/grid/tria_accessor.templates.h +++ b/include/deal.II/grid/tria_accessor.templates.h @@ -739,10 +739,8 @@ namespace internal * In 1d, face_flip is always false as there is no such concept as * "flipped" faces in 1d. * - * In 2d, we currently only support meshes where all faces are in - * standard orientation, so the result is also false. This also - * matches the fact that one can *always* orient faces in 2d in such a - * way that the don't need to be flipped + * In 2d the orientation is a single bit (i.e., the orientation of a + * line) encoded in face_orientation so face_flip is also always false. */ return false; } diff --git a/source/distributed/tria.cc b/source/distributed/tria.cc index adc93a90c6..7368fc0dbb 100644 --- a/source/distributed/tria.cc +++ b/source/distributed/tria.cc @@ -3748,13 +3748,13 @@ namespace parallel // p4est wants to know which corner the first corner on // the face with the lower id is mapped to on the face with // with the higher id. For d==2 there are only two possibilities - // that are determined by it->orientation[1]. - // For d==3 we have to use GridTools::OrientationLookupTable. + // that are determined by it->orientation[0] (the face_orientation + // flag). For d==3 we have to use GridTools::OrientationLookupTable. // The result is given below. unsigned int p4est_orientation = 0; if (dim == 2) - p4est_orientation = face_pair.orientation[1]; + p4est_orientation = face_pair.orientation[0] == true ? 0u : 1u; else { const unsigned int face_idx_list[] = {face_left, face_right}; diff --git a/source/dofs/dof_tools_constraints.cc b/source/dofs/dof_tools_constraints.cc index 9b0066cfc4..392eeff2f3 100644 --- a/source/dofs/dof_tools_constraints.cc +++ b/source/dofs/dof_tools_constraints.cc @@ -2355,7 +2355,7 @@ namespace DoFTools "(face_orientation, face_flip, face_rotation) " "is invalid for 1d")); - Assert((dim != 2) || (face_orientation == true && face_rotation == false), + Assert((dim != 2) || (face_flip == false && face_rotation == false), ExcMessage("The supplied orientation " "(face_orientation, face_flip, face_rotation) " "is invalid for 2d")); diff --git a/source/fe/fe_q_base.cc b/source/fe/fe_q_base.cc index da6d641e0a..890cd39dfd 100644 --- a/source/fe/fe_q_base.cc +++ b/source/fe/fe_q_base.cc @@ -1170,9 +1170,9 @@ FE_Q_Base::face_to_cell_index(const unsigned int face_index, break; case 2: - // in 2d, only face_flip has a meaning. if it is set, consider - // dofs in reverse order - if (face_flip == false) + // in 2d, only face_orientation has a meaning. if it is false (i.e., + // the non-default case), then consider dofs in reverse order + if (face_orientation == true) adjusted_dof_index_on_line = dof_index_on_line; else adjusted_dof_index_on_line = diff --git a/source/grid/grid_tools_dof_handlers.cc b/source/grid/grid_tools_dof_handlers.cc index 1b4308f301..ca9232c1d7 100644 --- a/source/grid/grid_tools_dof_handlers.cc +++ b/source/grid/grid_tools_dof_handlers.cc @@ -2472,17 +2472,9 @@ namespace GridTools // value so we have to do an additional translation step else if (dim == 2) { - // In 2D, calls in set_periodicity_constraints() ultimately require - // calling FiniteElement::face_to_cell_index(), which in turn (for - // hypercubes) calls GeometryInfo<2>::child_cell_on_face(). The - // final function assumes that orientation in 2D is encoded solely - // in face_flip (the second bit) whereas the orientation bit is - // ignored. Hence, the backwards orientation is 1 + 2 and the - // standard orientation is 1 + 0. - constexpr std::array translation{{3, 1}}; - AssertIndexRange(combined_orientation, translation.size()); - orientation = - translation[std::min(combined_orientation, 1u)]; + // In 2D only the first bit (orientation) is set + AssertIndexRange(combined_orientation, 2); + orientation = combined_orientation; } else { diff --git a/source/grid/tria.cc b/source/grid/tria.cc index 9d20313b0a..3db640cb6e 100644 --- a/source/grid/tria.cc +++ b/source/grid/tria.cc @@ -1665,7 +1665,7 @@ namespace "(face_orientation, face_flip, face_rotation) " "is invalid for 1d")); - Assert((dim != 2) || (face_orientation == true && face_rotation == false), + Assert((dim != 2) || (face_flip == false && face_rotation == false), ExcMessage("The supplied orientation " "(face_orientation, face_flip, face_rotation) " "is invalid for 2d")); @@ -1731,10 +1731,10 @@ namespace // see Documentation of GeometryInfo for details static const int lookup_table_2d[2][2] = - // flip: + // orientation: { - {0, 1}, // false - {1, 0} // true + {1, 0}, // false + {0, 1} // true }; static const int lookup_table_3d[2][2][2][4] = @@ -1779,7 +1779,7 @@ namespace switch (dim) { case 2: - j = lookup_table_2d[face_flip][i]; + j = lookup_table_2d[face_orientation][i]; break; case 3: j = lookup_table_3d[face_orientation][face_flip] diff --git a/tests/base/geometry_info_2.output b/tests/base/geometry_info_2.output index 31f4bf1a3a..1f6ed9947f 100644 --- a/tests/base/geometry_info_2.output +++ b/tests/base/geometry_info_2.output @@ -39,21 +39,21 @@ DEAL:2d::face normal1 +x0 DEAL:2d::face normal2 -x1 DEAL:2d::face normal3 +x1 DEAL:2d::face_children0[true ] 0 2 -DEAL:2d::face_children0[false] 0 2 +DEAL:2d::face_children0[false] 2 0 DEAL:2d::face_children1[true ] 1 3 -DEAL:2d::face_children1[false] 1 3 +DEAL:2d::face_children1[false] 3 1 DEAL:2d::face_children2[true ] 0 1 -DEAL:2d::face_children2[false] 0 1 +DEAL:2d::face_children2[false] 1 0 DEAL:2d::face_children3[true ] 2 3 -DEAL:2d::face_children3[false] 2 3 +DEAL:2d::face_children3[false] 3 2 DEAL:2d::face_vertices0[true ] 0 2 -DEAL:2d::face_vertices0[false] 0 2 +DEAL:2d::face_vertices0[false] 2 0 DEAL:2d::face_vertices1[true ] 1 3 -DEAL:2d::face_vertices1[false] 1 3 +DEAL:2d::face_vertices1[false] 3 1 DEAL:2d::face_vertices2[true ] 0 1 -DEAL:2d::face_vertices2[false] 0 1 +DEAL:2d::face_vertices2[false] 1 0 DEAL:2d::face_vertices3[true ] 2 3 -DEAL:2d::face_vertices3[false] 2 3 +DEAL:2d::face_vertices3[false] 3 2 DEAL:2d::face_lines0[true ] 0 DEAL:2d::face_lines0[false] 0 DEAL:2d::face_lines1[true ] 1 diff --git a/tests/base/geometry_info_8.output b/tests/base/geometry_info_8.output index 6fb781cb36..49b4c57f66 100644 --- a/tests/base/geometry_info_8.output +++ b/tests/base/geometry_info_8.output @@ -39,72 +39,72 @@ DEAL:: (0 -> 1 ) DEAL::2D: DEAL::face 0: DEAL::orientation 0, flip 0, rotation 0: -DEAL:: (0 -> 0 ) (1 -> 2 ) +DEAL:: (0 -> 2 ) (1 -> 0 ) DEAL::orientation 1, flip 0, rotation 0: DEAL:: (0 -> 0 ) (1 -> 2 ) DEAL::orientation 0, flip 1, rotation 0: DEAL:: (0 -> 2 ) (1 -> 0 ) DEAL::orientation 1, flip 1, rotation 0: -DEAL:: (0 -> 2 ) (1 -> 0 ) -DEAL::orientation 0, flip 0, rotation 1: DEAL:: (0 -> 0 ) (1 -> 2 ) +DEAL::orientation 0, flip 0, rotation 1: +DEAL:: (0 -> 2 ) (1 -> 0 ) DEAL::orientation 1, flip 0, rotation 1: DEAL:: (0 -> 0 ) (1 -> 2 ) DEAL::orientation 0, flip 1, rotation 1: DEAL:: (0 -> 2 ) (1 -> 0 ) DEAL::orientation 1, flip 1, rotation 1: -DEAL:: (0 -> 2 ) (1 -> 0 ) +DEAL:: (0 -> 0 ) (1 -> 2 ) DEAL::face 1: DEAL::orientation 0, flip 0, rotation 0: -DEAL:: (0 -> 1 ) (1 -> 3 ) +DEAL:: (0 -> 3 ) (1 -> 1 ) DEAL::orientation 1, flip 0, rotation 0: DEAL:: (0 -> 1 ) (1 -> 3 ) DEAL::orientation 0, flip 1, rotation 0: DEAL:: (0 -> 3 ) (1 -> 1 ) DEAL::orientation 1, flip 1, rotation 0: -DEAL:: (0 -> 3 ) (1 -> 1 ) -DEAL::orientation 0, flip 0, rotation 1: DEAL:: (0 -> 1 ) (1 -> 3 ) +DEAL::orientation 0, flip 0, rotation 1: +DEAL:: (0 -> 3 ) (1 -> 1 ) DEAL::orientation 1, flip 0, rotation 1: DEAL:: (0 -> 1 ) (1 -> 3 ) DEAL::orientation 0, flip 1, rotation 1: DEAL:: (0 -> 3 ) (1 -> 1 ) DEAL::orientation 1, flip 1, rotation 1: -DEAL:: (0 -> 3 ) (1 -> 1 ) +DEAL:: (0 -> 1 ) (1 -> 3 ) DEAL::face 2: DEAL::orientation 0, flip 0, rotation 0: -DEAL:: (0 -> 0 ) (1 -> 1 ) +DEAL:: (0 -> 1 ) (1 -> 0 ) DEAL::orientation 1, flip 0, rotation 0: DEAL:: (0 -> 0 ) (1 -> 1 ) DEAL::orientation 0, flip 1, rotation 0: DEAL:: (0 -> 1 ) (1 -> 0 ) DEAL::orientation 1, flip 1, rotation 0: -DEAL:: (0 -> 1 ) (1 -> 0 ) -DEAL::orientation 0, flip 0, rotation 1: DEAL:: (0 -> 0 ) (1 -> 1 ) +DEAL::orientation 0, flip 0, rotation 1: +DEAL:: (0 -> 1 ) (1 -> 0 ) DEAL::orientation 1, flip 0, rotation 1: DEAL:: (0 -> 0 ) (1 -> 1 ) DEAL::orientation 0, flip 1, rotation 1: DEAL:: (0 -> 1 ) (1 -> 0 ) DEAL::orientation 1, flip 1, rotation 1: -DEAL:: (0 -> 1 ) (1 -> 0 ) +DEAL:: (0 -> 0 ) (1 -> 1 ) DEAL::face 3: DEAL::orientation 0, flip 0, rotation 0: -DEAL:: (0 -> 2 ) (1 -> 3 ) +DEAL:: (0 -> 3 ) (1 -> 2 ) DEAL::orientation 1, flip 0, rotation 0: DEAL:: (0 -> 2 ) (1 -> 3 ) DEAL::orientation 0, flip 1, rotation 0: DEAL:: (0 -> 3 ) (1 -> 2 ) DEAL::orientation 1, flip 1, rotation 0: -DEAL:: (0 -> 3 ) (1 -> 2 ) -DEAL::orientation 0, flip 0, rotation 1: DEAL:: (0 -> 2 ) (1 -> 3 ) +DEAL::orientation 0, flip 0, rotation 1: +DEAL:: (0 -> 3 ) (1 -> 2 ) DEAL::orientation 1, flip 0, rotation 1: DEAL:: (0 -> 2 ) (1 -> 3 ) DEAL::orientation 0, flip 1, rotation 1: DEAL:: (0 -> 3 ) (1 -> 2 ) DEAL::orientation 1, flip 1, rotation 1: -DEAL:: (0 -> 3 ) (1 -> 2 ) +DEAL:: (0 -> 2 ) (1 -> 3 ) DEAL::3D: DEAL::face 0: DEAL::orientation 0, flip 0, rotation 0: diff --git a/tests/dofs/dof_tools_21_b.output b/tests/dofs/dof_tools_21_b.output index d4e21c2cdb..4e6366d90f 100644 --- a/tests/dofs/dof_tools_21_b.output +++ b/tests/dofs/dof_tools_21_b.output @@ -18,7 +18,7 @@ DEAL::DoFs of face_1: DEAL:: component 0: (4 - 1.000 3.000) (5 - -1.000 3.000) DEAL::DoFs of face_2: DEAL:: component 0: (0 - -1.000 -3.000) (1 - 1.000 -3.000) -DEAL::Orientation: 110 +DEAL::Orientation: 000 DEAL::Matching: 4 1: 1.000 5 0: 1.000 @@ -44,7 +44,7 @@ DEAL::DoFs of face_1: DEAL:: component 0: (4 - 1.000 3.000 0.000) (5 - -1.000 3.000 0.000) DEAL::DoFs of face_2: DEAL:: component 0: (0 - -1.000 -3.000 0.000) (1 - 1.000 -3.000 0.000) -DEAL::Orientation: 110 +DEAL::Orientation: 000 DEAL::Matching: 4 1: 1.000 5 0: 1.000 diff --git a/tests/dofs/dof_tools_21_b_x.cc b/tests/dofs/dof_tools_21_b_x.cc index 558718be1b..19693a4779 100644 --- a/tests/dofs/dof_tools_21_b_x.cc +++ b/tests/dofs/dof_tools_21_b_x.cc @@ -198,8 +198,8 @@ print_matching(DoFHandler &dof_handler) std::bitset<3> orientation; - orientation[0] = 1; - orientation[1] = 1; + orientation[0] = 0; + orientation[1] = 0; orientation[2] = 0; DoFTools::make_periodicity_constraints(face_1, diff --git a/tests/dofs/dof_tools_21_b_x_q3.cc b/tests/dofs/dof_tools_21_b_x_q3.cc index 841d489b16..60e3c9a65d 100644 --- a/tests/dofs/dof_tools_21_b_x_q3.cc +++ b/tests/dofs/dof_tools_21_b_x_q3.cc @@ -232,8 +232,8 @@ print_matching(DoFHandler &dof_handler) std::bitset<3> orientation; - orientation[0] = 1; - orientation[1] = 1; + orientation[0] = 0; + orientation[1] = 0; orientation[2] = 0; DoFTools::make_periodicity_constraints(face_1, diff --git a/tests/dofs/dof_tools_21_b_y.cc b/tests/dofs/dof_tools_21_b_y.cc index 6e68157f16..3bd1368a8d 100644 --- a/tests/dofs/dof_tools_21_b_y.cc +++ b/tests/dofs/dof_tools_21_b_y.cc @@ -198,8 +198,8 @@ print_matching(DoFHandler &dof_handler) std::bitset<3> orientation; - orientation[0] = 1; - orientation[1] = 1; + orientation[0] = 0; + orientation[1] = 0; orientation[2] = 0; DoFTools::make_periodicity_constraints(face_1, diff --git a/tests/dofs/dof_tools_21_c.output b/tests/dofs/dof_tools_21_c.output index badceef9da..02a9d34795 100644 --- a/tests/dofs/dof_tools_21_c.output +++ b/tests/dofs/dof_tools_21_c.output @@ -22,7 +22,7 @@ DEAL::DoFs of face_1: DEAL:: component 0: (0 - 1.000 3.000) (1 - -1.000 3.000) DEAL::DoFs of face_2: DEAL:: component 0: (4 - -1.000 -3.000) (8 - 1.000 -3.000) -DEAL::Orientation: 110 +DEAL::Orientation: 000 DEAL::Matching: 4 1: 1.000 5 1: 0.5000 @@ -56,7 +56,7 @@ DEAL::DoFs of face_1: DEAL:: component 0: (0 - 1.000 3.000 0.000) (1 - -1.000 3.000 0.000) DEAL::DoFs of face_2: DEAL:: component 0: (4 - -1.000 -3.000 0.000) (8 - 1.000 -3.000 0.000) -DEAL::Orientation: 110 +DEAL::Orientation: 000 DEAL::Matching: 4 1: 1.000 5 1: 0.5000 diff --git a/tests/fe/face_to_cell_q1_2d.cc b/tests/fe/face_to_cell_q1_2d.cc index 9cd5da9b93..d982f0b45b 100644 --- a/tests/fe/face_to_cell_q1_2d.cc +++ b/tests/fe/face_to_cell_q1_2d.cc @@ -14,8 +14,8 @@ // --------------------------------------------------------------------- -// it turns out that FE_Q::face_to_cell_index() had a bug for elements beyond -// Q2 when using the face flip flag. this test is for the 2d case for the Q1 +// it turns out that FE_Q::face_to_cell_index() had a bug for elements beyond Q2 +// when using the face orientation flag. this test is for the 2d case for the Q1 // case #include @@ -35,13 +35,17 @@ test() { deallog << "Face=" << face << std::endl; - for (int flip = 0; flip < 2; ++flip) + for (int orientation = 0; orientation < 2; ++orientation) { - deallog << " flip=" << (flip == 0 ? "false" : "true") << std::endl + deallog << " orientation=" << (orientation == 0 ? "false" : "true") + << std::endl << " "; for (unsigned int i = 0; i < dofs_per_face; ++i) - deallog << fe.face_to_cell_index( - i, face, true, (flip == 0 ? false : true), false) + deallog << fe.face_to_cell_index(i, + face, + (orientation == 0 ? false : true), + false, + false) << " - "; deallog << std::endl; } diff --git a/tests/fe/face_to_cell_q1_2d.output b/tests/fe/face_to_cell_q1_2d.output index af84b631db..3f8be74a0b 100644 --- a/tests/fe/face_to_cell_q1_2d.output +++ b/tests/fe/face_to_cell_q1_2d.output @@ -1,21 +1,21 @@ DEAL::Face=0 -DEAL:: flip=false -DEAL:: 0 - 2 - -DEAL:: flip=true +DEAL:: orientation=false DEAL:: 2 - 0 - +DEAL:: orientation=true +DEAL:: 0 - 2 - DEAL::Face=1 -DEAL:: flip=false -DEAL:: 1 - 3 - -DEAL:: flip=true +DEAL:: orientation=false DEAL:: 3 - 1 - +DEAL:: orientation=true +DEAL:: 1 - 3 - DEAL::Face=2 -DEAL:: flip=false -DEAL:: 0 - 1 - -DEAL:: flip=true +DEAL:: orientation=false DEAL:: 1 - 0 - +DEAL:: orientation=true +DEAL:: 0 - 1 - DEAL::Face=3 -DEAL:: flip=false -DEAL:: 2 - 3 - -DEAL:: flip=true +DEAL:: orientation=false DEAL:: 3 - 2 - +DEAL:: orientation=true +DEAL:: 2 - 3 - diff --git a/tests/fe/face_to_cell_q2_2d.cc b/tests/fe/face_to_cell_q2_2d.cc index 491bad8c68..4ae0330b4e 100644 --- a/tests/fe/face_to_cell_q2_2d.cc +++ b/tests/fe/face_to_cell_q2_2d.cc @@ -15,8 +15,8 @@ // it turns out that FE_Q::face_to_cell_index() had a bug for elements beyond -// Q2 when using the face flip flag. this test is for the 2d case for the Q2 -// case +// Q2 when using the face orientation flag. this test is for the 2d case for the +// Q2 case #include @@ -35,13 +35,17 @@ test() { deallog << "Face=" << face << std::endl; - for (int flip = 0; flip < 2; ++flip) + for (int orientation = 0; orientation < 2; ++orientation) { - deallog << " flip=" << (flip == 0 ? "false" : "true") << std::endl + deallog << " orientation=" << (orientation == 0 ? "false" : "true") + << std::endl << " "; for (unsigned int i = 0; i < dofs_per_face; ++i) - deallog << fe.face_to_cell_index( - i, face, true, (flip == 0 ? false : true), false) + deallog << fe.face_to_cell_index(i, + face, + (orientation == 0 ? false : true), + false, + false) << " - "; deallog << std::endl; } diff --git a/tests/fe/face_to_cell_q2_2d.output b/tests/fe/face_to_cell_q2_2d.output index fade4e0e77..fb9b682074 100644 --- a/tests/fe/face_to_cell_q2_2d.output +++ b/tests/fe/face_to_cell_q2_2d.output @@ -1,21 +1,21 @@ DEAL::Face=0 -DEAL:: flip=false -DEAL:: 0 - 2 - 4 - -DEAL:: flip=true +DEAL:: orientation=false DEAL:: 2 - 0 - 4 - +DEAL:: orientation=true +DEAL:: 0 - 2 - 4 - DEAL::Face=1 -DEAL:: flip=false -DEAL:: 1 - 3 - 5 - -DEAL:: flip=true +DEAL:: orientation=false DEAL:: 3 - 1 - 5 - +DEAL:: orientation=true +DEAL:: 1 - 3 - 5 - DEAL::Face=2 -DEAL:: flip=false -DEAL:: 0 - 1 - 6 - -DEAL:: flip=true +DEAL:: orientation=false DEAL:: 1 - 0 - 6 - +DEAL:: orientation=true +DEAL:: 0 - 1 - 6 - DEAL::Face=3 -DEAL:: flip=false -DEAL:: 2 - 3 - 7 - -DEAL:: flip=true +DEAL:: orientation=false DEAL:: 3 - 2 - 7 - +DEAL:: orientation=true +DEAL:: 2 - 3 - 7 - diff --git a/tests/fe/face_to_cell_q2xq2_2d.cc b/tests/fe/face_to_cell_q2xq2_2d.cc index cd6ff7ad9f..2dbb738228 100644 --- a/tests/fe/face_to_cell_q2xq2_2d.cc +++ b/tests/fe/face_to_cell_q2xq2_2d.cc @@ -15,8 +15,8 @@ // it turns out that FE_Q::face_to_cell_index() had a bug for elements beyond -// Q2XQ2 when using the face flip flag. this test is for the 2d case for the -// Q2XQ2 case +// Q2XQ2 when using the face orientation flag. this test is for the 2d case for +// the Q2XQ2 case #include #include @@ -36,13 +36,17 @@ test() { deallog << "Face=" << face << std::endl; - for (int flip = 0; flip < 2; ++flip) + for (int orientation = 0; orientation < 2; ++orientation) { - deallog << " flip=" << (flip == 0 ? "false" : "true") << std::endl + deallog << " orientation=" << (orientation == 0 ? "false" : "true") + << std::endl << " "; for (unsigned int i = 0; i < dofs_per_face; ++i) - deallog << fe.face_to_cell_index( - i, face, true, (flip == 0 ? false : true), false) + deallog << fe.face_to_cell_index(i, + face, + (orientation == 0 ? false : true), + false, + false) << " - "; deallog << std::endl; } diff --git a/tests/fe/face_to_cell_q2xq2_2d.output b/tests/fe/face_to_cell_q2xq2_2d.output index 133d8fde6c..f4097117a8 100644 --- a/tests/fe/face_to_cell_q2xq2_2d.output +++ b/tests/fe/face_to_cell_q2xq2_2d.output @@ -1,21 +1,21 @@ DEAL::Face=0 -DEAL:: flip=false -DEAL:: 0 - 1 - 4 - 5 - 8 - 9 - -DEAL:: flip=true +DEAL:: orientation=false DEAL:: 4 - 5 - 0 - 1 - 8 - 9 - +DEAL:: orientation=true +DEAL:: 0 - 1 - 4 - 5 - 8 - 9 - DEAL::Face=1 -DEAL:: flip=false -DEAL:: 2 - 3 - 6 - 7 - 10 - 11 - -DEAL:: flip=true +DEAL:: orientation=false DEAL:: 6 - 7 - 2 - 3 - 10 - 11 - +DEAL:: orientation=true +DEAL:: 2 - 3 - 6 - 7 - 10 - 11 - DEAL::Face=2 -DEAL:: flip=false -DEAL:: 0 - 1 - 2 - 3 - 12 - 13 - -DEAL:: flip=true +DEAL:: orientation=false DEAL:: 2 - 3 - 0 - 1 - 12 - 13 - +DEAL:: orientation=true +DEAL:: 0 - 1 - 2 - 3 - 12 - 13 - DEAL::Face=3 -DEAL:: flip=false -DEAL:: 4 - 5 - 6 - 7 - 14 - 15 - -DEAL:: flip=true +DEAL:: orientation=false DEAL:: 6 - 7 - 4 - 5 - 14 - 15 - +DEAL:: orientation=true +DEAL:: 4 - 5 - 6 - 7 - 14 - 15 - diff --git a/tests/fe/face_to_cell_q3_2d.cc b/tests/fe/face_to_cell_q3_2d.cc index 71392867ab..4357d1f31a 100644 --- a/tests/fe/face_to_cell_q3_2d.cc +++ b/tests/fe/face_to_cell_q3_2d.cc @@ -15,8 +15,8 @@ // it turns out that FE_Q::face_to_cell_index() had a bug for elements beyond -// Q2 when using the face flip flag. this test is for the 2d case for the Q3 -// case +// Q2 when using the face orientation flag. this test is for the 2d case for the +// Q3 case #include @@ -35,13 +35,17 @@ test() { deallog << "Face=" << face << std::endl; - for (int flip = 0; flip < 2; ++flip) + for (int orientation = 0; orientation < 2; ++orientation) { - deallog << " flip=" << (flip == 0 ? "false" : "true") << std::endl + deallog << " orientation=" << (orientation == 0 ? "false" : "true") + << std::endl << " "; for (unsigned int i = 0; i < dofs_per_face; ++i) - deallog << fe.face_to_cell_index( - i, face, true, (flip == 0 ? false : true), false) + deallog << fe.face_to_cell_index(i, + face, + (orientation == 0 ? false : true), + false, + false) << " - "; deallog << std::endl; } diff --git a/tests/fe/face_to_cell_q3_2d.output b/tests/fe/face_to_cell_q3_2d.output index ade889941d..d4d07f7880 100644 --- a/tests/fe/face_to_cell_q3_2d.output +++ b/tests/fe/face_to_cell_q3_2d.output @@ -1,21 +1,21 @@ DEAL::Face=0 -DEAL:: flip=false -DEAL:: 0 - 2 - 4 - 5 - -DEAL:: flip=true +DEAL:: orientation=false DEAL:: 2 - 0 - 5 - 4 - +DEAL:: orientation=true +DEAL:: 0 - 2 - 4 - 5 - DEAL::Face=1 -DEAL:: flip=false -DEAL:: 1 - 3 - 6 - 7 - -DEAL:: flip=true +DEAL:: orientation=false DEAL:: 3 - 1 - 7 - 6 - +DEAL:: orientation=true +DEAL:: 1 - 3 - 6 - 7 - DEAL::Face=2 -DEAL:: flip=false -DEAL:: 0 - 1 - 8 - 9 - -DEAL:: flip=true +DEAL:: orientation=false DEAL:: 1 - 0 - 9 - 8 - +DEAL:: orientation=true +DEAL:: 0 - 1 - 8 - 9 - DEAL::Face=3 -DEAL:: flip=false -DEAL:: 2 - 3 - 10 - 11 - -DEAL:: flip=true +DEAL:: orientation=false DEAL:: 3 - 2 - 11 - 10 - +DEAL:: orientation=true +DEAL:: 2 - 3 - 10 - 11 - diff --git a/tests/fe/face_to_cell_q3xq4_2d.cc b/tests/fe/face_to_cell_q3xq4_2d.cc index 910b0e5ce8..562934aaec 100644 --- a/tests/fe/face_to_cell_q3xq4_2d.cc +++ b/tests/fe/face_to_cell_q3xq4_2d.cc @@ -15,8 +15,8 @@ // it turns out that FE_Q::face_to_cell_index() had a bug for elements beyond -// Q2 when using the face flip flag. this test is for the 2d case for the Q3XQ4 -// case +// Q2 when using the face orientation flag. this test is for the 2d case for the +// Q3XQ4 case #include #include @@ -36,13 +36,17 @@ test() { deallog << "Face=" << face << std::endl; - for (int flip = 0; flip < 2; ++flip) + for (int orientation = 0; orientation < 2; ++orientation) { - deallog << " flip=" << (flip == 0 ? "false" : "true") << std::endl + deallog << " orientation=" << (orientation == 0 ? "false" : "true") + << std::endl << " "; for (unsigned int i = 0; i < dofs_per_face; ++i) - deallog << fe.face_to_cell_index( - i, face, true, (flip == 0 ? false : true), false) + deallog << fe.face_to_cell_index(i, + face, + (orientation == 0 ? false : true), + false, + false) << " - "; deallog << std::endl; } diff --git a/tests/fe/face_to_cell_q3xq4_2d.output b/tests/fe/face_to_cell_q3xq4_2d.output index fa15ff4d1d..c9a2f73134 100644 --- a/tests/fe/face_to_cell_q3xq4_2d.output +++ b/tests/fe/face_to_cell_q3xq4_2d.output @@ -1,21 +1,21 @@ DEAL::Face=0 -DEAL:: flip=false -DEAL:: 0 - 1 - 4 - 5 - 8 - 9 - 10 - 11 - 12 - -DEAL:: flip=true +DEAL:: orientation=false DEAL:: 4 - 5 - 0 - 1 - 9 - 8 - 12 - 11 - 10 - +DEAL:: orientation=true +DEAL:: 0 - 1 - 4 - 5 - 8 - 9 - 10 - 11 - 12 - DEAL::Face=1 -DEAL:: flip=false -DEAL:: 2 - 3 - 6 - 7 - 13 - 14 - 15 - 16 - 17 - -DEAL:: flip=true +DEAL:: orientation=false DEAL:: 6 - 7 - 2 - 3 - 14 - 13 - 17 - 16 - 15 - +DEAL:: orientation=true +DEAL:: 2 - 3 - 6 - 7 - 13 - 14 - 15 - 16 - 17 - DEAL::Face=2 -DEAL:: flip=false -DEAL:: 0 - 1 - 2 - 3 - 18 - 19 - 20 - 21 - 22 - -DEAL:: flip=true +DEAL:: orientation=false DEAL:: 2 - 3 - 0 - 1 - 19 - 18 - 22 - 21 - 20 - +DEAL:: orientation=true +DEAL:: 0 - 1 - 2 - 3 - 18 - 19 - 20 - 21 - 22 - DEAL::Face=3 -DEAL:: flip=false -DEAL:: 4 - 5 - 6 - 7 - 23 - 24 - 25 - 26 - 27 - -DEAL:: flip=true +DEAL:: orientation=false DEAL:: 6 - 7 - 4 - 5 - 24 - 23 - 27 - 26 - 25 - +DEAL:: orientation=true +DEAL:: 4 - 5 - 6 - 7 - 23 - 24 - 25 - 26 - 27 - diff --git a/tests/fe/face_to_cell_q4_2d.cc b/tests/fe/face_to_cell_q4_2d.cc index a2b5c95093..3bdde4eada 100644 --- a/tests/fe/face_to_cell_q4_2d.cc +++ b/tests/fe/face_to_cell_q4_2d.cc @@ -15,8 +15,8 @@ // it turns out that FE_Q::face_to_cell_index() had a bug for elements beyond -// Q2 when using the face flip flag. this test is for the 2d case for the Q4 -// case +// Q2 when using the face orientation flag. this test is for the 2d case for the +// Q4 case #include @@ -35,13 +35,17 @@ test() { deallog << "Face=" << face << std::endl; - for (int flip = 0; flip < 2; ++flip) + for (int orientation = 0; orientation < 2; ++orientation) { - deallog << " flip=" << (flip == 0 ? "false" : "true") << std::endl + deallog << " orientation=" << (orientation == 0 ? "false" : "true") + << std::endl << " "; for (unsigned int i = 0; i < dofs_per_face; ++i) - deallog << fe.face_to_cell_index( - i, face, true, (flip == 0 ? false : true), false) + deallog << fe.face_to_cell_index(i, + face, + (orientation == 0 ? false : true), + false, + false) << " - "; deallog << std::endl; } diff --git a/tests/fe/face_to_cell_q4_2d.output b/tests/fe/face_to_cell_q4_2d.output index 3cc09cb067..31fcf4a105 100644 --- a/tests/fe/face_to_cell_q4_2d.output +++ b/tests/fe/face_to_cell_q4_2d.output @@ -1,21 +1,21 @@ DEAL::Face=0 -DEAL:: flip=false -DEAL:: 0 - 2 - 4 - 5 - 6 - -DEAL:: flip=true +DEAL:: orientation=false DEAL:: 2 - 0 - 6 - 5 - 4 - +DEAL:: orientation=true +DEAL:: 0 - 2 - 4 - 5 - 6 - DEAL::Face=1 -DEAL:: flip=false -DEAL:: 1 - 3 - 7 - 8 - 9 - -DEAL:: flip=true +DEAL:: orientation=false DEAL:: 3 - 1 - 9 - 8 - 7 - +DEAL:: orientation=true +DEAL:: 1 - 3 - 7 - 8 - 9 - DEAL::Face=2 -DEAL:: flip=false -DEAL:: 0 - 1 - 10 - 11 - 12 - -DEAL:: flip=true +DEAL:: orientation=false DEAL:: 1 - 0 - 12 - 11 - 10 - +DEAL:: orientation=true +DEAL:: 0 - 1 - 10 - 11 - 12 - DEAL::Face=3 -DEAL:: flip=false -DEAL:: 2 - 3 - 13 - 14 - 15 - -DEAL:: flip=true +DEAL:: orientation=false DEAL:: 3 - 2 - 15 - 14 - 13 - +DEAL:: orientation=true +DEAL:: 2 - 3 - 13 - 14 - 15 - diff --git a/tests/grid/grid_tools_06.output b/tests/grid/grid_tools_06.output index 1fde49825d..9f16c53b44 100644 --- a/tests/grid/grid_tools_06.output +++ b/tests/grid/grid_tools_06.output @@ -9,7 +9,7 @@ DEAL:: DEAL::Triangulation: 1 DEAL::face 1 :: 1.000 3.000 :: -1.000 3.000 DEAL::face 2 :: -1.000 -3.000 :: 1.000 -3.000 -DEAL::orientation: 1 flip: 1 rotation: 0 +DEAL::orientation: 0 flip: 0 rotation: 0 DEAL:: DEAL::Test for 3D: Hypercube DEAL:: -- 2.39.5