]> https://gitweb.dealii.org/ - dealii.git/commitdiff
new functionality to compute mesh aspect ratio 8894/head
authorNiklas Fehn <fehn@lnm.mw.tum.de>
Tue, 8 Oct 2019 11:10:56 +0000 (13:10 +0200)
committerNiklas Fehn <fehn@lnm.mw.tum.de>
Mon, 14 Oct 2019 06:29:19 +0000 (08:29 +0200)
doc/news/changes/minor/20191009NiklasFehn [new file with mode: 0644]
include/deal.II/grid/grid_tools.h
source/grid/grid_tools.cc
source/grid/grid_tools.inst.in
tests/grid/grid_tools_aspect_ratio.cc [new file with mode: 0644]
tests/grid/grid_tools_aspect_ratio.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20191009NiklasFehn b/doc/news/changes/minor/20191009NiklasFehn
new file mode 100644 (file)
index 0000000..fe96e57
--- /dev/null
@@ -0,0 +1,4 @@
+New: The function GridTools::compute_maximum_aspect_ratio() computes the
+maximum aspect ratio of a mesh with arbitrarily deformed elements.
+<br>
+(Niklas Fehn, 2019/10/08)
index 307c305dbe73f3b43f999976680522926e6a31dc..8cb2bcae5bd2dd2015bc02377f7d7cf487de027b 100644 (file)
@@ -230,6 +230,55 @@ namespace GridTools
   double
   cell_measure(const T &, ...);
 
