]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Fix orientation
authorPeter Munch <peterrmuench@gmail.com>
Sun, 2 Jul 2023 01:22:27 +0000 (03:22 +0200)
committerPeter Munch <peterrmuench@gmail.com>
Sun, 2 Jul 2023 15:03:21 +0000 (17:03 +0200)
include/deal.II/grid/connectivity.h
include/deal.II/grid/reference_cell.h
tests/simplex/negative_measure_01.output
tests/simplex/orientation_02.cc
tests/simplex/orientation_02.output
tests/simplex/orientation_03.cc [new file with mode: 0644]
tests/simplex/orientation_03.output [new file with mode: 0644]
tests/simplex/orientation_04.cc [new file with mode: 0644]
tests/simplex/orientation_04.output [new file with mode: 0644]
tests/simplex/q_projection_01.output

index 93ebd597cb109a9cbacd390a274a6c944fc19d49..feff81c087131343bcd2b0f6a7518ba6e5550369 100644 (file)
@@ -1168,11 +1168,11 @@ namespace internal
                 offset_i,
                 ad_entity_types[offset_i]
                   .template get_combined_orientation<unsigned int>(
-                    make_array_view(ad_entity_vertices[offset_i].begin(),
-                                    ad_entity_vertices[offset_i].begin() +
-                                      ad_entity_types[offset_i].n_vertices()),
                     make_array_view(ref_indices.begin(),
                                     ref_indices.begin() +
+                                      ad_entity_types[offset_i].n_vertices()),
+                    make_array_view(ad_entity_vertices[offset_i].begin(),
+                                    ad_entity_vertices[offset_i].begin() +
                                       ad_entity_types[offset_i].n_vertices())));
             }
           col_d[offset_i] = counter;
