]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Use face_orientation to encode 2D orientation information. 16373/head
authorDavid Wells <drwells@email.unc.edu>
Wed, 20 Dec 2023 14:49:41 +0000 (09:49 -0500)
committerDavid Wells <drwells@email.unc.edu>
Wed, 20 Dec 2023 22:03:38 +0000 (17:03 -0500)
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.

29 files changed:
doc/news/changes/incompatibilities/20231220DavidWells [new file with mode: 0644]
include/deal.II/base/geometry_info.h
include/deal.II/grid/reference_cell.h
include/deal.II/grid/tria_accessor.templates.h
source/distributed/tria.cc
source/dofs/dof_tools_constraints.cc
source/fe/fe_q_base.cc
source/grid/grid_tools_dof_handlers.cc
source/grid/tria.cc
tests/base/geometry_info_2.output
tests/base/geometry_info_8.output
tests/dofs/dof_tools_21_b.output
tests/dofs/dof_tools_21_b_x.cc
tests/dofs/dof_tools_21_b_x_q3.cc
tests/dofs/dof_tools_21_b_y.cc
tests/dofs/dof_tools_21_c.output
tests/fe/face_to_cell_q1_2d.cc
tests/fe/face_to_cell_q1_2d.output
tests/fe/face_to_cell_q2_2d.cc
tests/fe/face_to_cell_q2_2d.output
tests/fe/face_to_cell_q2xq2_2d.cc
tests/fe/face_to_cell_q2xq2_2d.output
tests/fe/face_to_cell_q3_2d.cc
tests/fe/face_to_cell_q3_2d.output
tests/fe/face_to_cell_q3xq4_2d.cc
tests/fe/face_to_cell_q3xq4_2d.output
tests/fe/face_to_cell_q4_2d.cc
tests/fe/face_to_cell_q4_2d.output
tests/grid/grid_tools_06.output

diff --git a/doc/news/changes/incompatibilities/20231220DavidWells b/doc/news/changes/incompatibilities/20231220DavidWells
new file mode 100644 (file)
index 0000000..2ad2e24
--- /dev/null
@@ -0,0 +1,4 @@
+Changed: Meshes in 2D now use the face_orientation boolean (instead of the
+face_flip boolean) to represent their orientation.
+<br>
+(David Wells, 2023/12/20)
index 9b7f66d90759aa02588325f1b542e7e1250f4a2c..715adfaa002fc79da87b71b5b6f64adac29eef11 100644 (file)
@@ -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];
 }
 
 
index f7ab3ec133ef62f5127227c29a44d5a4f9fc0ba4..dbde9f41213b235b7ee18b9575f5163104e889e2 100644 (file)
@@ -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)
     {
index cd8f4d0ac107ba4e1f3f39669af5ad3133a8f04f..ddac77f8c7c5c4775526ab2403b231d9d77593bf 100644 (file)
@@ -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;
       }
