]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Provide binary CellId representation 3289/head
authorRene Gassmoeller <rene.gassmoeller@mailbox.org>
Mon, 19 Sep 2016 17:30:41 +0000 (11:30 -0600)
committerRene Gassmoeller <rene.gassmoeller@mailbox.org>
Fri, 28 Oct 2016 23:28:39 +0000 (17:28 -0600)
include/deal.II/distributed/p4est_wrappers.h
include/deal.II/grid/cell_id.h
include/deal.II/grid/tria_accessor.h
source/distributed/p4est_wrappers.cc
source/grid/cell_id.cc
source/grid/cell_id.inst.in
source/grid/tria_accessor.cc
tests/grid/cell_id_04.cc [new file with mode: 0644]
tests/grid/cell_id_04.output [new file with mode: 0644]

index 67181ec79e9f00e49d03f3ae9d175c09153dbbb7..2f703c009f82afc5887dacbb830475733c02851d 100644 (file)
@@ -272,7 +272,7 @@ namespace internal
       static
       size_t (&connectivity_memory_used) (types<2>::connectivity *p4est);
 
-      static const unsigned max_level;
+      static const unsigned int max_level = P4EST_MAXLEVEL;
     };
 
 
@@ -454,7 +454,7 @@ namespace internal
       static
       size_t (&connectivity_memory_used) (types<3>::connectivity *p4est);
 
-      static const unsigned max_level;
+      static const unsigned int max_level = P8EST_MAXLEVEL;
     };
 
 
index a4bce635b75e86705ef9c56d83607761d989902a..4c818fa45ec720e894f805ee608017be13fe929d 100644 (file)
 
 #include <deal.II/base/config.h>
 #include <deal.II/base/exceptions.h>
+#include <deal.II/base/std_cxx11/array.h>
 
 #include <vector>
 #include <iostream>
 
+#ifdef DEAL_II_WITH_P4EST
+#include <deal.II/distributed/p4est_wrappers.h>
+#endif
 
 DEAL_II_NAMESPACE_OPEN
 
@@ -48,30 +52,67 @@ template <int, int> class Triangulation;
  *
  * @note How this data is internally represented is not of importance (and not
  * exposed on purpose).