index 65c962fa11ca88b573550308d9d76c06f7947eb7..79c495f2344da80fed746b300f04f634a6870d60 100644 (file)
@@ -2752,11 +2752,11 @@ ReferenceCell::get_combined_orientation(
           return 1;
 
         // face_orientation=true, face_rotation=true, face_flip=false
-        if (v0_equals({vertices_1[1], vertices_1[2], vertices_1[0]}))
+        if (v0_equals({vertices_1[2], vertices_1[0], vertices_1[1]}))
           return 3;
 
         // face_orientation=true, face_rotation=false, face_flip=true
-        if (v0_equals({vertices_1[2], vertices_1[0], vertices_1[1]}))
+        if (v0_equals({vertices_1[1], vertices_1[2], vertices_1[0]}))
           return 5;
 
         // face_orientation=false, face_rotation=false, face_flip=false
@@ -2879,9 +2879,9 @@ ReferenceCell::permute_by_combined_orientation(
             case 1:
               return {vertices[0], vertices[1], vertices[2]};
             case 3:
-              return {vertices[1], vertices[2], vertices[0]};
-            case 5:
               return {vertices[2], vertices[0], vertices[1]};
+            case 5:
+              return {vertices[1], vertices[2], vertices[0]};
             case 0:
               return {vertices[0], vertices[2], vertices[1]};
             case 2:
index 20e21352e65c9df4619e2165c17ea5d3f17e2575..fc30f46757f37ddfcb740c27bbe7d355a21036d5 100644 (file)
@@ -23,7 +23,7 @@ POINTS 10 double
 CELLS 3 21
 8 0 1 3 2 4 5 7 6
 5 4 5 7 6 8
-5 5 1 3 7 9
+5 3 7 5 1 9
 
 CELL_TYPES 3
 12 14 14 
index aa48ff162f642261133f94e47cf863c897342040..3b72f582edd443daad1e599ca0223d882f577b5b 100644 (file)
@@ -47,11 +47,8 @@ test(const unsigned int orientation)
     const auto &face = dummy.begin()->face(face_no);
     const auto  permuted =
       ReferenceCell(ReferenceCells::Triangle)
-        .permute_according_orientation(
-          std::array<unsigned int, 3>{{face->vertex_index(0),
-                                       face->vertex_index(1),
-                                       face->vertex_index(2)}},
-          orientation);
+        .permute_according_orientation(std::array<unsigned int, 3>{{0, 1, 2}},
+                                       orientation);
 
     auto direction =
       cross_product_3d(vertices[permuted[1]] - vertices[permuted[0]],
@@ -61,7 +58,12 @@ test(const unsigned int orientation)
     vertices.push_back(face->center() + direction);
 
     CellData<3> cell;
-    cell.vertices = {permuted[0], permuted[1], permuted[2], 4u};
+    cell.vertices.resize(4);
+
+    cell.vertices[permuted[0]] = face->vertex_index(0);
+    cell.vertices[permuted[1]] = face->vertex_index(1);
+    cell.vertices[permuted[2]] = face->vertex_index(2);
+    cell.vertices[3]           = 4;
 
     cell.material_id = 1;
     cells.push_back(cell);
@@ -73,7 +75,8 @@ test(const unsigned int orientation)
   cell++;
 
   // check orientation
-  deallog << "face orientation: " << orientation << ' ' << std::endl;
+  deallog << "face orientation: " << orientation << ' '
+          << int(cell->combined_face_orientation(0)) << ' ' << std::endl;
   AssertDimension(orientation,
                   (cell->face_orientation(0) * 1 + cell->face_rotation(0) * 2 +
                    cell->face_flip(0) * 4));
index 4291a7ca460cd76be8c34dd65a905e392d4e7984..2d1a155dbe70c44c6fef68f836e0cca371de5870 100644 (file)
@@ -1,35 +1,35 @@
 
-DEAL::face orientation: 0 
+DEAL::face orientation: 0 
 DEAL::0 2 1 4  vs. 0 2 1 4 
 DEAL::0-2 vs. 0-2
 DEAL::1-2 vs. 1-2
 DEAL::0-1 vs. 0-1
 DEAL::
-DEAL::face orientation: 1 
+DEAL::face orientation: 1 
 DEAL::0 1 2 4  vs. 0 1 2 4 
 DEAL::0-1 vs. 0-1
 DEAL::1-2 vs. 1-2
 DEAL::0-2 vs. 0-2
 DEAL::
-DEAL::face orientation: 2 
+DEAL::face orientation: 2 
 DEAL::2 1 0 4  vs. 2 1 0 4 
 DEAL::1-2 vs. 1-2
 DEAL::0-1 vs. 0-1
 DEAL::0-2 vs. 0-2
 DEAL::
-DEAL::face orientation: 3 
+DEAL::face orientation: 3 
 DEAL::1 2 0 4  vs. 1 2 0 4 
 DEAL::1-2 vs. 1-2
 DEAL::0-2 vs. 0-2
 DEAL::0-1 vs. 0-1
 DEAL::
-DEAL::face orientation: 4 
+DEAL::face orientation: 4 
 DEAL::1 0 2 4  vs. 1 0 2 4 
 DEAL::0-1 vs. 0-1
 DEAL::0-2 vs. 0-2
 DEAL::1-2 vs. 1-2
 DEAL::
-DEAL::face orientation: 5 
+DEAL::face orientation: 5 
 DEAL::2 0 1 4  vs. 2 0 1 4 
 DEAL::0-2 vs. 0-2
 DEAL::0-1 vs. 0-1
diff --git a/tests/simplex/orientation_03.cc b/tests/simplex/orientation_03.cc
new file mode 100644 (file)
index 0000000..62bae39
--- /dev/null
@@ -0,0 +1,220 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2023 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+// Verify that the combination of QProjector, FEFaceValues, and various face
+// orientations place face quadrature points in the same locations on abutting
+// cells.
+
+#include <deal.II/fe/fe_nothing.h>
+#include <deal.II/fe/fe_simplex_p.h>
+#include <deal.II/fe/fe_values.h>
+#include <deal.II/fe/mapping.h>
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/tria.h>
+
+#include "../tests.h"
+
+template <int dim>
+void
+test()
+{
+  Triangulation<dim> tria;
+
+  const unsigned int orientation = 3; // 3, 5
+
+  if (false)
+    {
+      unsigned int ref = 2;
+      GridGenerator::subdivided_hyper_cube_with_simplices(
+        tria, Utilities::pow(2, ref), -1.0, +1.0);
+    }
+  else
+    {
+      const unsigned int face_no = 0;
+
+      Triangulation<3> dummy;
+      GridGenerator::reference_cell(dummy, ReferenceCells::Tetrahedron);
+
+      auto vertices = dummy.get_vertices();
+
+      std::vector<CellData<3>> cells;
+
+      {
+        CellData<3> cell;
+        cell.vertices    = {0, 1, 2, 3};
+        cell.material_id = 0;
+        cells.push_back(cell);
+      }
+
+      if (false)
+        {
+          const auto &face = dummy.begin()->face(face_no);
+          const auto  permuted =
+            ReferenceCell(ReferenceCells::Triangle)
+              .permute_according_orientation(
+                std::array<unsigned int, 3>{{face->vertex_index(0),
+                                             face->vertex_index(1),
+                                             face->vertex_index(2)}},
+                orientation);
+
+          for (const auto o : permuted)
+            std::cout << o << " ";
+          std::cout << std::endl;
+
+          auto direction =
+            cross_product_3d(vertices[permuted[1]] - vertices[permuted[0]],
+                             vertices[permuted[2]] - vertices[permuted[0]]);
+          direction = direction / direction.norm();
+
+          std::cout << direction << std::endl;
+
+          vertices.emplace_back(0.0, 0.0, direction[2]);
+
+          CellData<3> cell;
+          cell.vertices = {permuted[0], permuted[1], permuted[2], 4u};
+
+          cell.material_id = 1;
+          cells.push_back(cell);
+        }
+      else
+        {
+          const auto &face = dummy.begin()->face(face_no);
+          const auto  permuted =
+            ReferenceCell(ReferenceCells::Triangle)
+              .permute_according_orientation(
+                std::array<unsigned int, 3>{{0, 1, 2}}, orientation);
+
+          for (const auto o : permuted)
+            std::cout << o << " ";
+          std::cout << std::endl;
+
+          auto direction =
+            cross_product_3d(vertices[permuted[1]] - vertices[permuted[0]],
+                             vertices[permuted[2]] - vertices[permuted[0]]);
+          direction = direction / direction.norm();
+
+          std::cout << direction << std::endl;
+
+          vertices.emplace_back(0.0, 0.0, direction[2]);
+
+          CellData<3> cell;
+          cell.vertices.resize(4);
+
+          cell.vertices[permuted[0]] = face->vertex_index(0);
+          cell.vertices[permuted[1]] = face->vertex_index(1);
+          cell.vertices[permuted[2]] = face->vertex_index(2);
+          cell.vertices[3]           = 4;
+
+          cell.material_id = 1;
+          cells.push_back(cell);
+        }
+
+      tria.create_triangulation(vertices, cells, {});
+
+      for (const auto &cell : tria.active_cell_iterators())
+        {
+          for (const auto l : cell->line_indices())
+            std::cout << cell->line_orientation(l) << " ";
+          std::cout << std::endl;
+        }
+      std::cout << std::endl;
+    }
+
+  const ReferenceCell ref_cell = ReferenceCells::get_simplex<dim>();
+  FE_Nothing<dim>     fe(ref_cell);
+  const Mapping<dim> &mapping =
+    ref_cell.template get_default_linear_mapping<dim>();
+
+  // Quadrature<dim - 1> face_quad =
+  //  ref_cell.face_reference_cell(0).template get_gauss_type_quadrature<dim -
+  //  1>(
+  //    2u);
+
+  Quadrature<dim - 1> face_quad(
+    FE_SimplexP<dim>(1).get_unit_face_support_points());
+
+  FEFaceValues<dim> cell_face_values(mapping,
+                                     fe,
+                                     face_quad,
+                                     update_quadrature_points);
+  FEFaceValues<dim> ncell_face_values(mapping,
+                                      fe,
+                                      face_quad,
+                                      update_quadrature_points);
+
+  for (const auto &cell : tria.active_cell_iterators())
+    for (unsigned int f : cell->face_indices())
+      {
+        auto ncell = cell->neighbor(f);
+        if (ncell == tria.end())
+          continue;
+
+        auto face = cell->face(f);
+        cell_face_values.reinit(cell, face);
+        unsigned int nf = 0;
+        for (; nf < ref_cell.n_faces(); ++nf)
+          if (ncell->face(nf) == face)
+            break;
+        ncell_face_values.reinit(ncell, nf);
+
+        double max_distance = 0.0;
+        for (unsigned int q = 0; q < face_quad.size(); ++q)
+          max_distance =
+            std::max(max_distance,
+                     cell_face_values.get_quadrature_points()[q].distance(
+                       ncell_face_values.get_quadrature_points()[q]));
+
+        std::cout << orientation << " "
+                  << int(ncell->combined_face_orientation(nf)) << std::endl;
+
+        if (max_distance > 1e-12)
+          {
+            std::cout << "cell points =" << std::endl;
+            for (unsigned int q = 0; q < face_quad.size(); ++q)
+              {
+                std::cout << " " << cell_face_values.get_quadrature_points()[q]
+                          << std::endl;
+              }
+            std::cout << std::endl;
+            std::cout << "face orientation = "
+                      << int(cell->combined_face_orientation(f)) << std::endl;
+
+            std::cout << "ncell points =" << std::endl;
+            for (unsigned int q = 0; q < face_quad.size(); ++q)
+              {
+                std::cout << " " << ncell_face_values.get_quadrature_points()[q]
+                          << std::endl;
+              }
+            std::cout << std::endl;
+            std::cout << "nf orientation = "
+                      << int(ncell->combined_face_orientation(nf)) << std::endl;
+
+            AssertDimension(orientation,
+                            int(ncell->combined_face_orientation(nf)));
+            AssertThrow(false, ExcMessage("Should match"));
+          }
+      }
+}
+
+int
+main()
+{
+  initlog();
+  // test<2>();
+  test<3>();
+
+  deallog << "OK!" << std::endl;
+}
diff --git a/tests/simplex/orientation_03.output b/tests/simplex/orientation_03.output
new file mode 100644 (file)
index 0000000..5cfb783
--- /dev/null
@@ -0,0 +1,2 @@
+
+DEAL::OK!
diff --git a/tests/simplex/orientation_04.cc b/tests/simplex/orientation_04.cc
new file mode 100644 (file)
index 0000000..a7af82a
--- /dev/null
@@ -0,0 +1,78 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2021 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+
+// Test ReferenceCell::permute_according_orientation and compute_orientation
+// by reverting the permutation due to a given orientation.
+
+#include <deal.II/grid/reference_cell.h>
+
+#include "../tests.h"
+
+template <unsigned int n_points>
+void
+test(const ReferenceCell type, const unsigned int n_orientations)
+{
+  for (unsigned int o = 0; o < n_orientations; ++o)
+    {
+      std::array<unsigned int, n_points> origin;
+
+      for (unsigned int i = 0; i < n_points; ++i)
+        origin[i] = i;
+
+      const auto permuted = type.permute_according_orientation(origin, o);
+      const unsigned int origin_back =
+        type.compute_orientation(permuted, origin);
+
+      // AssertDimension(o, origin_back);
+      deallog << o << " " << origin_back << std::endl;
+    }
+
+  deallog << std::endl;
+
+  for (unsigned int o = 0; o < n_orientations; ++o)
+    {
+      std::array<unsigned int, n_points> origin;
+
+      for (unsigned int i = 0; i < n_points; ++i)
+        origin[i] = i;
+
+      const auto permutedi = type.permute_according_orientation(origin, o);
+
+      std::array<unsigned int, n_points> permuted;
+
+      for (unsigned int i = 0; i < n_points; ++i)
+        permuted[permutedi[i]] = origin[i];
+
+      const unsigned int origin_back =
+        type.compute_orientation(origin, permuted);
+
+      // AssertDimension(o, origin_back);
+      deallog << o << " " << origin_back << std::endl;
+    }
+}
+
+int
+main()
+{
+  initlog();
+
+  // test<2>(ReferenceCells::Line, 2);
+  test<3>(ReferenceCells::Triangle, 6);
+  // test<4>(ReferenceCells::Quadrilateral, 4);
+
+  deallog << "OK!" << std::endl;
+}
diff --git a/tests/simplex/orientation_04.output b/tests/simplex/orientation_04.output
new file mode 100644 (file)
index 0000000..2854a84
--- /dev/null
@@ -0,0 +1,15 @@
+
+DEAL::0 0
+DEAL::1 1
+DEAL::2 2
+DEAL::3 3
+DEAL::4 4
+DEAL::5 5
+DEAL::
+DEAL::0 0
+DEAL::1 1
+DEAL::2 2
+DEAL::3 3
+DEAL::4 4
+DEAL::5 5
+DEAL::OK!
index d6e0e3a2b0161f0781222f4e3986cd6c021f7a8d..18003c93c9c0762e7cb5ea6fa59bf2c709f9981e 100644 (file)
@@ -150,17 +150,17 @@ DEAL:3d-4::0.666390 0.155051 0.00000 0.159021
 DEAL:3d-4::0.280020 0.644949 0.00000 0.0909793 
 DEAL:3d-4::
 DEAL:3d-4::face_no=0 face_orientation=true face_flip=true face_rotation=false:
-DEAL:3d-4::0.155051 0.666390 0.00000 0.159021 
-DEAL:3d-4::0.644949 0.280020 0.00000 0.0909793 
-DEAL:3d-4::0.155051 0.178559 0.00000 0.159021 
-DEAL:3d-4::0.644949 0.0750311 0.00000 0.0909793 
-DEAL:3d-4::
-DEAL:3d-4::face_no=0 face_orientation=true face_flip=false face_rotation=true:
 DEAL:3d-4::0.666390 0.178559 0.00000 0.159021 
 DEAL:3d-4::0.280020 0.0750311 0.00000 0.0909793 
 DEAL:3d-4::0.178559 0.666390 0.00000 0.159021 
 DEAL:3d-4::0.0750311 0.280020 0.00000 0.0909793 
 DEAL:3d-4::
+DEAL:3d-4::face_no=0 face_orientation=true face_flip=false face_rotation=true:
+DEAL:3d-4::0.155051 0.666390 0.00000 0.159021 
+DEAL:3d-4::0.644949 0.280020 0.00000 0.0909793 
+DEAL:3d-4::0.155051 0.178559 0.00000 0.159021 
+DEAL:3d-4::0.644949 0.0750311 0.00000 0.0909793 
+DEAL:3d-4::
 DEAL:3d-4::face_no=0 face_orientation=false face_flip=false face_rotation=false:
 DEAL:3d-4::0.155051 0.178559 0.00000 0.159021 
 DEAL:3d-4::0.644949 0.0750311 0.00000 0.0909793 
@@ -186,17 +186,17 @@ DEAL:3d-4::0.178559 0.00000 0.155051 0.159021
 DEAL:3d-4::0.0750311 0.00000 0.644949 0.0909793 
 DEAL:3d-4::
 DEAL:3d-4::face_no=1 face_orientation=true face_flip=true face_rotation=false:
-DEAL:3d-4::0.178559 0.00000 0.666390 0.159021 
-DEAL:3d-4::0.0750311 0.00000 0.280020 0.0909793 
-DEAL:3d-4::0.666390 0.00000 0.178559 0.159021 
-DEAL:3d-4::0.280020 0.00000 0.0750311 0.0909793 
-DEAL:3d-4::
-DEAL:3d-4::face_no=1 face_orientation=true face_flip=false face_rotation=true:
 DEAL:3d-4::0.155051 0.00000 0.178559 0.159021 
 DEAL:3d-4::0.644949 0.00000 0.0750311 0.0909793 
 DEAL:3d-4::0.155051 0.00000 0.666390 0.159021 
 DEAL:3d-4::0.644949 0.00000 0.280020 0.0909793 
 DEAL:3d-4::
+DEAL:3d-4::face_no=1 face_orientation=true face_flip=false face_rotation=true:
+DEAL:3d-4::0.178559 0.00000 0.666390 0.159021 
+DEAL:3d-4::0.0750311 0.00000 0.280020 0.0909793 
+DEAL:3d-4::0.666390 0.00000 0.178559 0.159021 
+DEAL:3d-4::0.280020 0.00000 0.0750311 0.0909793 
+DEAL:3d-4::
 DEAL:3d-4::face_no=1 face_orientation=false face_flip=false face_rotation=false:
 DEAL:3d-4::0.666390 0.00000 0.178559 0.159021 
 DEAL:3d-4::0.280020 0.00000 0.0750311 0.0909793 
@@ -222,17 +222,17 @@ DEAL:3d-4::0.00000 0.666390 0.155051 0.159021
 DEAL:3d-4::0.00000 0.280020 0.644949 0.0909793 
 DEAL:3d-4::
 DEAL:3d-4::face_no=2 face_orientation=true face_flip=true face_rotation=false:
-DEAL:3d-4::0.00000 0.155051 0.666390 0.159021 
-DEAL:3d-4::0.00000 0.644949 0.280020 0.0909793 
-DEAL:3d-4::0.00000 0.155051 0.178559 0.159021 
-DEAL:3d-4::0.00000 0.644949 0.0750311 0.0909793 
-DEAL:3d-4::
-DEAL:3d-4::face_no=2 face_orientation=true face_flip=false face_rotation=true:
 DEAL:3d-4::0.00000 0.666390 0.178559 0.159021 
 DEAL:3d-4::0.00000 0.280020 0.0750311 0.0909793 
 DEAL:3d-4::0.00000 0.178559 0.666390 0.159021 
 DEAL:3d-4::0.00000 0.0750311 0.280020 0.0909793 
 DEAL:3d-4::
+DEAL:3d-4::face_no=2 face_orientation=true face_flip=false face_rotation=true:
+DEAL:3d-4::0.00000 0.155051 0.666390 0.159021 
+DEAL:3d-4::0.00000 0.644949 0.280020 0.0909793 
+DEAL:3d-4::0.00000 0.155051 0.178559 0.159021 
+DEAL:3d-4::0.00000 0.644949 0.0750311 0.0909793 
+DEAL:3d-4::
 DEAL:3d-4::face_no=2 face_orientation=false face_flip=false face_rotation=false:
 DEAL:3d-4::0.00000 0.155051 0.178559 0.159021 
 DEAL:3d-4::0.00000 0.644949 0.0750311 0.0909793 
@@ -258,17 +258,17 @@ DEAL:3d-4::0.666390 0.178559 0.155051 0.275432
 DEAL:3d-4::0.280020 0.0750311 0.644949 0.157581 
 DEAL:3d-4::
 DEAL:3d-4::face_no=3 face_orientation=true face_flip=true face_rotation=false:
-DEAL:3d-4::0.155051 0.178559 0.666390 0.275432 
-DEAL:3d-4::0.644949 0.0750311 0.280020 0.157581 
-DEAL:3d-4::0.155051 0.666390 0.178559 0.275432 
-DEAL:3d-4::0.644949 0.280020 0.0750311 0.157581 
-DEAL:3d-4::
-DEAL:3d-4::face_no=3 face_orientation=true face_flip=false face_rotation=true:
 DEAL:3d-4::0.666390 0.155051 0.178559 0.275432 
 DEAL:3d-4::0.280020 0.644949 0.0750311 0.157581 
 DEAL:3d-4::0.178559 0.155051 0.666390 0.275432 
 DEAL:3d-4::0.0750311 0.644949 0.280020 0.157581 
 DEAL:3d-4::
+DEAL:3d-4::face_no=3 face_orientation=true face_flip=false face_rotation=true:
+DEAL:3d-4::0.155051 0.178559 0.666390 0.275432 
+DEAL:3d-4::0.644949 0.0750311 0.280020 0.157581 
+DEAL:3d-4::0.155051 0.666390 0.178559 0.275432 
+DEAL:3d-4::0.644949 0.280020 0.0750311 0.157581 
+DEAL:3d-4::
 DEAL:3d-4::face_no=3 face_orientation=false face_flip=false face_rotation=false:
 DEAL:3d-4::0.155051 0.666390 0.178559 0.275432 
 DEAL:3d-4::0.644949 0.280020 0.0750311 0.157581 
@@ -298,21 +298,21 @@ DEAL:3d-10::0.470142 0.470142 0.00000 0.0661971
 DEAL:3d-10::
 DEAL:3d-10::face_no=0 face_orientation=true face_flip=true face_rotation=false:
 DEAL:3d-10::0.333333 0.333333 0.00000 0.112500 
+DEAL:3d-10::0.101287 0.797427 0.00000 0.0629696 
 DEAL:3d-10::0.101287 0.101287 0.00000 0.0629696 
 DEAL:3d-10::0.797427 0.101287 0.00000 0.0629696 
-DEAL:3d-10::0.101287 0.797427 0.00000 0.0629696 
+DEAL:3d-10::0.470142 0.0597159 0.00000 0.0661971 
 DEAL:3d-10::0.470142 0.470142 0.00000 0.0661971 
 DEAL:3d-10::0.0597159 0.470142 0.00000 0.0661971 
-DEAL:3d-10::0.470142 0.0597159 0.00000 0.0661971 
 DEAL:3d-10::
 DEAL:3d-10::face_no=0 face_orientation=true face_flip=false face_rotation=true:
 DEAL:3d-10::0.333333 0.333333 0.00000 0.112500 
-DEAL:3d-10::0.101287 0.797427 0.00000 0.0629696 
 DEAL:3d-10::0.101287 0.101287 0.00000 0.0629696 
 DEAL:3d-10::0.797427 0.101287 0.00000 0.0629696 
-DEAL:3d-10::0.470142 0.0597159 0.00000 0.0661971 
+DEAL:3d-10::0.101287 0.797427 0.00000 0.0629696 
 DEAL:3d-10::0.470142 0.470142 0.00000 0.0661971 
 DEAL:3d-10::0.0597159 0.470142 0.00000 0.0661971 
+DEAL:3d-10::0.470142 0.0597159 0.00000 0.0661971 
 DEAL:3d-10::
 DEAL:3d-10::face_no=0 face_orientation=false face_flip=false face_rotation=false:
 DEAL:3d-10::0.333333 0.333333 0.00000 0.112500 
@@ -352,21 +352,21 @@ DEAL:3d-10::0.0597159 0.00000 0.470142 0.0661971
 DEAL:3d-10::
 DEAL:3d-10::face_no=1 face_orientation=true face_flip=true face_rotation=false:
 DEAL:3d-10::0.333333 0.00000 0.333333 0.112500 
+DEAL:3d-10::0.101287 0.00000 0.797427 0.0629696 
 DEAL:3d-10::0.797427 0.00000 0.101287 0.0629696 
 DEAL:3d-10::0.101287 0.00000 0.101287 0.0629696 
-DEAL:3d-10::0.101287 0.00000 0.797427 0.0629696 
+DEAL:3d-10::0.470142 0.00000 0.0597159 0.0661971 
 DEAL:3d-10::0.0597159 0.00000 0.470142 0.0661971 
 DEAL:3d-10::0.470142 0.00000 0.470142 0.0661971 
-DEAL:3d-10::0.470142 0.00000 0.0597159 0.0661971 
 DEAL:3d-10::
 DEAL:3d-10::face_no=1 face_orientation=true face_flip=false face_rotation=true:
 DEAL:3d-10::0.333333 0.00000 0.333333 0.112500 
-DEAL:3d-10::0.101287 0.00000 0.797427 0.0629696 
 DEAL:3d-10::0.797427 0.00000 0.101287 0.0629696 
 DEAL:3d-10::0.101287 0.00000 0.101287 0.0629696 
-DEAL:3d-10::0.470142 0.00000 0.0597159 0.0661971 
+DEAL:3d-10::0.101287 0.00000 0.797427 0.0629696 
 DEAL:3d-10::0.0597159 0.00000 0.470142 0.0661971 
 DEAL:3d-10::0.470142 0.00000 0.470142 0.0661971 
+DEAL:3d-10::0.470142 0.00000 0.0597159 0.0661971 
 DEAL:3d-10::
 DEAL:3d-10::face_no=1 face_orientation=false face_flip=false face_rotation=false:
 DEAL:3d-10::0.333333 0.00000 0.333333 0.112500 
@@ -406,21 +406,21 @@ DEAL:3d-10::0.00000 0.470142 0.470142 0.0661971
 DEAL:3d-10::
 DEAL:3d-10::face_no=2 face_orientation=true face_flip=true face_rotation=false:
 DEAL:3d-10::0.00000 0.333333 0.333333 0.112500 
+DEAL:3d-10::0.00000 0.101287 0.797427 0.0629696 
 DEAL:3d-10::0.00000 0.101287 0.101287 0.0629696 
 DEAL:3d-10::0.00000 0.797427 0.101287 0.0629696 
-DEAL:3d-10::0.00000 0.101287 0.797427 0.0629696 
+DEAL:3d-10::0.00000 0.470142 0.0597159 0.0661971 
 DEAL:3d-10::0.00000 0.470142 0.470142 0.0661971 
 DEAL:3d-10::0.00000 0.0597159 0.470142 0.0661971 
-DEAL:3d-10::0.00000 0.470142 0.0597159 0.0661971 
 DEAL:3d-10::
 DEAL:3d-10::face_no=2 face_orientation=true face_flip=false face_rotation=true:
 DEAL:3d-10::0.00000 0.333333 0.333333 0.112500 
-DEAL:3d-10::0.00000 0.101287 0.797427 0.0629696 
 DEAL:3d-10::0.00000 0.101287 0.101287 0.0629696 
 DEAL:3d-10::0.00000 0.797427 0.101287 0.0629696 
-DEAL:3d-10::0.00000 0.470142 0.0597159 0.0661971 
+DEAL:3d-10::0.00000 0.101287 0.797427 0.0629696 
 DEAL:3d-10::0.00000 0.470142 0.470142 0.0661971 
 DEAL:3d-10::0.00000 0.0597159 0.470142 0.0661971 
+DEAL:3d-10::0.00000 0.470142 0.0597159 0.0661971 
 DEAL:3d-10::
 DEAL:3d-10::face_no=2 face_orientation=false face_flip=false face_rotation=false:
 DEAL:3d-10::0.00000 0.333333 0.333333 0.112500 
@@ -460,21 +460,21 @@ DEAL:3d-10::0.470142 0.0597159 0.470142 0.114657
 DEAL:3d-10::
 DEAL:3d-10::face_no=3 face_orientation=true face_flip=true face_rotation=false:
 DEAL:3d-10::0.333333 0.333333 0.333333 0.194856 
+DEAL:3d-10::0.101287 0.101287 0.797427 0.109067 
 DEAL:3d-10::0.101287 0.797427 0.101287 0.109067 
 DEAL:3d-10::0.797427 0.101287 0.101287 0.109067 
-DEAL:3d-10::0.101287 0.101287 0.797427 0.10906
+DEAL:3d-10::0.470142 0.470142 0.0597159 0.11465
 DEAL:3d-10::0.470142 0.0597159 0.470142 0.114657 
 DEAL:3d-10::0.0597159 0.470142 0.470142 0.114657 
-DEAL:3d-10::0.470142 0.470142 0.0597159 0.114657 
 DEAL:3d-10::
 DEAL:3d-10::face_no=3 face_orientation=true face_flip=false face_rotation=true:
 DEAL:3d-10::0.333333 0.333333 0.333333 0.194856 
-DEAL:3d-10::0.101287 0.101287 0.797427 0.109067 
 DEAL:3d-10::0.101287 0.797427 0.101287 0.109067 
 DEAL:3d-10::0.797427 0.101287 0.101287 0.109067 
-DEAL:3d-10::0.470142 0.470142 0.0597159 0.11465
+DEAL:3d-10::0.101287 0.101287 0.797427 0.10906
 DEAL:3d-10::0.470142 0.0597159 0.470142 0.114657 
 DEAL:3d-10::0.0597159 0.470142 0.470142 0.114657 
+DEAL:3d-10::0.470142 0.470142 0.0597159 0.114657 
 DEAL:3d-10::
 DEAL:3d-10::face_no=3 face_orientation=false face_flip=false face_rotation=false:
 DEAL:3d-10::0.333333 0.333333 0.333333 0.194856 

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.