index adc93a90c692e5264b37380edd4771de1756e15a..7368fc0dbb49b2a742032b5596540298e4541de7 100644 (file)
@@ -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};
index 9b0066cfc405c75c8497ec3736cb078691c9740d..392eeff2f34ad044b852c74fb581e5ee82690f23 100644 (file)
@@ -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"));
index da6d641e0a96633cae19a541bcd49c3a15bfa9bc..890cd39dfd952884626249a8f1eef944ebeea469 100644 (file)
@@ -1170,9 +1170,9 @@ FE_Q_Base<dim, spacedim>::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 =
index 1b4308f30171cacb645df1f1a35155ac0a8fa170..ca9232c1d7a7c874c0818326c3c539a6ca865500 100644 (file)
@@ -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<unsigned int, 2> translation{{3, 1}};
-            AssertIndexRange(combined_orientation, translation.size());
-            orientation =
-              translation[std::min<unsigned int>(combined_orientation, 1u)];
+            // In 2D only the first bit (orientation) is set
+            AssertIndexRange(combined_orientation, 2);
+            orientation = combined_orientation;
           }
         else
           {
index 9d20313b0a575e7be6eff274c1b0f1485aa1be46..3db640cb6e06344847ed7b323640430cffbfd537 100644 (file)
@@ -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]
index 31f4bf1a3ac9e2ce62220938713e1891cf3eaffa..1f6ed9947f215d95d417b7d6cf184756eaefa8e5 100644 (file)
@@ -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
index 6fb781cb3632da3334b397f5d1366afa656f1d05..49b4c57f6660871236260dd6340ff75c087f82e6 100644 (file)
@@ -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:
index d4e21c2cdbf8623a810e3a6b6ce3ce6efc644f99..4e6366d90f664ec4c8086d0cee3ea6a3ef7ad842 100644 (file)
@@ -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
index 558718be1bce8a3011facbd03623941d26e0482b..19693a47790248bcdca372ec6300d4a8a0db3650 100644 (file)
@@ -198,8 +198,8 @@ print_matching(DoFHandler<dim, spacedim> &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,
index 841d489b161b0281e2ecf68d1c80c0183815f027..60e3c9a65dac10af585f7ea6ae4bb49dd2d22cd6 100644 (file)
@@ -232,8 +232,8 @@ print_matching(DoFHandler<dim, spacedim> &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,
index 6e68157f16958b87c9cfafc328ed9872a54b2e04..3bd1368a8df520411db8a7680872c4b0e2fa88dd 100644 (file)
@@ -198,8 +198,8 @@ print_matching(DoFHandler<dim, spacedim> &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,
index badceef9dac9a4a75e6bc381b0d1aa64ef58a977..02a9d34795d4aed1831f3d15484667454d4c220f 100644 (file)
@@ -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
index 9cd5da9b93070e96b5ace0d5be002c6b187402b3..d982f0b45bd6df2fed3818a417e35051c3073f79 100644 (file)
@@ -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 <deal.II/fe/fe_q.h>
@@ -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;
         }
index af84b631dba2545021db10d4d99dd185ca7d2f67..3f8be74a0b4b5c434df0ccc4b38f9d75f8fbeafd 100644 (file)
@@ -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 -
index 491bad8c68f52c550350fee8efeadf1375b0017d..4ae0330b4e6e97b55634f4a8855c9c65b712ae6c 100644 (file)
@@ -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 <deal.II/fe/fe_q.h>
 
@@ -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;
         }
index fade4e0e77d80fcb2eef59011ed5ed5a47553939..fb9b682074d09a75a5706ddb7a7b852fe849829d 100644 (file)
@@ -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 -
index cd6ff7ad9f4b7c3969a6db817543114e3fb05b6a..2dbb7382288eaa2b0b71ac0795acf22746800ac5 100644 (file)
@@ -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 <deal.II/fe/fe_q.h>
 #include <deal.II/fe/fe_system.h>
@@ -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;
         }
index 133d8fde6c3b7f253f088750459d4b06568dd483..f4097117a8efdf2c10a95866a1236f6614b67d30 100644 (file)
@@ -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 -
index 71392867ab561548c3838d68da4e7a5965aa4efb..4357d1f31af5a029b7976495d00363584af54a96 100644 (file)
@@ -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 <deal.II/fe/fe_q.h>
 
@@ -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;
         }
index ade889941d013c6b0e63923b8a7b5d9a8d5b8bcd..d4d07f788032aae546c1b1b90af1ff1bbec771cc 100644 (file)
@@ -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 -
index 910b0e5ce85d5776b1d33fe091ed45f8fb70a2d0..562934aaec0715959d2b3399254eb8c7317d1c63 100644 (file)
@@ -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 <deal.II/fe/fe_q.h>
 #include <deal.II/fe/fe_system.h>
@@ -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;
         }
index fa15ff4d1d02d8adcb85df90b4e2463ebe2c366a..c9a2f731348fa251a02b664baa149e506060f8c4 100644 (file)
@@ -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 -
index a2b5c9509378356a03c3b33cc302e444d3dd82e7..3bdde4eadaf76498e0f2a6d9a7de3e405b19a805 100644 (file)
@@ -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 <deal.II/fe/fe_q.h>
 
@@ -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;
         }
index 3cc09cb06785b6e0690cf1845ea04445429bf121..31fcf4a10561d25d22413971968f75735d03d791 100644 (file)
@@ -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 -
index 1fde49825d1b407be3c8839b6fd009b9f898715b..9f16c53b448ba3ac95bf7beb69d0a4a2c852cc89 100644 (file)
@@ -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::

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.