]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Implement various functions for 1d.
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Sun, 29 Apr 2001 16:03:55 +0000 (16:03 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Sun, 29 Apr 2001 16:03:55 +0000 (16:03 +0000)
git-svn-id: https://svn.dealii.org/trunk@4501 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/dofs/dof_tools.h
deal.II/deal.II/include/numerics/vectors.h
deal.II/deal.II/source/dofs/dof_handler.cc
deal.II/deal.II/source/dofs/dof_tools.cc
deal.II/deal.II/source/numerics/vectors.cc
tests/deal.II/boundaries.cc [new file with mode: 0644]
tests/deal.II/boundaries.checked [new file with mode: 0644]

index 7a872b17bf29e177250fd5a118ed012fc202387b..8dc9c6442b41372307a6935bd7c4f5dd38fdefa9 100644 (file)
@@ -88,7 +88,7 @@ template <int dim> class Mapping;
  * should not rely on it since it may be changed sometimes: we loop
  * over all faces, check whether it is on the boundary, if so get the
  * global numbers of the degrees of freedom on that face, and for each
- * of these we give a subsequent boundary number if none has already
+ * of these we give a sequential boundary number if none has already
  * been given to this dof. But it should be emphasized again that you
  * should not try to use this internal knowledge about the used
  * algorithm, you are better off if you just accept the mapping `as
@@ -119,7 +119,7 @@ template <int dim> class Mapping;
  * easy.  It should really be a degree of freedom of which the
  * respective basis function has nonzero values on the boundary. At
  * least for Lagrange elements this definition is equal to the
- * statement that the off-point of the trial function, i.e. the point
+ * statement that the off-point of the shape function, i.e. the point
  * where the function assumes its nominal value (for Lagrange elements
  * this is the point where it has the function value @p{1}), is
  * located on the boundary. We do not check this directly, the
@@ -821,7 +821,7 @@ class DoFTools
                                      * index of that degree of
                                      * freedom on the boundary. After
                                      * this operation, @p{mapping[dof]}
-                                     * gives the index of the the
+                                     * gives the index of the
                                      * degree of freedom with global
                                      * number @p{dof} in the list of
                                      * degrees of freedom on the
@@ -841,13 +841,6 @@ class DoFTools
                                      *
                                      * Prior content of @p{mapping}
                                      * is deleted.
-                                     *
-                                     * This function is not
-                                     * implemented for one
-                                     * dimension. See the general doc
-                                     * of this class for more
-                                     * information on boundary
-                                     * treatment.
                                      */
     template <int dim>
     static void
index 0c556f34ad64fd5ed8e533b47ce7ea98f481d264..341fb804a6f5212c811532f211acd956e17384cd 100644 (file)
@@ -542,7 +542,6 @@ class VectorTools
                                             const std::vector<bool>       &component_mask = std::vector<bool>());
 
     
-//TODO:[WB] Update project_boundary_values for more components
                                     /**
                                      * Project @p{function} to the boundary
                                      * of the domain, using the given quadrature
index 31c0d1a34ead292d65b89fa016feb845fd71c78f..f263da4fa831ce47bc6ccadbcbc562a6e8b8c86f 100644 (file)
@@ -1154,8 +1154,7 @@ template <>
 unsigned int DoFHandler<1>::n_boundary_dofs (const FunctionMap &) const
 {
   Assert (selected_fe != 0, ExcNoFESelected());
-  Assert (false, ExcNotImplemented());
-  return 0;
+  return 2*selected_fe->dofs_per_vertex;
 };
 
 
@@ -1659,8 +1658,7 @@ template <>
 unsigned int DoFHandler<1>::max_couplings_between_boundary_dofs () const
 {
   Assert (selected_fe != 0, ExcNoFESelected());
-  Assert (false, ExcInternalError());
-  return 0;
+  return selected_fe->dofs_per_vertex;
 };
 
 #endif
index d92f5bc05a3e128c2e0e43377ecf9c3520932440..14ccdeacf84338ffecaee98f24180e69a8137f9d 100644 (file)
@@ -118,23 +118,65 @@ DoFTools::make_sparsity_pattern (const DoFHandler<dim>                 &dof,
 #if deal_II_dimension == 1
 
 template <>
-void DoFTools::make_boundary_sparsity_pattern (const DoFHandler<1>&,
-                                              const std::vector<unsigned int>  &,
-                                              SparsityPattern    &)
+void
+DoFTools::make_boundary_sparsity_pattern (const DoFHandler<1> &dof_handler,
+                                         const DoFHandler<1>::FunctionMap  &function_map,
+                                         const std::vector<unsigned int>  &dof_to_boundary_mapping,
+                                         SparsityPattern    &sparsity)
 {
-  Assert (false, ExcInternalError());
+  const unsigned int dofs_per_vertex = dof_handler.get_fe().dofs_per_vertex;
+  std::vector<unsigned int> boundary_dof_boundary_indices (dofs_per_vertex);
+  
+                                  // first check left, the right
+                                  // boundary point
+  for (unsigned int direction=0; direction<2; ++direction)
+    {
+                                      // if this boundary is not
+                                      // requested, then go on with next one
+      if (function_map.find(direction) ==
+         function_map.end())
+       continue;
+
+                                      // find active cell at that
+                                      // boundary: first go to
+                                      // left/right, then to children
+      DoFHandler<1>::cell_iterator cell = dof_handler.begin(0);
+      while (!cell->at_boundary(direction))
+       cell = cell->neighbor(direction);
+      while (!cell->active())
+       cell = cell->child(direction);
+
+                                      // next get boundary mapped dof
+                                      // indices of boundary dofs
+      for (unsigned int i=0; i<dofs_per_vertex; ++i)
+       boundary_dof_boundary_indices[i]
+         = dof_to_boundary_mapping[cell->vertex_dof_index(direction,i)];
+
+      for (unsigned int i=0; i<dofs_per_vertex; ++i)
+       for (unsigned int j=0; j<dofs_per_vertex; ++j)
+         sparsity.add (boundary_dof_boundary_indices[i],
+                       boundary_dof_boundary_indices[j]);
+    };
 };
 
 
+
 template <>
-void
-DoFTools::make_boundary_sparsity_pattern (const DoFHandler<1>&,
-                                         const DoFHandler<1>::FunctionMap  &,
-                                         const std::vector<unsigned int>  &,
-                                         SparsityPattern    &)
+void DoFTools::make_boundary_sparsity_pattern (const DoFHandler<1> &dof_handler,
+                                              const std::vector<unsigned int>  &dof_to_boundary_mapping,
+                                              SparsityPattern    &sparsity)
 {
-  Assert (false, ExcInternalError());
-}
+                                  // there are only 2 boundary
+                                  // indicators in 1d, so it is no
+                                  // performance problem to call the
+                                  // other function
+  DoFHandler<1>::FunctionMap boundary_indicators;
+  boundary_indicators[0] = 0;
+  boundary_indicators[1] = 0;
+  make_boundary_sparsity_pattern (dof_handler, boundary_indicators,
+                                 dof_to_boundary_mapping, sparsity);
+};
+
 
 #endif
 
@@ -222,8 +264,8 @@ void DoFTools::make_boundary_sparsity_pattern (const DoFHandler<dim>& dof,
 
   const unsigned int dofs_per_face = dof.get_fe().dofs_per_face;
   std::vector<unsigned int> dofs_on_this_face(dofs_per_face);
-  DoFHandler<dim>::active_face_iterator face = dof.begin_active_face(),
-                                       endf = dof.end_face();
+  typename DoFHandler<dim>::active_face_iterator face = dof.begin_active_face(),
+                                                endf = dof.end_face();
   for (; face!=endf; ++face)
     if (boundary_indicators.find(face->boundary_indicator()) !=
        boundary_indicators.end())
@@ -1757,23 +1799,60 @@ DoFTools::compute_intergrid_weights_1 (const DoFHandler<dim>              &coars
 
 template <>
 void DoFTools::map_dof_to_boundary_indices (const DoFHandler<1> &dof_handler,
-                                           std::vector<unsigned int> &)
+                                           const std::set<unsigned char> &boundary_indicators,
+                                           std::vector<unsigned int> &mapping)
 {
   Assert (&dof_handler.get_fe() != 0, ExcNoFESelected());
-  Assert (false, ExcNotImplemented());
+
+  mapping.clear ();
+  mapping.insert (mapping.end(), dof_handler.n_dofs(),
+                 DoFHandler<1>::invalid_dof_index);
+
+  unsigned int next_free_index = 0;
+  
+                                  // first check left, the right
+                                  // boundary point
+  for (unsigned int direction=0; direction<2; ++direction)
+    {
+                                      // if this boundary is not
+                                      // requested, then go on with next one
+      if (boundary_indicators.find(direction) ==
+         boundary_indicators.end())
+       continue;
+
+                                      // find active cell at that
+                                      // boundary: first go to
+                                      // left/right, then to children
+      DoFHandler<1>::cell_iterator cell = dof_handler.begin(0);
+      while (!cell->at_boundary(direction))
+       cell = cell->neighbor(direction);
+      while (!cell->active())
+       cell = cell->child(direction);
+
+                                      // next enumerate these degrees
+                                      // of freedom
+      for (unsigned int i=0; i<dof_handler.get_fe().dofs_per_vertex; ++i)
+       mapping[cell->vertex_dof_index(direction,i)] = next_free_index++;
+    };
 };
 
 
 
 template <>
 void DoFTools::map_dof_to_boundary_indices (const DoFHandler<1> &dof_handler,
-                                           const std::set<unsigned char> &,
-                                           std::vector<unsigned int> &)
+                                           std::vector<unsigned int> &mapping)
 {
   Assert (&dof_handler.get_fe() != 0, ExcNoFESelected());
-  Assert (false, ExcNotImplemented());
-};
 
+                                  // in 1d, there are only 2 boundary
+                                  // indicators, so enumerate them
+                                  // and pass on to the other
+                                  // function
+  std::set<unsigned char> boundary_indicators;
+  boundary_indicators.insert (0U);
+  boundary_indicators.insert (1U);
+  map_dof_to_boundary_indices (dof_handler, boundary_indicators, mapping);
+};
 
 #else
 
index b21614cd81285382f92a683c88dd8524f6c886b9..2fc2165f46b9043423356a8f73fbe1d62ccc1013 100644 (file)
@@ -804,6 +804,24 @@ VectorTools::interpolate_boundary_values (const DoFHandler<dim>         &dof,
 }
 
 
+#if deal_II_dimension == 1
+
+template <>
+void
+VectorTools::project_boundary_values (const Mapping<1>       &mapping,
+                                     const DoFHandler<1>    &dof,
+                                     const FunctionMap<1>::type &boundary_functions,
+                                     const Quadrature<0>  &,
+                                     std::map<unsigned int,double> &boundary_values)
+{
+                                  // projection in 1d is equivalent
+                                  // to interpolation
+  interpolate_boundary_values (mapping, dof, boundary_functions,
+                              boundary_values, std::vector<bool>());
+};
+
+#endif
+
 
 template <int dim>
 void
@@ -1371,12 +1389,15 @@ void VectorTools::integrate_difference (const DoFHandler<deal_II_dimension> &,
                                        const Quadrature<deal_II_dimension> &,
                                        const NormType                      &,
                                        const Function<deal_II_dimension>   *);
+#if deal_II_dimension != 1
 template
 void VectorTools::project_boundary_values (const Mapping<deal_II_dimension>     &,
                                           const DoFHandler<deal_II_dimension>  &,
                                           const FunctionMap<deal_II_dimension>::FunctionMap &,
                                           const Quadrature<deal_II_dimension-1>&,
                                           std::map<unsigned int,double>        &);
+#endif
+
 template
 void VectorTools::project_boundary_values (const DoFHandler<deal_II_dimension>  &,
                                           const FunctionMap<deal_II_dimension>::FunctionMap &,
diff --git a/tests/deal.II/boundaries.cc b/tests/deal.II/boundaries.cc
new file mode 100644 (file)
index 0000000..8f488f4
--- /dev/null
@@ -0,0 +1,183 @@
+//----------------------------  boundaries.cc  ---------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2000, 2001 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.
+//
+//----------------------------  boundaries.cc  ---------------------------
+
+
+/* Author: Wolfgang Bangerth, University of Heidelberg, 2001 */
+/* Purpose: check interpolation and projection of boundary values. */
+
+
+
+#include <base/logstream.h>
+#include <base/function_lib.h>
+#include <base/quadrature_lib.h>
+#include <lac/vector.h>
+#include <grid/tria.h>
+#include <grid/tria_iterator.h>
+#include <grid/tria_boundary_lib.h>
+#include <dofs/dof_accessor.h>
+#include <grid/grid_generator.h>
+#include <dofs/dof_handler.h>
+#include <dofs/dof_tools.h>
+#include <fe/fe_q.h>
+#include <fe/fe_system.h>
+#include <fe/mapping_q.h>
+#include <numerics/vectors.h>
+
+#include <fstream>
+
+
+template<int dim>
+class MySquareFunction : public Function<dim>
+{
+  public:
+    MySquareFunction () : Function<dim>(2) {};
+    
+    virtual double value (const Point<dim>   &p,
+                         const unsigned int  component) const
+      {        return 100*(component+1)*p.square()*sin(p.square()); };
+    
+    virtual void   vector_value (const Point<dim>   &p,
+                                Vector<double>     &values) const
+      { for (unsigned int d=0; d<n_components; ++d) values(d) = value(p,d); };
+};
+
+
+template <int dim>
+const Quadrature<dim-1> &
+boundary_q (const DoFHandler<dim> &)
+{
+  static const QGauss4<dim-1> q;
+  return q;
+};
+
+
+const Quadrature<0> &
+boundary_q (const DoFHandler<1> &)
+{
+  return *static_cast<const Quadrature<0>*>(0);
+};
+
+
+void write_map (const std::map<unsigned int,double> &bv)
+{
+  for (std::map<unsigned int,double>::const_iterator
+        i=bv.begin(); i!=bv.end(); ++i)
+    deallog << i->first << ' ' << i->second <<std::endl;
+};
+
+      
+
+
+template <int dim>
+void
+check ()
+{
+  Triangulation<dim> tr;  
+  if (dim==2)
+    {
+      GridGenerator::hyper_ball(tr);
+    }
+  else
+    GridGenerator::hyper_cube(tr, -1./sqrt(dim),1./sqrt(dim));
+  if (dim != 1)
+    {
+      static const HyperBallBoundary<dim> boundary;
+      tr.set_boundary (0, boundary);
+    };
+  tr.refine_global (1);
+  tr.begin_active()->set_refine_flag ();
+  tr.execute_coarsening_and_refinement ();
+  if (dim==1)
+    tr.refine_global(2);
+
+                                  // use a cubic mapping to make
+                                  // things a little more complicated
+  MappingQ<dim> mapping(3);
+
+  
+                                  // list of finite elements for
+                                  // which we want check, and
+                                  // associated list of boundary
+                                  // value functions
+  std::vector<const FiniteElement<dim>*> fe_list;
+  std::vector<const Function<dim>*> function_list;
+
+                                  // FE1: a system of a quadratic and
+                                  // a linear element
+  fe_list.push_back (new FESystem<dim> (FE_Q<dim>(2), 1, FE_Q<dim>(1), 1));
+  function_list.push_back (new MySquareFunction<dim>());
+
+                                  // FE2: a linear element, to make
+                                  // things simple
+  fe_list.push_back (new FE_Q<dim> (1));
+  function_list.push_back (new SquareFunction<dim>());
+  
+                                  // check all of them
+  for (unsigned int i=0; i<fe_list.size(); ++i)
+    {
+      const FiniteElement<dim> &fe = *fe_list[i];
+      
+      DoFHandler<dim> dof(tr);
+      dof.distribute_dofs(fe);
+
+      typename FunctionMap<dim>::type function_map;
+      function_map[0] = function_list[i];
+
+                                      // interpolate boundary values
+      deallog << "Interpolated boundary values" << std::endl;
+      std::map<unsigned int,double> interpolated_bv;
+      VectorTools::interpolate_boundary_values (mapping, dof, function_map,
+                                               interpolated_bv);
+      write_map (interpolated_bv);
+
+                                      // project boundary values
+                                      // presently this is not
+                                      // implemented for 3d
+      if (dim != 3)
+       {
+         deallog << "Projected boundary values" << std::endl;
+         std::map<unsigned int,double> projected_bv;
+         VectorTools::project_boundary_values (mapping, dof, function_map,
+                                               boundary_q(dof), projected_bv);
+         write_map (projected_bv);
+       };
+    };
+  
+      
+                                  // delete objects now no more needed
+  for (unsigned int i=0; i<fe_list.size(); ++i)
+    {
+      delete fe_list[i];
+      delete function_list[i];
+    };
+}
+
+
+int main ()
+{
+  ofstream logfile ("boundaries.output");
+  logfile.precision (2);
+  logfile.setf(ios::fixed);  
+  deallog.attach(logfile);
+  deallog.depth_console (0);
+
+  deallog.push ("1d");
+  check<1> ();
+  deallog.pop ();
+  deallog.push ("2d");
+  check<2> ();
+  deallog.pop ();
+  deallog.push ("3d");
+  check<3> ();
+  deallog.pop ();
+}
diff --git a/tests/deal.II/boundaries.checked b/tests/deal.II/boundaries.checked
new file mode 100644 (file)
index 0000000..9fb2ca9
--- /dev/null
@@ -0,0 +1,315 @@
+
+DEAL:1d::Interpolated boundary 
+DEAL:1d::14 84.15
+DEAL:1d::15 168.29
+DEAL:1d::Projected boundary 
+DEAL:1d::14 84.15
+DEAL:1d::15 168.29
+DEAL:1d::Interpolated boundary 
+DEAL:1d::5 1.00
+DEAL:1d::Projected boundary 
+DEAL:1d::5 1.00
+DEAL:2d::Interpolated boundary 
+DEAL:2d::0 84.15
+DEAL:2d::1 168.29
+DEAL:2d::2 84.15
+DEAL:2d::3 168.29
+DEAL:2d::8 84.12
+DEAL:2d::29 84.15
+DEAL:2d::30 168.29
+DEAL:2d::33 84.15
+DEAL:2d::34 168.29
+DEAL:2d::38 84.12
+DEAL:2d::53 84.15
+DEAL:2d::54 168.29
+DEAL:2d::56 84.12
+DEAL:2d::78 84.15
+DEAL:2d::79 168.29
+DEAL:2d::82 84.12
+DEAL:2d::86 84.15
+DEAL:2d::87 168.29
+DEAL:2d::90 84.12
+DEAL:2d::100 84.15
+DEAL:2d::101 168.29
+DEAL:2d::104 84.12
+DEAL:2d::110 84.12
+DEAL:2d::112 84.15
+DEAL:2d::113 168.29
+DEAL:2d::118 84.15
+DEAL:2d::125 84.15
+DEAL:2d::Projected boundary 
+DEAL:2d:cg::Starting 
+DEAL:2d:cg::Convergence step 15 
+DEAL:2d::0 84.19
+DEAL:2d::1 168.32
+DEAL:2d::2 84.20
+DEAL:2d::3 168.33
+DEAL:2d::8 84.14
+DEAL:2d::29 84.19
+DEAL:2d::30 168.32
+DEAL:2d::33 84.20
+DEAL:2d::34 168.33
+DEAL:2d::38 84.14
+DEAL:2d::53 84.21
+DEAL:2d::54 168.32
+DEAL:2d::56 84.14
+DEAL:2d::78 84.21
+DEAL:2d::79 168.32
+DEAL:2d::82 84.14
+DEAL:2d::86 84.21
+DEAL:2d::87 168.33
+DEAL:2d::90 84.14
+DEAL:2d::100 84.21
+DEAL:2d::101 168.33
+DEAL:2d::104 84.14
+DEAL:2d::110 84.14
+DEAL:2d::112 84.16
+DEAL:2d::113 168.29
+DEAL:2d::118 84.14
+DEAL:2d::125 84.14
+DEAL:2d::Interpolated boundary 
+DEAL:2d::0 1.00
+DEAL:2d::1 1.00
+DEAL:2d::8 1.00
+DEAL:2d::10 1.00
+DEAL:2d::14 1.00
+DEAL:2d::19 1.00
+DEAL:2d::21 1.00
+DEAL:2d::24 1.00
+DEAL:2d::25 1.00
+DEAL:2d::Projected boundary 
+DEAL:2d:cg::Starting 
+DEAL:2d:cg::Convergence step 9 
+DEAL:2d::0 1.00
+DEAL:2d::1 1.00
+DEAL:2d::8 1.00
+DEAL:2d::10 1.00
+DEAL:2d::14 1.00
+DEAL:2d::19 1.00
+DEAL:2d::21 1.00
+DEAL:2d::24 1.00
+DEAL:2d::25 1.00
+DEAL:3d::Interpolated boundary 
+DEAL:3d::0 84.15
+DEAL:3d::1 168.29
+DEAL:3d::2 84.15
+DEAL:3d::3 168.29
+DEAL:3d::4 84.15
+DEAL:3d::5 168.29
+DEAL:3d::6 84.15
+DEAL:3d::7 168.29
+DEAL:3d::8 84.15
+DEAL:3d::9 168.29
+DEAL:3d::10 84.15
+DEAL:3d::11 168.29
+DEAL:3d::12 84.15
+DEAL:3d::13 168.29
+DEAL:3d::16 84.14
+DEAL:3d::17 84.14
+DEAL:3d::18 84.12
+DEAL:3d::19 84.12
+DEAL:3d::20 84.12
+DEAL:3d::21 84.12
+DEAL:3d::24 84.12
+DEAL:3d::25 84.14
+DEAL:3d::26 84.12
+DEAL:3d::28 84.08
+DEAL:3d::30 84.08
+DEAL:3d::31 84.08
+DEAL:3d::35 84.15
+DEAL:3d::36 168.29
+DEAL:3d::37 84.15
+DEAL:3d::38 168.29
+DEAL:3d::39 84.15
+DEAL:3d::40 168.29
+DEAL:3d::41 84.15
+DEAL:3d::42 168.29
+DEAL:3d::43 84.14
+DEAL:3d::44 84.14
+DEAL:3d::45 84.12
+DEAL:3d::46 84.12
+DEAL:3d::47 84.12
+DEAL:3d::49 84.14
+DEAL:3d::50 84.12
+DEAL:3d::51 84.08
+DEAL:3d::53 84.08
+DEAL:3d::54 84.08
+DEAL:3d::57 84.15
+DEAL:3d::58 168.29
+DEAL:3d::59 84.15
+DEAL:3d::60 168.29
+DEAL:3d::61 84.15
+DEAL:3d::62 168.29
+DEAL:3d::63 84.15
+DEAL:3d::64 168.29
+DEAL:3d::65 84.12
+DEAL:3d::66 84.14
+DEAL:3d::67 84.14
+DEAL:3d::69 84.12
+DEAL:3d::70 84.12
+DEAL:3d::71 84.12
+DEAL:3d::72 84.14
+DEAL:3d::73 84.08
+DEAL:3d::76 84.08
+DEAL:3d::77 84.08
+DEAL:3d::79 84.15
+DEAL:3d::80 168.29
+DEAL:3d::81 84.15
+DEAL:3d::82 168.29
+DEAL:3d::83 84.15
+DEAL:3d::84 168.29
+DEAL:3d::85 84.15
+DEAL:3d::86 168.29
+DEAL:3d::87 84.15
+DEAL:3d::88 168.29
+DEAL:3d::89 84.12
+DEAL:3d::90 84.12
+DEAL:3d::91 84.14
+DEAL:3d::92 84.12
+DEAL:3d::93 84.12
+DEAL:3d::94 84.14
+DEAL:3d::95 84.14
+DEAL:3d::96 84.12
+DEAL:3d::98 84.12
+DEAL:3d::100 84.08
+DEAL:3d::101 84.08
+DEAL:3d::104 84.08
+DEAL:3d::106 84.15
+DEAL:3d::107 168.29
+DEAL:3d::108 84.15
+DEAL:3d::109 168.29
+DEAL:3d::110 84.14
+DEAL:3d::111 84.14
+DEAL:3d::112 84.12
+DEAL:3d::113 84.14
+DEAL:3d::114 84.12
+DEAL:3d::115 84.08
+DEAL:3d::116 84.08
+DEAL:3d::117 84.08
+DEAL:3d::120 84.15
+DEAL:3d::121 168.29
+DEAL:3d::122 84.15
+DEAL:3d::123 168.29
+DEAL:3d::124 84.14
+DEAL:3d::125 84.14
+DEAL:3d::126 84.12
+DEAL:3d::127 84.14
+DEAL:3d::128 84.12
+DEAL:3d::129 84.08
+DEAL:3d::130 84.08
+DEAL:3d::131 84.08
+DEAL:3d::134 84.15
+DEAL:3d::135 168.29
+DEAL:3d::136 84.14
+DEAL:3d::137 84.14
+DEAL:3d::138 84.14
+DEAL:3d::139 84.08
+DEAL:3d::140 84.08
+DEAL:3d::141 84.08
+DEAL:3d::143 84.15
+DEAL:3d::144 168.29
+DEAL:3d::145 84.15
+DEAL:3d::146 168.29
+DEAL:3d::147 84.15
+DEAL:3d::148 168.29
+DEAL:3d::149 84.15
+DEAL:3d::150 168.29
+DEAL:3d::151 84.15
+DEAL:3d::152 168.29
+DEAL:3d::153 84.15
+DEAL:3d::154 168.29
+DEAL:3d::157 84.15
+DEAL:3d::158 168.29
+DEAL:3d::159 84.15
+DEAL:3d::160 84.15
+DEAL:3d::161 84.15
+DEAL:3d::162 84.15
+DEAL:3d::163 84.15
+DEAL:3d::166 84.15
+DEAL:3d::167 84.15
+DEAL:3d::168 84.15
+DEAL:3d::170 84.15
+DEAL:3d::171 84.14
+DEAL:3d::173 84.14
+DEAL:3d::176 84.14
+DEAL:3d::178 84.15
+DEAL:3d::179 168.29
+DEAL:3d::180 84.15
+DEAL:3d::181 168.29
+DEAL:3d::184 84.15
+DEAL:3d::185 84.15
+DEAL:3d::186 84.15
+DEAL:3d::187 84.15
+DEAL:3d::190 84.15
+DEAL:3d::192 84.14
+DEAL:3d::194 84.14
+DEAL:3d::198 84.15
+DEAL:3d::199 168.29
+DEAL:3d::204 84.15
+DEAL:3d::205 84.15
+DEAL:3d::206 84.15
+DEAL:3d::212 84.14
+DEAL:3d::218 84.15
+DEAL:3d::219 168.29
+DEAL:3d::220 84.15
+DEAL:3d::221 84.15
+DEAL:3d::223 84.15
+DEAL:3d::224 84.15
+DEAL:3d::225 84.14
+DEAL:3d::228 84.14
+DEAL:3d::230 84.15
+DEAL:3d::231 168.29
+DEAL:3d::234 84.15
+DEAL:3d::235 168.29
+DEAL:3d::236 84.15
+DEAL:3d::239 84.15
+DEAL:3d::240 84.15
+DEAL:3d::241 84.15
+DEAL:3d::243 84.15
+DEAL:3d::245 84.14
+DEAL:3d::248 84.14
+DEAL:3d::252 84.15
+DEAL:3d::255 84.15
+DEAL:3d::258 84.14
+DEAL:3d::275 84.15
+DEAL:3d::276 84.15
+DEAL:3d::279 84.14
+DEAL:3d::Interpolated boundary 
+DEAL:3d::0 1.00
+DEAL:3d::1 1.00
+DEAL:3d::2 1.00
+DEAL:3d::3 1.00
+DEAL:3d::4 1.00
+DEAL:3d::5 1.00
+DEAL:3d::6 1.00
+DEAL:3d::8 1.00
+DEAL:3d::9 1.00
+DEAL:3d::10 1.00
+DEAL:3d::11 1.00
+DEAL:3d::12 1.00
+DEAL:3d::13 1.00
+DEAL:3d::14 1.00
+DEAL:3d::15 1.00
+DEAL:3d::16 1.00
+DEAL:3d::17 1.00
+DEAL:3d::18 1.00
+DEAL:3d::19 1.00
+DEAL:3d::20 1.00
+DEAL:3d::21 1.00
+DEAL:3d::22 1.00
+DEAL:3d::23 1.00
+DEAL:3d::24 1.00
+DEAL:3d::25 1.00
+DEAL:3d::26 1.00
+DEAL:3d::27 1.00
+DEAL:3d::28 1.00
+DEAL:3d::29 1.00
+DEAL:3d::30 1.00
+DEAL:3d::31 1.00
+DEAL:3d::33 1.00
+DEAL:3d::34 1.00
+DEAL:3d::35 1.00
+DEAL:3d::37 1.00
+DEAL:3d::40 1.00
+DEAL:3d::41 1.00
+DEAL:3d::43 1.00

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.