]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Throw an exception if we create cells that are distorted during mesh refinement....
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Mon, 13 Jul 2009 03:07:03 +0000 (03:07 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Mon, 13 Jul 2009 03:07:03 +0000 (03:07 +0000)
git-svn-id: https://svn.dealii.org/trunk@19061 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/grid/tria.h
deal.II/deal.II/source/grid/persistent_tria.cc
deal.II/deal.II/source/grid/tria.cc
deal.II/doc/doxygen/headers/glossary.h
deal.II/doc/news/changes.h
tests/bits/distorted_cells_03.cc [new file with mode: 0644]
tests/bits/distorted_cells_03/2d [new file with mode: 0644]
tests/bits/distorted_cells_03/3d [new file with mode: 0644]
tests/bits/distorted_cells_03/cmp/generic [new file with mode: 0644]

index 0c9e26fbcf7943016f36a13950c7881677173dcc..c2ce29771ab31f6d90845a977e7eb08dec2b4827 100644 (file)
@@ -1780,6 +1780,8 @@ class Triangulation : public Subscriptor
                                      *  overload the
                                      *  @p execute_coarsening_and_refinement
                                      *  function.
+                                     *
+                                     * 
                                      */
     void refine_global (const unsigned int times);
 
@@ -1794,6 +1796,28 @@ class Triangulation : public Subscriptor
                                      * will therefore be overwritten by
                                      * undefined content.
                                      *
+                                     * If the boundary description is
+                                     * sufficiently irregular, it can
+                                     * happen that some of the
+                                     * children produced by mesh
+                                     * refinement are distorted (see
+                                     * the extensive discussion on
+                                     * @ref GlossDistorted "distorted cells").
+                                     *
+                                     * To allow user programs to fix
+                                     * up these cells if that is
+                                     * desired, this function after
+                                     * completing all other work may
+                                     * throw an exception of type
+                                     * DistortedCellList that
+                                     * contains a list of those cells
+                                     * that have been refined and
+                                     * have at least one child that
+                                     * is distorted. The function
+                                     * does not create such an
+                                     * exception if no cells have
+                                     * created distorted children.
+                                     *
                                       * See the general docs for more
                                       * information.
                                      *
@@ -3161,8 +3185,14 @@ class Triangulation : public Subscriptor
                                      *  for <tt>dim=2,3</tt> and the
                                      *  <tt>quad->user_flags</tt> for
                                      *  <tt>dim=3</tt>.
+                                     *
+                                     *  The function returns a list
+                                     *  of cells that have produced
+                                     *  children that satisfy the
+                                     *  criteria of
+                                     *  @ref GlossDistorted "distorted cells".
                                      */ 
-    void execute_refinement ();
+    DistortedCellList execute_refinement ();
 
                                     /**
                                      * Coarsen all cells which were flagged for
index f4747d391f81ca52d4bcd908c7c3f7c309996500..36ec9cee1cc67c95891caa8593a4cb70c1d977e2 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 1999, 2000, 2001, 2002, 2003, 2005, 2006, 2007 by the deal.II authors
+//    Copyright (C) 1999, 2000, 2001, 2002, 2003, 2005, 2006, 2007, 2009 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -67,7 +67,10 @@ PersistentTriangulation<dim>::execute_coarsening_and_refinement ()
   this->save_refine_flags (refine_flags.back());
   this->save_coarsen_flags (coarsen_flags.back());
 
-                                  // then refine triangulation
+                                  // then refine triangulation. if
+                                  // this function throws an
+                                  // exception, that's fine since it
+                                  // is the last call here
   Triangulation<dim>::execute_coarsening_and_refinement ();
 }
 
@@ -107,7 +110,7 @@ PersistentTriangulation<dim>::restore (const unsigned int step)
       this->load_coarsen_flags (coarsen_flags[step-1]);
 
       Triangulation<dim>::execute_coarsening_and_refinement ();
-    };
+    }
 }
 
 
index 3e58529a2ac3ed807fa943cd4c69290b3a8d57be..04514d8c9e2e81d5b378ae8ea3221d0c83dde896 100644 (file)
@@ -1189,6 +1189,59 @@ namespace
   {
     return typename Triangulation<dim,spacedim>::DistortedCellList();
   }
+
+
+
+                                  /**
+                                   * Return whether any of the
+                                   * children of the given cell is
+                                   * distorted or not. This is the
+                                   * function for dim==spacedim.
+                                   */
+  template <int dim>
+  bool
+  has_distorted_children (const typename Triangulation<dim,dim>::cell_iterator &cell,
+                         internal::int2type<dim>,
+                         internal::int2type<dim>)
+  {
+    Assert (cell->has_children(), ExcInternalError());
+
+    for (unsigned int c=0; c<cell->n_children(); ++c)
+      {
+       Point<dim> vertices[GeometryInfo<dim>::vertices_per_cell];
+       for (unsigned int i=0; i<GeometryInfo<dim>::vertices_per_cell; ++i)
+         vertices[i] = cell->child(c)->vertex(i);
+
+       double determinants[GeometryInfo<dim>::vertices_per_cell];
+       GeometryInfo<dim>::jacobian_determinants_at_vertices (vertices,
+                                                             determinants);
+
+       for (unsigned int i=0; i<GeometryInfo<dim>::vertices_per_cell; ++i)
+         if (determinants[i] <= 1e-9 * std::pow (cell->child(c)->diameter(),
+                                                 1.*dim))
+           return true;
+      }
+
+    return false;
+  }
+
+
+                                  /**
+                                   * Function for dim!=spacedim. As
+                                   * for
+                                   * collect_distorted_coarse_cells,
+                                   * there is nothing that we can do
+                                   * in this case.
+                                   */
+  template <int dim, int spacedim>
+  bool
+  has_distorted_children (const typename Triangulation<dim,spacedim>::cell_iterator &,
+                         internal::int2type<dim>,
+                         internal::int2type<spacedim>)
+  {
+    return false;
+  }
+  
   
 }// end of anonymous namespace
 