- *
- * @todo Does it make sense to implement a more efficient representation
- * (internally and/or as a string)? If yes, something like a 64bit int as in
- * p4est would be a good option.
  */
 class CellId
 {
 public:
   /**
-   * Construct a CellId object with a given @p coarse_cell_index and list of child indices.
+   * A type that is used to encode the CellId data in a compact and fast way
+   * (e.g. for MPI transfer to other processes). Note that it limits the
+   * number of children that can be transferred to 20 in 3D and 30 in 2D
+   * (using 2 times 32 bit for storage), a limitation that is identical to
+   * the one used by p4est.
+   */
+  typedef std_cxx11::array<unsigned int, 4> binary_type;
+
+  /**
+   * Construct a CellId object with a given @p coarse_cell_id and vector of
+   * child indices. @p child_indices is
+   * interpreted identical to the member variable with the same name, namely
+   * each entry denotes which child to pick from one refinement level to the
+   * next, starting with the coarse cell, until we get to the cell represented
+   * by the current object. Therefore, each entry should be a number between 0
+   * and the number of children of a cell in the current space dimension.
    */
   CellId(const unsigned int coarse_cell_id,
          const std::vector<unsigned char> &child_indices);
 
+  /**
+   * Construct a CellId object with a given @p coarse_cell_id and array of
+   * child indices provided in @p child_indices. @p child_indices is
+   * interpreted identical to the member variable with the same name, namely
+   * each entry denotes which child to pick from one refinement level to the
+   * next, starting with the coarse cell, until we get to the cell represented
+   * by the current object. Therefore, each entry should be a number between 0
+   * and the number of children of a cell in the current space dimension.
+   * @p child_indices shall have at least @p n_child_indices valid entries.
+   */
+  CellId(const unsigned int coarse_cell_id,
+         const unsigned int n_child_indices,
+         const unsigned char *child_indices);
+
+  /**
+   * Construct a CellId object with a given binary representation that was
+   * previously constructed by CellId::to_binary.
+   */
+  CellId(const binary_type &binary_representation);
+
   /**
    * Construct an invalid CellId.
    */
   CellId();
 
   /**
-   * Return a string representation of this CellId.
+   * Return a human readable string representation of this CellId.
    */
   std::string to_string() const;
 
+  /**
+   * Return a compact and fast binary representation of this CellId.
+   */
+  template<int dim>
+  binary_type to_binary() const;
+
   /**
    * Return a cell_iterator to the cell represented by this CellId.
    */
@@ -104,11 +145,25 @@ private:
   unsigned int coarse_cell_id;
 
   /**
-   * A list of integers that denote which child to pick from one
+   * The number of child indices stored in the child_indices array. This is
+   * equivalent to (level-1) of the current cell.
+   */
+  unsigned int n_child_indices;
+
+  /**
+   * An array of integers that denotes which child to pick from one
    * refinement level to the next, starting with the coarse cell,
    * until we get to the cell represented by the current object.
+   * Only the first n_child_indices entries are used, but we use a statically
+   * allocated array instead of a vector of size n_child_indices to speed up
+   * creation of this object. If the given dimensions ever become a limitation
+   * the array can be extended.
    */
-  std::vector<unsigned char> child_indices;
+#ifdef DEAL_II_WITH_P4EST
+  std_cxx11::array<char,internal::p4est::functions<2>::max_level> child_indices;
+#else
+  std_cxx11::array<char,30> child_indices;
+#endif
 
   friend std::istream &operator>> (std::istream &is, CellId &cid);
   friend std::ostream &operator<< (std::ostream &os, const CellId &cid);
@@ -124,8 +179,8 @@ inline
 std::ostream &operator<< (std::ostream &os,
                           const CellId &cid)
 {
-  os << cid.coarse_cell_id << '_' << cid.child_indices.size() << ':';
-  for (unsigned int i=0; i<cid.child_indices.size(); ++i)
+  os << cid.coarse_cell_id << '_' << cid.n_child_indices << ':';
+  for (unsigned int i=0; i<cid.n_child_indices; ++i)
     // write the child indices. because they are between 0 and 2^dim-1, they all
     // just have one digit, so we could write them as integers. it's
     // probably clearer to write them as one-digit characters starting
@@ -152,47 +207,36 @@ std::istream &operator>> (std::istream &is,
   char dummy;
   is >> dummy;
   Assert(dummy=='_', ExcMessage("invalid CellId"));
-  unsigned int idsize;
-  is >> idsize;
+  is >> cid.n_child_indices;
   is >> dummy;
   Assert(dummy==':', ExcMessage("invalid CellId"));
 
   char value;
-  cid.child_indices.clear();
-  for (unsigned int i=0; i<idsize; ++i)
+  for (unsigned int i=0; i<cid.n_child_indices; ++i)
     {
       // read the one-digit child index (as an integer number) and
       // convert it back into unsigned char
       is >> value;
-      cid.child_indices.push_back(value-'0');
+      cid.child_indices[i] = value-'0';
     }
   return is;
 }
 
 
-inline
-CellId::CellId(const unsigned int coarse_cell_id,
-               const std::vector<unsigned char> &id)
-  :
-  coarse_cell_id(coarse_cell_id),
-  child_indices(id)
-{}
-
-
-inline
-CellId::CellId()
-  :
-  coarse_cell_id(static_cast<unsigned int>(-1))
-{}
-
-
 
 inline bool
 CellId::operator== (const CellId &other) const
 {
   if (this->coarse_cell_id != other.coarse_cell_id)
     return false;
-  return child_indices == other.child_indices;
+  if (n_child_indices != other.n_child_indices)
+    return false;
+
+  for (unsigned int i=0; i<n_child_indices; ++i)
+    if (child_indices[i] != other.child_indices[i])
+      return false;
+
+  return true;
 }
 
 
@@ -212,9 +256,9 @@ bool CellId::operator<(const CellId &other) const
     return this->coarse_cell_id < other.coarse_cell_id;
 
   unsigned int idx = 0;
-  while (idx < child_indices.size())
+  while (idx < n_child_indices)
     {
-      if (idx>=other.child_indices.size())
+      if (idx>=other.n_child_indices)
         return false;
 
       if (child_indices[idx] != other.child_indices[idx])
@@ -223,7 +267,7 @@ bool CellId::operator<(const CellId &other) const
       ++idx;
     }
 
-  if (child_indices.size() == other.child_indices.size())
+  if (n_child_indices == other.n_child_indices)
     return false;
   return true; // other.id is longer
 }
index c9822c4fcf42280562e81f26939e89e04f3b60eb..8bb37682bc45c21ab8e6ff39aa8d0ab9d2b3c0e7 100644 (file)
@@ -3234,39 +3234,6 @@ CellAccessor (const TriaAccessor<structdim2,dim2,spacedim2> &)
                       "the conversion is not valid in the current context."));
 }
 