+  /**
+   * Computes an aspect ratio measure for all locally-owned active cells and
+   * fills a vector with one entry per cell, given a @p triangulation and
+   * @p mapping. The size of the vector that is returned equals the number of
+   * active cells. The vector contains zero for non locally-owned cells. The
+   * aspect ratio of a cell is defined as the ratio of the maximum to minimum
+   * singular value of the Jacobian, taking the maximum over all quadrature
+   * points of a quadrature rule specified via @p quadrature. For example, for
+   * the special case of rectangular elements in 2d with dimensions $a$ and $b$
+   * ($a \geq b$), this function returns the usual aspect ratio definition
+   * $a/b$. The above definition using singular values is a generalization to
+   * arbitrarily deformed elements. This function is intended to be used for
+   * $d=2,3$ space dimensions, but it can also be used for $d=1$ returning a
+   * value of 1.
+   *
+   * @note Inverted elements do not throw an exception. Instead, a value of inf
+   * is written into the vector in case of inverted elements.
+   *
+   * @note Make sure to use enough quadrature points for a precise calculation
+   * of the aspect ratio in case of deformed elements.
+   *
+   * @note In parallel computations the return value will have the length
+   * n_active_cells but the aspect ratio is only computed for the cells that
+   * are locally owned and placed at index CellAccessor::active_cell_index(),
+   * respectively. All other values are set to 0.
+   *
+   * @author Niklas Fehn, Martin Kronbichler, 2019
+   */
+  template <int dim>
+  Vector<double>
+  compute_aspect_ratio_of_cells(const Triangulation<dim> &triangulation,
+                                const Mapping<dim> &      mapping,
+                                const Quadrature<dim> &   quadrature);
+
+  /**
+   * Computes the maximum aspect ratio by taking the maximum over all cells.
+   *
+   * @note When running in parallel with a Triangulation that supports MPI,
+   * this is a collective call and the return value is the maximum over all
+   * processors.
+   *
+   * @author Niklas Fehn, Martin Kronbichler, 2019
+   */
+  template <int dim>
+  double
+  compute_maximum_aspect_ratio(const Triangulation<dim> &triangulation,
+                               const Mapping<dim> &      mapping,
+                               const Quadrature<dim> &   quadrature);
+
   /**
    * Compute the smallest box containing the entire triangulation.
    *
index bbbaeecc2099330cbbc29c08aa9aca8fd6e5774c..b824ba319290282840f8ed470ff5bca2075c9c0b 100644 (file)
@@ -52,6 +52,7 @@
 #include <deal.II/lac/vector_memory.h>
 
 #include <deal.II/numerics/matrix_tools.h>
+#include <deal.II/numerics/vector_tools.h>
 
 #include <boost/random/mersenne_twister.hpp>
 #include <boost/random/uniform_real_distribution.hpp>
@@ -414,6 +415,79 @@ namespace GridTools
   }
 
 
+  template <int dim>
+  Vector<double>
+  compute_aspect_ratio_of_cells(const Triangulation<dim> &triangulation,
+                                const Mapping<dim> &      mapping,
+                                const Quadrature<dim> &   quadrature)
+  {
+    FE_Nothing<dim> fe;
+    FEValues<dim>   fe_values(mapping, fe, quadrature, update_jacobians);
+
+    Vector<double> aspect_ratio_vector(triangulation.n_active_cells());
+
+    // loop over cells of processor
+    for (const auto &cell : triangulation.active_cell_iterators())
+      {
+        if (cell->is_locally_owned())
+          {
+            double aspect_ratio_cell = 0.0;
+
+            fe_values.reinit(cell);
+
+            // loop over quadrature points
+            for (unsigned int q = 0; q < quadrature.size(); ++q)
+              {
+                Tensor<2, dim, double> jacobian =
+                  Tensor<2, dim, double>(fe_values.jacobian(q));
+
+                LAPACKFullMatrix<double> J = LAPACKFullMatrix<double>(dim);
+                for (unsigned int i = 0; i < dim; i++)
+                  for (unsigned int j = 0; j < dim; j++)
+                    J(i, j) = jacobian[i][j];
+
+                // We intentionally do not want to throw an exception in case of
+                // inverted elements since this is not the task of this
+                // function. Instead, inf is written into the vector in case of
+                // inverted elements.
+                if (determinant(jacobian) <= 0)
+                  {
+                    aspect_ratio_cell = std::numeric_limits<double>::infinity();
+                  }
+                else
+                  {
+                    J.compute_svd();
+
+                    double const max_sv = J.singular_value(0);
+                    double const min_sv = J.singular_value(dim - 1);
+                    double const ar     = max_sv / min_sv;
+
+                    aspect_ratio_cell = std::max(aspect_ratio_cell, ar);
+                  }
+              }
+
+            // fill vector
+            aspect_ratio_vector(cell->active_cell_index()) = aspect_ratio_cell;
+          }
+      }
+
+    return aspect_ratio_vector;
+  }
+
+  template <int dim>
+  double
+  compute_maximum_aspect_ratio(const Triangulation<dim> &triangulation,
+                               const Mapping<dim> &      mapping,
+                               const Quadrature<dim> &   quadrature)
+  {
+    Vector<double> aspect_ratio_vector =
+      compute_aspect_ratio_of_cells(triangulation, mapping, quadrature);
+
+    return VectorTools::compute_global_error(triangulation,
+                                             aspect_ratio_vector,
+                                             VectorTools::Linfty_norm);
+  }
+
 
   template <int dim, int spacedim>
   BoundingBox<spacedim>
index 900a1caf993c23eed49061d8c288afe21741c47c..edf8ab5280d61fe17c4eaf43aaaa59e3739151c2 100644 (file)
@@ -175,6 +175,16 @@ for (deal_II_space_dimension : SPACE_DIMENSIONS)
       std::pair<BoundingBox<deal_II_space_dimension>, unsigned int>>
     GridTools::build_global_description_tree(
       const std::vector<BoundingBox<deal_II_space_dimension>> &, MPI_Comm);
+
+    template Vector<double> GridTools::compute_aspect_ratio_of_cells(
+      const Triangulation<deal_II_space_dimension> &,
+      const Mapping<deal_II_space_dimension> &,
+      const Quadrature<deal_II_space_dimension> &);
+
+    template double GridTools::compute_maximum_aspect_ratio(
+      const Triangulation<deal_II_space_dimension> &,
+      const Mapping<deal_II_space_dimension> &,
+      const Quadrature<deal_II_space_dimension> &);
   }
 
 
diff --git a/tests/grid/grid_tools_aspect_ratio.cc b/tests/grid/grid_tools_aspect_ratio.cc
new file mode 100644 (file)
index 0000000..54507c9
--- /dev/null
@@ -0,0 +1,234 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2019 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.
+//
+// ---------------------------------------------------------------------
+
+// Short test to validate GridTools::compute_maximum_aspect_ratio()
+
+#include <deal.II/base/exceptions.h>
+#include <deal.II/base/mpi.h>
+#include <deal.II/base/partitioner.h>
+#include <deal.II/base/quadrature_lib.h>
+#include <deal.II/base/utilities.h>
+
+#include <deal.II/distributed/tria.h>
+
+#include <deal.II/dofs/dof_handler.h>
+
+#include <deal.II/fe/fe_dgq.h>
+#include <deal.II/fe/mapping_q_generic.h>
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_tools.h>
+
+#include "../tests.h"
+
+using namespace dealii;
+
+/*
+ * The parameters left, right, refinements follow the hyper-rectangle
+ * notation.
+ *
+ * degree denotes the polynomial degree of the mapping.
+ *
+ * n_q_points denotes the number of 1D quadrature points used during
+ * the computation of the aspect ratio.
+ *
+ * If deform is true, the first vertex of the first element of the
+ * hyper_rectangle is shifted, with factor being a scaling factor to
+ * be able to generate inverted elements.
+ */
+template <int dim>
+double
+compute_aspect_ratio_hyper_rectangle(
+  Point<dim> const &               left,
+  Point<dim> const &               right,
+  std::vector<unsigned int> const &refinements,
+  unsigned int                     degree     = 1,
+  unsigned int                     n_q_points = 2,
+  bool                             deform     = false,
+  double                           factor     = 1.0)
+{
+  Triangulation<dim> tria;
+  GridGenerator::subdivided_hyper_rectangle(
+    tria, refinements, left, right, false);
+
+  if (deform)
+    {
+      Tensor<1, dim> diag = right - left;
+      double         l    = diag.norm();
+      Point<dim>     shift;
+      for (unsigned int d = 0; d < dim; ++d)
+        shift[d] = l * factor * (0.05 + d * 0.01);
+      tria.begin_active()->vertex(0) += shift;
+    }
+
+  MappingQGeneric<dim> const mapping(degree);
+  QGauss<dim> const          gauss(n_q_points);
+
+  Vector<double> ratios =
+    GridTools::compute_aspect_ratio_of_cells(tria, mapping, gauss);
+
+  deallog << std::endl
+          << "Parameters:"
+          << " d = " << dim << ", degree = " << degree
+          << ", n_q_points = " << n_q_points << ", deform = " << deform
+          << ", factor = " << factor << std::endl;
+
+  deallog << "Aspect ratio vector = " << std::endl;
+  ratios.print(deallog.get_file_stream());
+
+  return GridTools::compute_maximum_aspect_ratio(tria, mapping, gauss);
+}
+
+int
+main(int argc, char **argv)
+{
+  initlog();
+  deallog << std::setprecision(6);
+
+  try
+    {
+      Utilities::MPI::MPI_InitFinalize mpi(argc, argv, 1);
+
+      // 1D
+      {
+        Point<1> left  = Point<1>(0.0);
+        Point<1> right = Point<1>(1.0);
+
+        double ar = 0.0;
+
+        std::vector<unsigned int> refine(1, 2);
+
+        ar = compute_aspect_ratio_hyper_rectangle(left, right, refine, 1, 2);
+        deallog << "aspect ratio max    = " << ar << std::endl << std::endl;
+
+        refine[0] = 5;
+
+        ar = compute_aspect_ratio_hyper_rectangle(left, right, refine, 1, 2);
+        deallog << "aspect ratio max    = " << ar << std::endl << std::endl;
+      }
+
+      // 2D
+      {
+        Point<2> left  = Point<2>(0.0, 0.0);
+        Point<2> right = Point<2>(1.0, 1.0);
+
+        double ar = 0.0;
+
+        std::vector<unsigned int> refine(2, 2);
+
+        ar = compute_aspect_ratio_hyper_rectangle(left, right, refine, 1, 2);
+        deallog << "aspect ratio max    = " << ar << std::endl << std::endl;
+
+        ar =
+          compute_aspect_ratio_hyper_rectangle(left, right, refine, 1, 2, true);
+        deallog << "aspect ratio max    = " << ar << std::endl << std::endl;
+
+        ar = compute_aspect_ratio_hyper_rectangle(
+          left, right, refine, 1, 2, true, 10);
+        deallog << "aspect ratio max    = " << ar << std::endl << std::endl;
+
+        refine[0] = 5;
+
+        ar = compute_aspect_ratio_hyper_rectangle(left, right, refine, 1, 2);
+        deallog << "aspect ratio max    = " << ar << std::endl << std::endl;
+
+        ar =
+          compute_aspect_ratio_hyper_rectangle(left, right, refine, 1, 2, true);
+        deallog << "aspect ratio max    = " << ar << std::endl << std::endl;
+
+        ar =
+          compute_aspect_ratio_hyper_rectangle(left, right, refine, 1, 5, true);
+        deallog << "aspect ratio max    = " << ar << std::endl << std::endl;
+
+        ar = compute_aspect_ratio_hyper_rectangle(
+          left, right, refine, 1, 15, true);
+        deallog << "aspect ratio max    = " << ar << std::endl << std::endl;
+
+        ar = compute_aspect_ratio_hyper_rectangle(
+          left, right, refine, 1, 2, true, 10);
+        deallog << "aspect ratio max    = " << ar << std::endl << std::endl;
+      }
+
+      // 3D
+      {
+        Point<3> left  = Point<3>(0.0, 0.0, 0.0);
+        Point<3> right = Point<3>(1.0, 1.0, 1.0);
+
+        double ar = 0.0;
+
+        std::vector<unsigned int> refine(3, 2);
+
+        ar = compute_aspect_ratio_hyper_rectangle(left, right, refine, 1, 2);
+        deallog << "aspect ratio max    = " << ar << std::endl << std::endl;
+
+        ar =
+          compute_aspect_ratio_hyper_rectangle(left, right, refine, 1, 2, true);
+        deallog << "aspect ratio max    = " << ar << std::endl << std::endl;
+
+        ar = compute_aspect_ratio_hyper_rectangle(
+          left, right, refine, 1, 2, true, 10);
+        deallog << "aspect ratio max    = " << ar << std::endl << std::endl;
+
+        refine[0] = 5;
+        refine[1] = 3;
+
+        ar = compute_aspect_ratio_hyper_rectangle(left, right, refine, 1, 2);
+        deallog << "aspect ratio max    = " << ar << std::endl << std::endl;
+
+        ar =
+          compute_aspect_ratio_hyper_rectangle(left, right, refine, 1, 2, true);
+        deallog << "aspect ratio max    = " << ar << std::endl << std::endl;
+
+        ar =
+          compute_aspect_ratio_hyper_rectangle(left, right, refine, 1, 5, true);
+        deallog << "aspect ratio max    = " << ar << std::endl << std::endl;
+
+        ar = compute_aspect_ratio_hyper_rectangle(
+          left, right, refine, 1, 15, true);
+        deallog << "aspect ratio max    = " << ar << std::endl << std::endl;
+
+        ar = compute_aspect_ratio_hyper_rectangle(
+          left, right, refine, 1, 2, true, 10);
+        deallog << "aspect ratio max    = " << ar << std::endl << std::endl;
+      }
+    }
+  catch (std::exception &exc)
+    {
+      std::cerr << std::endl
+                << std::endl
+                << "----------------------------------------------------"
+                << std::endl;
+      std::cerr << "Exception on processing: " << std::endl
+                << exc.what() << std::endl
+                << "Aborting!" << std::endl
+                << "----------------------------------------------------"
+                << std::endl;
+      return 1;
+    }
+  catch (...)
+    {
+      std::cerr << std::endl
+                << std::endl
+                << "----------------------------------------------------"
+                << std::endl;
+      std::cerr << "Unknown exception!" << std::endl
+                << "Aborting!" << std::endl
+                << "----------------------------------------------------"
+                << std::endl;
+      return 1;
+    }
+
+  return 0;
+}
diff --git a/tests/grid/grid_tools_aspect_ratio.output b/tests/grid/grid_tools_aspect_ratio.output
new file mode 100644 (file)
index 0000000..395806a
--- /dev/null
@@ -0,0 +1,109 @@
+
+DEAL::
+DEAL::Parameters: d = 1, degree = 1, n_q_points = 2, deform = 0, factor = 1.00000
+DEAL::Aspect ratio vector = 
+1.000e+00 1.000e+00 
+DEAL::aspect ratio max    = 1.00000
+DEAL::
+DEAL::
+DEAL::Parameters: d = 1, degree = 1, n_q_points = 2, deform = 0, factor = 1.00000
+DEAL::Aspect ratio vector = 
+1.000e+00 1.000e+00 1.000e+00 1.000e+00 1.000e+00 
+DEAL::aspect ratio max    = 1.00000
+DEAL::
+DEAL::
+DEAL::Parameters: d = 2, degree = 1, n_q_points = 2, deform = 0, factor = 1.00000
+DEAL::Aspect ratio vector = 
+1.000e+00 1.000e+00 1.000e+00 1.000e+00 
+DEAL::aspect ratio max    = 1.00000
+DEAL::
+DEAL::
+DEAL::Parameters: d = 2, degree = 1, n_q_points = 2, deform = 1, factor = 1.00000
+DEAL::Aspect ratio vector = 
+1.327e+00 1.000e+00 1.000e+00 1.000e+00 
+DEAL::aspect ratio max    = 1.32670
+DEAL::
+DEAL::
+DEAL::Parameters: d = 2, degree = 1, n_q_points = 2, deform = 1, factor = 10.0000
+DEAL::Aspect ratio vector = 
+inf 1.000e+00 1.000e+00 1.000e+00 
+DEAL::aspect ratio max    = inf
+DEAL::
+DEAL::
+DEAL::Parameters: d = 2, degree = 1, n_q_points = 2, deform = 0, factor = 1.00000
+DEAL::Aspect ratio vector = 
+2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 
+DEAL::aspect ratio max    = 2.50000
+DEAL::
+DEAL::
+DEAL::Parameters: d = 2, degree = 1, n_q_points = 2, deform = 1, factor = 1.00000
+DEAL::Aspect ratio vector = 
+3.476e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 
+DEAL::aspect ratio max    = 3.47552
+DEAL::
+DEAL::
+DEAL::Parameters: d = 2, degree = 1, n_q_points = 5, deform = 1, factor = 1.00000
+DEAL::Aspect ratio vector = 
+3.866e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 
+DEAL::aspect ratio max    = 3.86561
+DEAL::
+DEAL::
+DEAL::Parameters: d = 2, degree = 1, n_q_points = 15, deform = 1, factor = 1.00000
+DEAL::Aspect ratio vector = 
+3.971e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 
+DEAL::aspect ratio max    = 3.97069
+DEAL::
+DEAL::
+DEAL::Parameters: d = 2, degree = 1, n_q_points = 2, deform = 1, factor = 10.0000
+DEAL::Aspect ratio vector = 
+inf 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 
+DEAL::aspect ratio max    = inf
+DEAL::
+DEAL::
+DEAL::Parameters: d = 3, degree = 1, n_q_points = 2, deform = 0, factor = 1.00000
+DEAL::Aspect ratio vector = 
+1.000e+00 1.000e+00 1.000e+00 1.000e+00 1.000e+00 1.000e+00 1.000e+00 1.000e+00 
+DEAL::aspect ratio max    = 1.00000
+DEAL::
+DEAL::
+DEAL::Parameters: d = 3, degree = 1, n_q_points = 2, deform = 1, factor = 1.00000
+DEAL::Aspect ratio vector = 
+1.641e+00 1.000e+00 1.000e+00 1.000e+00 1.000e+00 1.000e+00 1.000e+00 1.000e+00 
+DEAL::aspect ratio max    = 1.64083
+DEAL::
+DEAL::
+DEAL::Parameters: d = 3, degree = 1, n_q_points = 2, deform = 1, factor = 10.0000
+DEAL::Aspect ratio vector = 
+inf 1.000e+00 1.000e+00 1.000e+00 1.000e+00 1.000e+00 1.000e+00 1.000e+00 
+DEAL::aspect ratio max    = inf
+DEAL::
+DEAL::
+DEAL::Parameters: d = 3, degree = 1, n_q_points = 2, deform = 0, factor = 1.00000
+DEAL::Aspect ratio vector = 
+2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 
+DEAL::aspect ratio max    = 2.50000
+DEAL::
+DEAL::
+DEAL::Parameters: d = 3, degree = 1, n_q_points = 2, deform = 1, factor = 1.00000
+DEAL::Aspect ratio vector = 
+4.611e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 
+DEAL::aspect ratio max    = 4.61134
+DEAL::
+DEAL::
+DEAL::Parameters: d = 3, degree = 1, n_q_points = 5, deform = 1, factor = 1.00000
+DEAL::Aspect ratio vector = 
+1.613e+01 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 
+DEAL::aspect ratio max    = 16.1349
+DEAL::
+DEAL::
+DEAL::Parameters: d = 3, degree = 1, n_q_points = 15, deform = 1, factor = 1.00000
+DEAL::Aspect ratio vector = 
+6.696e+01 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 
+DEAL::aspect ratio max    = 66.9573
+DEAL::
+DEAL::
+DEAL::Parameters: d = 3, degree = 1, n_q_points = 2, deform = 1, factor = 10.0000
+DEAL::Aspect ratio vector = 
+inf 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 2.500e+00 
+DEAL::aspect ratio max    = inf
+DEAL::
\ No newline at end of file

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.