@@ -4564,7 +4617,7 @@ namespace internal
                                          */
        template <int spacedim>
        static
-       void
+       typename Triangulation<1,spacedim>::DistortedCellList
        execute_refinement (Triangulation<1,spacedim> &triangulation)
          {
            const unsigned int dim = 1;
@@ -4814,6 +4867,17 @@ namespace internal
                          }
                    }      
              }
+
+                                            // in 1d, we can not have
+                                            // distorted children
+                                            // unless the parent was
+                                            // already distorted
+                                            // (that is because we
+                                            // don't use boundary
+                                            // information for 1d
+                                            // triangulations). so
+                                            // return an empty list
+           return typename Triangulation<1,spacedim>::DistortedCellList();
          }
 
 
@@ -4823,7 +4887,7 @@ namespace internal
                                          */
        template <int spacedim>
        static
-       void
+       typename Triangulation<2,spacedim>::DistortedCellList
        execute_refinement (Triangulation<2,spacedim> &triangulation)
          {
            const unsigned int dim = 2;
@@ -5073,7 +5137,8 @@ namespace internal
                                                           // the two vertices:
                          if (line->at_boundary())
                            triangulation.vertices[next_unused_vertex]
-                             = triangulation.boundary[line->boundary_indicator()]->get_new_point_on_line (line);
+                             = triangulation.boundary[line->boundary_indicator()]
+                             ->get_new_point_on_line (line);
                          else
                            triangulation.vertices[next_unused_vertex]
                              = (line->vertex(0) + line->vertex(1)) / 2;
@@ -5163,6 +5228,9 @@ namespace internal
            triangulation.faces->lines.
              reserve_space (0,n_single_lines);
 
+           typename Triangulation<2,spacedim>::DistortedCellList
+             cells_with_distorted_children;
+           
                                             // reset next_unused_line, as
                                             // now also single empty places
                                             // in the vector can be used
@@ -5203,8 +5271,15 @@ namespace internal
                                       next_unused_line,
                                       next_unused_cell,
                                       cell);