-template <int dim, int spacedim>
-CellId
-CellAccessor<dim,spacedim>::id() const
-{
-  std::vector<unsigned char> id(this->level(), -1);
-  unsigned int coarse_index;
-
-  CellAccessor<dim,spacedim> ptr = *this;
-  while (ptr.level()>0)
-    {
-      // determine which child we are
-      unsigned char v = static_cast<unsigned char>(-1);
-      for (unsigned int c=0; c<ptr.parent()->n_children(); ++c)
-        {
-          if (ptr.parent()->child_index(c)==ptr.index())
-            {
-              v = c;
-              break;
-            }
-        }
-
-      Assert(v != (unsigned char)-1, ExcInternalError());
-      id[ptr.level()-1] = v;
-
-      ptr.copy_from( *(ptr.parent()));
-    }
-
-  Assert(ptr.level()==0, ExcInternalError());
-  coarse_index = ptr.index();
-
-  return CellId(coarse_index, id);
-}
-
 
 #ifndef DOXYGEN
 
index 77abbfcf4088dd83479a94c7ffd1e6342f56f5ed..d2e6d4c614cb0971d56ccab43f4b0c37e8884db2 100644 (file)
@@ -191,7 +191,7 @@ namespace internal
     size_t (&functions<2>::connectivity_memory_used) (types<2>::connectivity *p4est)
       = p4est_connectivity_memory_used;
 
-    const unsigned int functions<2>::max_level = P4EST_MAXLEVEL;
+    const unsigned int functions<2>::max_level;
 
     int (&functions<3>::quadrant_compare) (const void *v1, const void *v2)
       = p8est_quadrant_compare;
@@ -379,7 +379,7 @@ namespace internal
     size_t (&functions<3>::connectivity_memory_used) (types<3>::connectivity *p4est)
       = p8est_connectivity_memory_used;
 
-    const unsigned int functions<3>::max_level = P8EST_MAXLEVEL;
+    const unsigned int functions<3>::max_level;
 
     template <int dim>
     void
index b50878c0b32ca57bb41cadf722cb1a4ebfe7cab1..fa42f4952abc8d77fe31abbd0310681271715b8e 100644 (file)
 
 DEAL_II_NAMESPACE_OPEN
 