+
+                     if (has_distorted_children (cell,
+                                                 internal::int2type<dim>(),
+                                                 internal::int2type<spacedim>()))
+                       cells_with_distorted_children.distorted_cells.push_back (cell);
                    }
              }
+
+           return cells_with_distorted_children;
          }
 
 
@@ -5214,7 +5289,7 @@ namespace internal
                                          */
        template <int spacedim>
        static
-       void
+       typename Triangulation<3,spacedim>::DistortedCellList
        execute_refinement (Triangulation<3,spacedim> &triangulation)
          {
            const unsigned int dim = 3;
@@ -6520,6 +6595,10 @@ namespace internal
                                             // Now, finally, set up the new
                                             // cells
                                             ///////////////////////////////////
+
+           typename Triangulation<3,spacedim>::DistortedCellList
+             cells_with_distorted_children;
+
            for (unsigned int level=0; level!=triangulation.levels.size()-1; ++level)
              {
                                                 // only active objects can be
@@ -9256,6 +9335,20 @@ namespace internal
                                                       // separate function
                      update_neighbors<spacedim> (hex, true);
 
+                                                      // now see if
+                                                      // we have
+                                                      // created
+                                                      // cells that
+                                                      // are
+                                                      // distorted
+                                                      // and if so
+                                                      // add them to
+                                                      // our list
+                     if (has_distorted_children (hex,
+                                                 internal::int2type<dim>(),
+                                                 internal::int2type<spacedim>()))
+                       cells_with_distorted_children.distorted_cells.push_back (hex);
+                     
                                                       // note that the
                                                       // refinement flag was
                                                       // already cleared at the
@@ -9270,7 +9363,10 @@ namespace internal
                                             // whether we used indices or pointers is
                                             // still present. reset it now to enable the
                                             // user to use whichever he likes later on.
-           triangulation.faces->quads.clear_user_data();       
+           triangulation.faces->quads.clear_user_data();
+
+                                            // return the list with distorted children
+           return cells_with_distorted_children;
          }
 
 
@@ -12624,7 +12720,10 @@ Triangulation<dim, spacedim>::execute_coarsening_and_refinement ()
     (*ref_listener)->pre_refinement_notification (*this);
 
   execute_coarsening();
-  execute_refinement();
+
+  const DistortedCellList
+    cells_with_distorted_children
+    = execute_refinement();
 
                                   // verify a case with which we have had
                                   // some difficulty in the past (see the
@@ -12638,6 +12737,9 @@ Triangulation<dim, spacedim>::execute_coarsening_and_refinement ()
   for (ref_listener = refinement_listeners.begin ();
        ref_listener != end_listener; ++ref_listener)
     (*ref_listener)->post_refinement_notification (*this);
+
+  AssertThrow (cells_with_distorted_children.distorted_cells.size() == 0,
+              cells_with_distorted_children);
 }
 
 
@@ -12666,10 +12768,12 @@ Triangulation<dim, spacedim>::clear_despite_subscriptions()
 
 
 template <int dim, int spacedim>
-void
+typename Triangulation<dim,spacedim>::DistortedCellList
 Triangulation<dim,spacedim>::execute_refinement ()
 {
-  internal::Triangulation::Implementation::execute_refinement (*this);
+  const DistortedCellList
+    cells_with_distorted_children
+    = internal::Triangulation::Implementation::execute_refinement (*this);
 
 
 
@@ -12693,7 +12797,9 @@ Triangulation<dim,spacedim>::execute_refinement ()
                endc = end();
   while (cell != endc)
     Assert (!(cell++)->refine_flag_set(), ExcInternalError ());  
-#endif  
+#endif
+
+  return cells_with_distorted_children;
 }
 
 
index 42927b0b921c362a849d4e4306f58eae44f0d208..ce3fd1b44505934fef9a6ac486aabba007e0c914 100644 (file)
  * non-positive somewhere in the cell. Typically, we only check the sign
  * of this determinant at the vertices of the cell. The function
  * GeometryInfo::jacobian_determinants_at_vertices computes these
- * determinants at the vertices, and the function
- * Triangulation::create_triangulation (along with the various grid
- * readers in the GridIn class) will reject meshes with distorted cells
- * by default.
+ * determinants at the vertices.
  *
  * By way of example, if all of the determinants are of roughly equal value
  * and on the order of $h^\text{dim}$ then the cell is well-shaped. For
  * @image html distorted_3d.png "A well-formed, a pinched, and a twisted cell in 3d."
  * </dd>
  *
+ * Distorted cells can appear in two different ways: The original
+ * coarse mesh can already contain such cells, or they can be created
+ * as the result of mesh refinement if the boundary description in use
+ * is sufficiently irregular.
+ *
+ * The function Triangulation::create_triangulation, which is called
+ * by the various functions in GridGenerator and GridIn (but can also
+ * be called from user code, see @ref step_14 "step-14" will signal
+ * the creation of coarse meshes with distorted cells by throwing an
+ * exception of type Triangulation::DistortedCellList. There are
+ * legitimate cases for creating meshes with distorted cells (in
+ * particular collapsed/pinched cells) if you don't intend to assemble
+ * anything on these cells. For example, consider a case where one
+ * would like to simulate the behavior of an elastic material with a
+ * fluid-filled crack such as an oil reservoir. If the pressure
+ * becomes too large, the crack is closed -- and the cells that
+ * discretize the crack volume are collapsed to zero volume. As long
+ * as you don't integrate over these cells to simulate the behavior of
+ * the fluid (of which there isn't any if the crack has zero volume),
+ * such meshes are perfectly legitimate. As a consequence,
+ * Triangulation::create_triangulation does not simply abort the
+ * program, but throws an exception that contains a list of cells that
+ * are distorted; this exception can be caught and, if you believe
+ * that you can ignore this condition, you can react by doing nothing
+ * with the caught exception.
+ *
+ * The second case in which distorted cells can appear is through mesh
+ * refinement when we have curved boundaries. Consider, for example,
+ * the following case where he dashed line shows the exact boundary
+ * that the lower edge of the cell is supposed to approximate (let's
+ * assume that the left, top and right edges are interior edges and
+ * therefore will be considered as straight):
+ *
+ * @image html distorted_2d_refinement_01.png "One cell with a edge approximating a curved boundary"
+ *
+ * Now, if this cell is refined, we first split all edges and place
+ * new mid-points on them. For the left, top and right edge, this is
+ * trivial: because they are considered straight, we just take the
+ * point in the middle between the two vertices. For the lower edge,
+ * the Triangulation class asks the Boundary object associated with
+ * this boundary (and in particular the Boundary::new_point_on_line
+ * function) where the new point should lie. The four old vertices and
+ * the four new points are shown here:
+ *
+ * @image html distorted_2d_refinement_02.png "Cell after edge refinement"
+ *
+ * The last step is to compute the location of the new point in the
+ * interior of the cell. By default, it is chosen as the average
+ * location (arithmetic mean of the coordinates) of these 8 points:
+ *
+ * @image html distorted_2d_refinement_03.png "Cell after edge refinement"
+ *
+ * The problem with that is, of course, that the bottom two child cells are
+ * twisted, whereas the top two children are well-shaped. While such
+ * meshes can happen with sufficiently irregular boundary descriptions
+ * (and if the coarse mesh is entirely inadequate to resolve the
+ * complexity of the boundary), the Triangulation class does not know
+ * what to do in such situations. Consequently, the
+ * Triangulation::execute_coarsening_and_refinement function does
+ * create such meshes, but it keeps a list of cells whose children are
+ * distorted. If this list is non-empty at the end of a refinement
+ * step, it will throw an exception of type
+ * Triangualtion::DistortedCellList that contains those cells that
+ * have distorted children. The caller of
+ * Triangulation::execute_coarsening_and_refinement can then decide
+ * what to do with this situation.
+ *
  *
  * <dt class="glossary">@anchor GlossFaceOrientation Face orientation</dt>
  * <dd>In a triangulation, the normal vector to a face
index 11c0a1666a795af80ee0c95a6c1c107edbf55bf1..142dee9b73433f70c8c34e2d5a23499a6831a0d9 100644 (file)
@@ -32,16 +32,24 @@ inconvenience this causes.
 <ol>
   <li>
   <p>
-  Changed: Previously, the Triangulation::create_triangulation function silently
-  accepted input meshes with inverted cells (i.e. cells with a zero or negative
-  determinant of the Jacobian of the mapping from the reference cell to the
-  real cell). This has been changed now: The function checks whether cells
-  are distorted or inverted (see the entry on @ref GlossDistorted "distorted cells"
-  in the glossary), and may throw an exception containing a list of cells for which this
-  is the case. If you know that this is harmless, for example
-  if you have cells with collapsed vertices in your mesh but you do not intend
-  to integrate on them, then you can catch and ignore this message. In all
-  other cases, the output of your computations are likely to be wrong anyway.
+  Changed: Previously, the Triangulation::create_triangulation
+  function silently accepted input meshes with inverted cells
+  (i.e. cells with a zero or negative determinant of the Jacobian of
+  the mapping from the reference cell to the real cell). This has been
+  changed now: The function checks whether cells are distorted or
+  inverted, and may throw an exception containing a list of cells
+  for which this is the case. If you know that this is harmless, for
+  example if you have cells with collapsed vertices in your mesh but
+  you do not intend to integrate on them, then you can catch and
+  ignore this message. In all other cases, the output of your
+  computations are likely to be wrong anyway.
+  <br>
+  The same is true for the Triangulation::execute_coarsening_and_refinement
+  function: if it creates cells that are distorted, it throws a list of cells
+  whose children are distorted.
+  <br>
+  The whole issue is described in some detail in the entry on
+  @ref GlossDistorted "distorted cells" in the glossary.
   <br>
   (WB 2009/06/29)
   </p>
diff --git a/tests/bits/distorted_cells_03.cc b/tests/bits/distorted_cells_03.cc
new file mode 100644 (file)
index 0000000..0f79c32
--- /dev/null
@@ -0,0 +1,128 @@
+//----------------------------  distorted_cells_03.cc  ---------------------------
+//    $Id$
+//    Version: $Name$ 
+//
+//    Copyright (C) 2003, 2004, 2005, 2009 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//----------------------------  distorted_cells_03.cc  ---------------------------
+
+
+// check that indeed Triangulation::execute_coarsening_and_refinement
+// throws an exception if a distorted child cell is created
+//
+// the 2d case is as described in the example in the Glossary on
+// distorted cells resulting from refinement
+
+#include "../tests.h"
+#include <base/logstream.h>
+#include <base/quadrature_lib.h>
+#include <grid/tria.h>
+#include <grid/tria_accessor.h>
+#include <grid/tria_iterator.h>
+#include <grid/grid_reordering.h>
+#include <grid/grid_generator.h>
+#include <grid/tria_boundary.h>
+#include <grid/grid_out.h>
+#include <dofs/dof_handler.h>
+#include <fe/fe_q.h>
+#include <fe/fe_values.h>
+
+#include <fstream>
+
+
+template <int dim>
+class MyBoundary : public Boundary<dim>
+{
+    virtual Point<dim>
+    get_new_point_on_line (const typename Triangulation<dim>::line_iterator &line) const
+      {
+       deallog << "Finding point between "
+               << line->vertex(0) << " and "
+               << line->vertex(1) << std::endl;
+
+                                        // in 2d, find a point that
+                                        // lies on the opposite side
+                                        // of the quad. in 3d, choose
+                                        // the midpoint of the edge
+       if (dim == 2)
+         return Point<dim>(0,0.75);
+       else
+         return (line->vertex(0) + line->vertex(1)) / 2;
+      }
+
+    virtual Point<dim>
+    get_new_point_on_quad (const typename Triangulation<dim>::quad_iterator &quad) const
+      {
+       deallog << "Finding point between "
+               << quad->vertex(0) << " and "
+               << quad->vertex(1) << " and "
+               << quad->vertex(2) << " and "
+               << quad->vertex(3) << std::endl;
+
+       return Point<dim>(0,0,.75);
+      }
+};
+
+
+
+template <int dim>
+void check ()
+{
+  MyBoundary<dim> my_boundary;
+  
+                                  // create a single square/cube
+  Triangulation<dim> coarse_grid;
+  GridGenerator::hyper_cube (coarse_grid, -1, 1);
+
+                                  // set bottom face to use MyBoundary
+  for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
+    if (coarse_grid.begin_active()->face(f)->center()[dim-1] == -1)
+      coarse_grid.begin_active()->face(f)->set_boundary_indicator (1);
+  coarse_grid.set_boundary (1, my_boundary);
+
+                                  // now try to refine this one
+                                  // cell. we should get an exception
+  bool flag = false;
+  try
+    {
+      coarse_grid.begin_active()->set_refine_flag ();
+      coarse_grid.execute_coarsening_and_refinement ();
+    }
+  catch (typename Triangulation<dim>::DistortedCellList &dcv)
+    {
+      flag = true;
+      
+      deallog << dcv.distorted_cells.size() << " distorted cells"
+             << std::endl;
+      Assert (dcv.distorted_cells.front() == coarse_grid.begin(0),
+             ExcInternalError());
+    }
+
+  Assert (flag == true, ExcInternalError());
+  Assert (coarse_grid.n_levels() == 2, ExcInternalError());
+  Assert (coarse_grid.n_active_cells() == 1<<dim, ExcInternalError());
+
+                                  // output the coordinates of the
+                                  // child cells
+  GridOut().write_gnuplot (coarse_grid, deallog.get_file_stream());
+}
+
+
+int main () 
+{
+  std::ofstream logfile("distorted_cells_03/output");
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-10);
+
+  check<2> ();
+  check<3> ();
+}
+
+  
+  
diff --git a/tests/bits/distorted_cells_03/2d b/tests/bits/distorted_cells_03/2d
new file mode 100644 (file)
index 0000000..1a53889
--- /dev/null
@@ -0,0 +1,6 @@
+4 1 0 0 0 
+1 0 0 0
+2 1 0 0 
+3 0 1 0
+4 1 1 0
+1 0 quad 1 2 3 4
diff --git a/tests/bits/distorted_cells_03/3d b/tests/bits/distorted_cells_03/3d
new file mode 100644 (file)
index 0000000..dd37099
--- /dev/null
@@ -0,0 +1,10 @@
+8 1 0 0 0
+1 0 0 0
+2 1 0 0 
+3 0 1 0
+4 1 1 0
+5 0 0 1
+6 1 0 1 
+7 0 1 1
+8 1 1 1
+1 0 hex 1 2 3 4 5 6 7 8
diff --git a/tests/bits/distorted_cells_03/cmp/generic b/tests/bits/distorted_cells_03/cmp/generic
new file mode 100644 (file)
index 0000000..3564f85
--- /dev/null
@@ -0,0 +1,225 @@
+
+DEAL::Finding point between -1.00000 -1.00000 and 1.00000 -1.00000
+DEAL::1 distorted cells
+-1.00000 -1.00000 1 0
+0.00000 0.750000 1 0
+0.00000 0.218750 1 0
+-1.00000 0.00000 1 0
+-1.00000 -1.00000 1 0
+
+
+0.00000 0.750000 1 0
+1.00000 -1.00000 1 0
+1.00000 0.00000 1 0
+0.00000 0.218750 1 0
+0.00000 0.750000 1 0
+
+
+-1.00000 0.00000 1 0
+0.00000 0.218750 1 0
+0.00000 1.00000 1 0
+-1.00000 1.00000 1 0
+-1.00000 0.00000 1 0
+
+
+0.00000 0.218750 1 0
+1.00000 0.00000 1 0
+1.00000 1.00000 1 0
+0.00000 1.00000 1 0
+0.00000 0.218750 1 0
+
+
+DEAL::Finding point between -1.00000 -1.00000 -1.00000 and 1.00000 -1.00000 -1.00000 and -1.00000 1.00000 -1.00000 and 1.00000 1.00000 -1.00000
+DEAL::1 distorted cells
+-1.00000 -1.00000 -1.00000 1 0
+0.00000 -1.00000 -1.00000 1 0
+0.00000 -1.00000 0.00000 1 0
+-1.00000 -1.00000 0.00000 1 0
+-1.00000 -1.00000 -1.00000 1 0
+
+-1.00000 0.00000 -1.00000 1 0
+0.00000 0.00000 0.750000 1 0
+0.00000 0.00000 0.145833 1 0
+-1.00000 0.00000 0.00000 1 0
+-1.00000 0.00000 -1.00000 1 0
+
+-1.00000 -1.00000 -1.00000 1 0
+-1.00000 0.00000 -1.00000 1 0
+
+0.00000 -1.00000 -1.00000 1 0
+0.00000 0.00000 0.750000 1 0
+
+0.00000 -1.00000 0.00000 1 0
+0.00000 0.00000 0.145833 1 0
+
+-1.00000 -1.00000 0.00000 1 0
+-1.00000 0.00000 0.00000 1 0
+
+0.00000 -1.00000 -1.00000 1 0
+1.00000 -1.00000 -1.00000 1 0
+1.00000 -1.00000 0.00000 1 0
+0.00000 -1.00000 0.00000 1 0
+0.00000 -1.00000 -1.00000 1 0
+
+0.00000 0.00000 0.750000 1 0
+1.00000 0.00000 -1.00000 1 0
+1.00000 0.00000 0.00000 1 0
+0.00000 0.00000 0.145833 1 0
+0.00000 0.00000 0.750000 1 0
+
+0.00000 -1.00000 -1.00000 1 0
+0.00000 0.00000 0.750000 1 0
+
+1.00000 -1.00000 -1.00000 1 0
+1.00000 0.00000 -1.00000 1 0
+
+1.00000 -1.00000 0.00000 1 0
+1.00000 0.00000 0.00000 1 0
+
+0.00000 -1.00000 0.00000 1 0
+0.00000 0.00000 0.145833 1 0
+
+-1.00000 0.00000 -1.00000 1 0
+0.00000 0.00000 0.750000 1 0
+0.00000 0.00000 0.145833 1 0
+-1.00000 0.00000 0.00000 1 0
+-1.00000 0.00000 -1.00000 1 0
+
+-1.00000 1.00000 -1.00000 1 0
+0.00000 1.00000 -1.00000 1 0
+0.00000 1.00000 0.00000 1 0
+-1.00000 1.00000 0.00000 1 0
+-1.00000 1.00000 -1.00000 1 0
+
+-1.00000 0.00000 -1.00000 1 0
+-1.00000 1.00000 -1.00000 1 0
+
+0.00000 0.00000 0.750000 1 0
+0.00000 1.00000 -1.00000 1 0
+
+0.00000 0.00000 0.145833 1 0
+0.00000 1.00000 0.00000 1 0
+
+-1.00000 0.00000 0.00000 1 0
+-1.00000 1.00000 0.00000 1 0
+
+0.00000 0.00000 0.750000 1 0
+1.00000 0.00000 -1.00000 1 0
+1.00000 0.00000 0.00000 1 0
+0.00000 0.00000 0.145833 1 0
+0.00000 0.00000 0.750000 1 0
+
+0.00000 1.00000 -1.00000 1 0
+1.00000 1.00000 -1.00000 1 0
+1.00000 1.00000 0.00000 1 0
+0.00000 1.00000 0.00000 1 0
+0.00000 1.00000 -1.00000 1 0
+
+0.00000 0.00000 0.750000 1 0
+0.00000 1.00000 -1.00000 1 0
+
+1.00000 0.00000 -1.00000 1 0
+1.00000 1.00000 -1.00000 1 0
+
+1.00000 0.00000 0.00000 1 0
+1.00000 1.00000 0.00000 1 0
+
+0.00000 0.00000 0.145833 1 0
+0.00000 1.00000 0.00000 1 0
+
+-1.00000 -1.00000 0.00000 1 0
+0.00000 -1.00000 0.00000 1 0
+0.00000 -1.00000 1.00000 1 0
+-1.00000 -1.00000 1.00000 1 0
+-1.00000 -1.00000 0.00000 1 0
+
+-1.00000 0.00000 0.00000 1 0
+0.00000 0.00000 0.145833 1 0
+0.00000 0.00000 1.00000 1 0
+-1.00000 0.00000 1.00000 1 0
+-1.00000 0.00000 0.00000 1 0
+
+-1.00000 -1.00000 0.00000 1 0
+-1.00000 0.00000 0.00000 1 0
+
+0.00000 -1.00000 0.00000 1 0
+0.00000 0.00000 0.145833 1 0
+
+0.00000 -1.00000 1.00000 1 0
+0.00000 0.00000 1.00000 1 0
+
+-1.00000 -1.00000 1.00000 1 0
+-1.00000 0.00000 1.00000 1 0
+
+0.00000 -1.00000 0.00000 1 0
+1.00000 -1.00000 0.00000 1 0
+1.00000 -1.00000 1.00000 1 0
+0.00000 -1.00000 1.00000 1 0
+0.00000 -1.00000 0.00000 1 0
+
+0.00000 0.00000 0.145833 1 0
+1.00000 0.00000 0.00000 1 0
+1.00000 0.00000 1.00000 1 0
+0.00000 0.00000 1.00000 1 0
+0.00000 0.00000 0.145833 1 0
+
+0.00000 -1.00000 0.00000 1 0
+0.00000 0.00000 0.145833 1 0
+
+1.00000 -1.00000 0.00000 1 0
+1.00000 0.00000 0.00000 1 0
+
+1.00000 -1.00000 1.00000 1 0
+1.00000 0.00000 1.00000 1 0
+
+0.00000 -1.00000 1.00000 1 0
+0.00000 0.00000 1.00000 1 0
+
+-1.00000 0.00000 0.00000 1 0
+0.00000 0.00000 0.145833 1 0
+0.00000 0.00000 1.00000 1 0
+-1.00000 0.00000 1.00000 1 0
+-1.00000 0.00000 0.00000 1 0
+
+-1.00000 1.00000 0.00000 1 0
+0.00000 1.00000 0.00000 1 0
+0.00000 1.00000 1.00000 1 0
+-1.00000 1.00000 1.00000 1 0
+-1.00000 1.00000 0.00000 1 0
+
+-1.00000 0.00000 0.00000 1 0
+-1.00000 1.00000 0.00000 1 0
+
+0.00000 0.00000 0.145833 1 0
+0.00000 1.00000 0.00000 1 0
+
+0.00000 0.00000 1.00000 1 0
+0.00000 1.00000 1.00000 1 0
+
+-1.00000 0.00000 1.00000 1 0
+-1.00000 1.00000 1.00000 1 0
+
+0.00000 0.00000 0.145833 1 0
+1.00000 0.00000 0.00000 1 0
+1.00000 0.00000 1.00000 1 0
+0.00000 0.00000 1.00000 1 0
+0.00000 0.00000 0.145833 1 0
+
+0.00000 1.00000 0.00000 1 0
+1.00000 1.00000 0.00000 1 0
+1.00000 1.00000 1.00000 1 0
+0.00000 1.00000 1.00000 1 0
+0.00000 1.00000 0.00000 1 0
+
+0.00000 0.00000 0.145833 1 0
+0.00000 1.00000 0.00000 1 0
+
+1.00000 0.00000 0.00000 1 0
+1.00000 1.00000 0.00000 1 0
+
+1.00000 0.00000 1.00000 1 0
+1.00000 1.00000 1.00000 1 0
+
+0.00000 0.00000 1.00000 1 0
+0.00000 1.00000 1.00000 1 0
+

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.