+
+CellId::CellId()
+  :
+  coarse_cell_id(numbers::invalid_unsigned_int),
+  n_child_indices(numbers::invalid_unsigned_int)
+{}
+
+
+
+CellId::CellId(const unsigned int coarse_cell_id,
+               const std::vector<unsigned char> &id)
+  :
+  coarse_cell_id(coarse_cell_id),
+  n_child_indices(id.size())
+{
+  Assert(n_child_indices < child_indices.size(),
+         ExcInternalError());
+  std::copy(id.begin(),id.end(),child_indices.begin());
+}
+
+
+
+CellId::CellId(const unsigned int coarse_cell_id,
+               const unsigned int n_child_indices,
+               const unsigned char *id)
+  :
+  coarse_cell_id(coarse_cell_id),
+  n_child_indices(n_child_indices)
+{
+  Assert(n_child_indices < child_indices.size(),
+         ExcInternalError());
+  memcpy(&(child_indices[0]),id,n_child_indices);
+}
+
+
+
+CellId::CellId(const CellId::binary_type &binary_representation)
+{
+  // The first entry stores the coarse cell id
+  coarse_cell_id = binary_representation[0];
+
+  // The rightmost two bits of the second entry store the dimension,
+  // the rest stores the number of child indices.
+  const unsigned int two_bit_mask = (1<<2)-1;
+  const unsigned int dim = binary_representation[1] & two_bit_mask;
+  n_child_indices = (binary_representation[1] >> 2);
+
+  Assert(n_child_indices < child_indices.size(),
+         ExcInternalError());
+
+  // Each child requires 'dim' bits to store its index
+  const unsigned int children_per_value = sizeof(binary_type::value_type) * 8 / dim;
+  const unsigned int child_mask = (1<<dim) - 1;
+
+  // Loop until all child indices have been read
+  unsigned int child_level=0;
+  unsigned int binary_entry = 2;
+  while (child_level < n_child_indices)
+    {
+      for (unsigned int j=0; j<children_per_value; ++j)
+        {
+          // Read the current child index by shifting to the current
+          // index's position and doing a bitwise-and with the child_mask.
+          child_indices[child_level] = (binary_representation[binary_entry] >> (dim*j)) & child_mask;
+          ++child_level;
+          if (child_level == n_child_indices)
+            break;
+        }
+      ++binary_entry;
+    }
+}
+
+
+
+template <int dim>
+CellId::binary_type
+CellId::to_binary() const
+{
+  CellId::binary_type binary_representation;
+  binary_representation.fill(0);
+
+  Assert(n_child_indices < child_indices.size(),
+         ExcInternalError());
+
+  // The first entry stores the coarse cell id
+  binary_representation[0] = coarse_cell_id;
+
+  // The rightmost two bits of the second entry store the dimension,
+  // the rest stores the number of child indices.
+  binary_representation[1] = (n_child_indices << 2);
+  binary_representation[1] |= dim;
+
+  // Each child requires 'dim' bits to store its index
+  const unsigned int children_per_value = sizeof(binary_type::value_type) * 8 / dim;
+  unsigned int child_level=0;
+  unsigned int binary_entry = 2;
+
+  // Loop until all child indices have been written
+  while (child_level < n_child_indices)
+    {
+      Assert(binary_entry < binary_representation.size(),
+             ExcInternalError());
+
+      for (unsigned int j=0; j<children_per_value; ++j)
+        {
+          const unsigned int child_index = static_cast<unsigned int>(child_indices[child_level]);
+          // Shift the child index to its position in the unsigned int and store it
+          binary_representation[binary_entry] |= (child_index << (j*dim));
+          ++child_level;
+          if (child_level == n_child_indices)
+            break;
+        }
+      ++binary_entry;
+    }
+
+  return binary_representation;
+}
+
+
+
 std::string
 CellId::to_string() const
 {
@@ -29,13 +149,15 @@ CellId::to_string() const
   return ss.str();
 }
 
-template<int dim, int spacedim>
+
+
+template <int dim, int spacedim>
 typename Triangulation<dim,spacedim>::cell_iterator
 CellId::to_cell(const Triangulation<dim,spacedim> &tria) const
 {
   typename Triangulation<dim,spacedim>::cell_iterator cell (&tria,0,coarse_cell_id);
 
-  for (unsigned int i = 0; i < child_indices.size(); ++i)
+  for (unsigned int i = 0; i < n_child_indices; ++i)
     cell = cell->child(static_cast<unsigned int> (child_indices[i]));
 
   return cell;
index 137c1a5faf8a4e54c7f18d5bf779d93f149e638a..61ff77fe9901646189bcebd4542399d4f61f25a7 100644 (file)
@@ -21,3 +21,8 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension :  SPACE_DIMENSIONS
     template Triangulation<deal_II_dimension,deal_II_space_dimension>::cell_iterator CellId::to_cell(const Triangulation<deal_II_dimension,deal_II_space_dimension> &tria) const;
 #endif
 }
+
+for (deal_II_dimension : DIMENSIONS)
+{
+    template CellId::binary_type CellId::to_binary<deal_II_dimension>() const;
+}
index 4f3b53a30cf8de68aec68eb8f04f0ba3fdf036f2..cfcc932819d83ff24768738156801ad90e4da589 100644 (file)
@@ -19,6 +19,7 @@
 #include <deal.II/grid/tria_accessor.h>
 #include <deal.II/grid/tria_accessor.templates.h>
 #include <deal.II/grid/tria_iterator.templates.h>
+#include <deal.II/base/std_cxx11/array.h>
 #include <deal.II/base/geometry_info.h>
 #include <deal.II/base/quadrature.h>
 #include <deal.II/grid/grid_tools.h>
@@ -1635,6 +1636,45 @@ void CellAccessor<dim, spacedim>::set_neighbor (const unsigned int i,
 
 
 
+template <int dim, int spacedim>
+CellId
+CellAccessor<dim,spacedim>::id() const
+{
+  std_cxx11::array<unsigned char,30> id;
+
+  CellAccessor<dim,spacedim> ptr = *this;
+  const unsigned int n_child_indices = ptr.level();
+
+  while (ptr.level()>0)
+    {
+      const TriaIterator<CellAccessor<dim,spacedim> > parent = ptr.parent();
+      const unsigned int n_children = parent->n_children();
+
+      // determine which child we are
+      unsigned char v = static_cast<unsigned char> (-1);
+      for (unsigned int c=0; c<n_children; ++c)
+        {
+          if (parent->child_index(c)==ptr.index())
+            {
+              v = c;
+              break;
+            }
+        }
+
+      Assert(v != static_cast<unsigned char> (-1), ExcInternalError());
+      id[ptr.level()-1] = v;
+
+      ptr.copy_from(*parent);
+    }
+
+  Assert(ptr.level()==0, ExcInternalError());
+  const unsigned int coarse_index = ptr.index();
+
+  return CellId(coarse_index, n_child_indices, &(id[0]));
+}
+
+
+
 template <int dim, int spacedim>
 unsigned int CellAccessor<dim,spacedim>::neighbor_of_neighbor_internal (const unsigned int neighbor) const
 {
diff --git a/tests/grid/cell_id_04.cc b/tests/grid/cell_id_04.cc
new file mode 100644 (file)
index 0000000..340e38c
--- /dev/null
@@ -0,0 +1,94 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2016 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 at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+
+
+// check CellId
+
+#include "../tests.h"
+
+#include <deal.II/base/geometry_info.h>
+#include <deal.II/base/logstream.h>
+#include <deal.II/base/quadrature_lib.h>
+#include <deal.II/dofs/dof_handler.h>
+#include <deal.II/dofs/dof_accessor.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/distributed/tria.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_refinement.h>
+
+#include <fstream>
+#include <sstream>
+
+using namespace dealii;
+
+template <int dim>
+void check (Triangulation<dim> &tr)
+{
+  typename Triangulation<dim>::cell_iterator cell = tr.begin(),
+                                             endc = tr.end();
+
+
+  for (; cell!=endc; ++cell)
+    {
+      deallog << cell->level() << " " << cell->index() << std::endl;
+
+      // Store the CellId, convert it to a binary representation,
+      // create a new CellId from that, and create a cell iterator
+      // pointing to the same cell
+
+      const CellId cid = cell->id();
+
+      const CellId::binary_type cids = cid.to_binary<dim>();
+
+      const CellId cid2(cids);
+
+      typename Triangulation<dim>::cell_iterator cell2 = cid2.to_cell(tr);
+
+      Assert(cid2 == cid, ExcInternalError());
+      Assert(cell2 == cell, ExcInternalError());
+
+      deallog << cell2->level() << " " << cell2->index() << std::endl;
+    }
+
+  deallog << "OK" << std::endl;
+}
+
+
+int main (int argc, char *argv[])
+{
+  //Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv, testing_max_num_threads());
+
+  initlog();
+  deal_II_exceptions::disable_abort_on_exception();
+
+  Triangulation<2> tria;
+  GridGenerator::hyper_cube (tria);
+  tria.refine_global (3);
+  tria.begin_active()->set_refine_flag();
+  tria.execute_coarsening_and_refinement();
+  check(tria);
+
+  Triangulation<3> tria2;
+  GridGenerator::hyper_cube (tria2);
+  tria2.refine_global (1);
+  tria2.begin_active()->set_refine_flag();
+  tria2.execute_coarsening_and_refinement();
+  check(tria2);
+}
+
+
+
diff --git a/tests/grid/cell_id_04.output b/tests/grid/cell_id_04.output
new file mode 100644 (file)
index 0000000..3e722c6
--- /dev/null
@@ -0,0 +1,215 @@
+
+DEAL::0 0
+DEAL::0 0
+DEAL::1 0
+DEAL::1 0
+DEAL::1 1
+DEAL::1 1
+DEAL::1 2
+DEAL::1 2
+DEAL::1 3
+DEAL::1 3
+DEAL::2 0
+DEAL::2 0
+DEAL::2 1
+DEAL::2 1
+DEAL::2 2
+DEAL::2 2
+DEAL::2 3
+DEAL::2 3
+DEAL::2 4
+DEAL::2 4
+DEAL::2 5
+DEAL::2 5
+DEAL::2 6
+DEAL::2 6
+DEAL::2 7
+DEAL::2 7
+DEAL::2 8
+DEAL::2 8
+DEAL::2 9
+DEAL::2 9
+DEAL::2 10
+DEAL::2 10
+DEAL::2 11
+DEAL::2 11
+DEAL::2 12
+DEAL::2 12
+DEAL::2 13
+DEAL::2 13
+DEAL::2 14
+DEAL::2 14
+DEAL::2 15
+DEAL::2 15
+DEAL::3 0
+DEAL::3 0
+DEAL::3 1
+DEAL::3 1
+DEAL::3 2
+DEAL::3 2
+DEAL::3 3
+DEAL::3 3
+DEAL::3 4
+DEAL::3 4
+DEAL::3 5
+DEAL::3 5
+DEAL::3 6
+DEAL::3 6
+DEAL::3 7
+DEAL::3 7
+DEAL::3 8
+DEAL::3 8
+DEAL::3 9
+DEAL::3 9
+DEAL::3 10
+DEAL::3 10
+DEAL::3 11
+DEAL::3 11
+DEAL::3 12
+DEAL::3 12
+DEAL::3 13
+DEAL::3 13
+DEAL::3 14
+DEAL::3 14
+DEAL::3 15
+DEAL::3 15
+DEAL::3 16
+DEAL::3 16
+DEAL::3 17
+DEAL::3 17
+DEAL::3 18
+DEAL::3 18
+DEAL::3 19
+DEAL::3 19
+DEAL::3 20
+DEAL::3 20
+DEAL::3 21
+DEAL::3 21
+DEAL::3 22
+DEAL::3 22
+DEAL::3 23
+DEAL::3 23
+DEAL::3 24
+DEAL::3 24
+DEAL::3 25
+DEAL::3 25
+DEAL::3 26
+DEAL::3 26
+DEAL::3 27
+DEAL::3 27
+DEAL::3 28
+DEAL::3 28
+DEAL::3 29
+DEAL::3 29
+DEAL::3 30
+DEAL::3 30
+DEAL::3 31
+DEAL::3 31
+DEAL::3 32
+DEAL::3 32
+DEAL::3 33
+DEAL::3 33
+DEAL::3 34
+DEAL::3 34
+DEAL::3 35
+DEAL::3 35
+DEAL::3 36
+DEAL::3 36
+DEAL::3 37
+DEAL::3 37
+DEAL::3 38
+DEAL::3 38
+DEAL::3 39
+DEAL::3 39
+DEAL::3 40
+DEAL::3 40
+DEAL::3 41
+DEAL::3 41
+DEAL::3 42
+DEAL::3 42
+DEAL::3 43
+DEAL::3 43
+DEAL::3 44
+DEAL::3 44
+DEAL::3 45
+DEAL::3 45
+DEAL::3 46
+DEAL::3 46
+DEAL::3 47
+DEAL::3 47
+DEAL::3 48
+DEAL::3 48
+DEAL::3 49
+DEAL::3 49
+DEAL::3 50
+DEAL::3 50
+DEAL::3 51
+DEAL::3 51
+DEAL::3 52
+DEAL::3 52
+DEAL::3 53
+DEAL::3 53
+DEAL::3 54
+DEAL::3 54
+DEAL::3 55
+DEAL::3 55
+DEAL::3 56
+DEAL::3 56
+DEAL::3 57
+DEAL::3 57
+DEAL::3 58
+DEAL::3 58
+DEAL::3 59
+DEAL::3 59
+DEAL::3 60
+DEAL::3 60
+DEAL::3 61
+DEAL::3 61
+DEAL::3 62
+DEAL::3 62
+DEAL::3 63
+DEAL::3 63
+DEAL::4 0
+DEAL::4 0
+DEAL::4 1
+DEAL::4 1
+DEAL::4 2
+DEAL::4 2
+DEAL::4 3
+DEAL::4 3
+DEAL::OK
+DEAL::0 0
+DEAL::0 0
+DEAL::1 0
+DEAL::1 0
+DEAL::1 1
+DEAL::1 1
+DEAL::1 2
+DEAL::1 2
+DEAL::1 3
+DEAL::1 3
+DEAL::1 4
+DEAL::1 4
+DEAL::1 5
+DEAL::1 5
+DEAL::1 6
+DEAL::1 6
+DEAL::1 7
+DEAL::1 7
+DEAL::2 0
+DEAL::2 0
+DEAL::2 1
+DEAL::2 1
+DEAL::2 2
+DEAL::2 2
+DEAL::2 3
+DEAL::2 3
+DEAL::2 4
+DEAL::2 4
+DEAL::2 5
+DEAL::2 5
+DEAL::2 6
+DEAL::2 6
+DEAL::2 7
+DEAL::2 7
+DEAL::OK

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.