]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Provide a second function that computes the support point locations even for distribu...
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 26 Apr 2012 16:11:48 +0000 (16:11 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 26 Apr 2012 16:11:48 +0000 (16:11 +0000)
git-svn-id: https://svn.dealii.org/trunk@25455 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/doc/news/changes.h
deal.II/include/deal.II/dofs/dof_tools.h
deal.II/source/dofs/dof_tools.cc
deal.II/source/dofs/dof_tools.inst.in
tests/mpi/map_dofs_to_support_points.cc [new file with mode: 0644]
tests/mpi/map_dofs_to_support_points/ncpu_10/cmp/generic [new file with mode: 0644]
tests/mpi/map_dofs_to_support_points/ncpu_4/cmp/generic [new file with mode: 0644]

index fbbccfd87ebfdabc4555507376414b40fa1bb624..eb9b381013bbc674194a3a6023e203f0860e3609 100644 (file)
@@ -175,6 +175,12 @@ enabled due to a missing include file in file
 <h3>Specific improvements</h3>
 
 <ol>
+<li> New: There is now a second DoFTools::map_dofs_to_support_points
+function that also works for parallel::distributed::Triangulation
+triangulations.
+<br>
+(Wolfgang Bangerth, 2012/04/26)
+
 <li> New: There is now a second DoFTools::extract_boundary_dofs
 function that also works for parallel::distributed::Triangulation
 triangulations.
index 27014eeb90648f38a7ae757eb5ffcef87a8f4c32..3db11e7135fd8be027c2e63f81eeff27943189f5 100644 (file)
@@ -2251,9 +2251,19 @@ namespace DoFTools
                                     * or the like. Otherwise, an
                                     * exception is thrown.
                                     *
-                                    * The given array must have a
+                                    * @pre The given array must have a
                                     * length of as many elements as
                                     * there are degrees of freedom.
+                                    *
+                                    * @note The precondition to this function
+                                    * that the output argument needs to have
+                                    * size equal to the total number of degrees
+                                    * of freedom makes this function
+                                    * unsuitable for the case that the given
+                                    * DoFHandler object derives from a
+                                    * parallel::distributed::Triangulation object.
+                                    * Consequently, this function will produce an
+                                    * error if called with such a DoFHandler.
                                     */
   template <int dim, int spacedim>
   void
@@ -2262,15 +2272,55 @@ namespace DoFTools
                               std::vector<Point<spacedim> >     &support_points);
 
             /**
-             * Same as above for the hp case.
+             * Same as the previous function but for the hp case.
              */
-
   template <int dim, int spacedim>
   void
   map_dofs_to_support_points (const dealii::hp::MappingCollection<dim,spacedim>   &mapping,
         const hp::DoFHandler<dim,spacedim>    &dof_handler,
         std::vector<Point<spacedim> > &support_points);
 
+  /**
+   * This function is a version of the above map_dofs_to_support_points
+   * function that doesn't simply return a vector of support_points with one
+   * entry for each global degree of freedom, but instead a map that
+   * maps from the DoFs index to its location. The point of this
+   * function is that it is also usable in cases where the DoFHandler
+   * is based on a parallel::distributed::Triangulation object. In such cases,
+   * each processor will not be able to determine the support point location
+   * of all DoFs, and worse no processor may be able to hold a vector that
+   * would contain the locations of all DoFs even if they were known. As
+   * a consequence, this function constructs a map from those DoFs for which
+   * we can know the locations (namely, those DoFs that are
+   * locally relevant (see @ref GlossLocallyRelevantDof "locally relevant DoFs")
+   * to their locations.
+   *
+   * For non-distributed triangulations, the map returned as @p support_points
+   * is of course dense, i.e., every DoF is to be found in it.
+   *
+   * @param mapping The mapping from the reference cell to the real cell on
+   *        which DoFs are defined.
+   * @param dof_handler The object that describes which DoF indices live on
+   *        which cell of the triangulation.
+   * @param support_points A map that for every locally relevant DoF index
+   *        contains the corresponding location in real space coordinates.
+   *        Previous content of this object is deleted in this function.
+   */
+  template <int dim, int spacedim>
+  void
+  map_dofs_to_support_points (const Mapping<dim,spacedim>       &mapping,
+                              const DoFHandler<dim,spacedim>    &dof_handler,
+                              std::map<unsigned int, Point<spacedim> >     &support_points);
+
+            /**
+             * Same as the previous function but for the hp case.
+             */
+  template <int dim, int spacedim>
+  void
+  map_dofs_to_support_points (const dealii::hp::MappingCollection<dim,spacedim>   &mapping,
+                              const hp::DoFHandler<dim,spacedim>    &dof_handler,
+                              std::map<unsigned int, Point<spacedim> > &support_points);
+
                                    /**
                                     * This is the opposite function
                                     * to the one above. It generates
index b5eb53fe81ed0bc46808492bfced7cc88fa31dd4..b6efd20146b87cea6012d84b8e93a22851d2f5f0 100644 (file)
@@ -5743,14 +5743,12 @@ namespace DoFTools
   {
     namespace
     {
-      template<class DH>
-      void  map_dofs_to_support_points_internal(
-          const hp::MappingCollection<DH::dimension, DH::space_dimension> & mapping,
-          const DH  &dof_handler,
-          std::vector<Point<DH::space_dimension> >  &support_points)
+      template <class DH>
+      void
+      map_dofs_to_support_points(const hp::MappingCollection<DH::dimension, DH::space_dimension> & mapping,
+                                 const DH  &dof_handler,
+                                 std::map<unsigned int,Point<DH::space_dimension> >  &support_points)
       {
-          AssertDimension(support_points.size(), dof_handler.n_dofs());
-
           const unsigned int dim = DH::dimension;
           const unsigned int spacedim = DH::space_dimension;
 
@@ -5762,7 +5760,7 @@ namespace DoFTools
               // check whether every fe in the collection
               // has support points
               Assert(fe_collection[fe_index].has_support_points(),
-                  typename FiniteElement<dim>::ExcFEHasNoSupportPoints());
+                     typename FiniteElement<dim>::ExcFEHasNoSupportPoints());
               q_coll_dummy.push_back(
                   Quadrature<dim> (
                       fe_collection[fe_index].get_unit_support_points()));
@@ -5781,25 +5779,49 @@ namespace DoFTools
           // rule have been set to invalid values
           // by the used constructor.
           hp::FEValues<dim, spacedim> hp_fe_values(mapping, fe_collection,
-              q_coll_dummy, update_quadrature_points);
+                                                   q_coll_dummy, update_quadrature_points);
           typename DH::active_cell_iterator cell =
               dof_handler.begin_active(), endc = dof_handler.end();
 
           std::vector<unsigned int> local_dof_indices;
           for (; cell != endc; ++cell)
-            {
-              hp_fe_values.reinit(cell);
-              const FEValues<dim, spacedim> &fe_values = hp_fe_values.get_present_fe_values();
+            // only work on locally relevant cells
+            if (cell->is_artificial() == false)
+              {
+                hp_fe_values.reinit(cell);
+                const FEValues<dim, spacedim> &fe_values = hp_fe_values.get_present_fe_values();
 
-              local_dof_indices.resize(cell->get_fe().dofs_per_cell);
-              cell->get_dof_indices(local_dof_indices);
+                local_dof_indices.resize(cell->get_fe().dofs_per_cell);
+                cell->get_dof_indices(local_dof_indices);
 
-              const std::vector<Point<spacedim> > & points =
-                  fe_values.get_quadrature_points();
-              for (unsigned int i = 0; i < cell->get_fe().dofs_per_cell; ++i)
-                support_points[local_dof_indices[i]] = points[i];
-            };
+                const std::vector<Point<spacedim> > & points =
+                    fe_values.get_quadrature_points();
+                for (unsigned int i = 0; i < cell->get_fe().dofs_per_cell; ++i)
+                  // insert the values into the map
+                  support_points[local_dof_indices[i]] = points[i];
+              }
       }
+
+
+      template <class DH>
+      void
+      map_dofs_to_support_points(const hp::MappingCollection<DH::dimension, DH::space_dimension> & mapping,
+                                 const DH  &dof_handler,
+                                 std::vector<Point<DH::space_dimension> >  &support_points)
+        {
+          // get the data in the form of the map as above
+          std::map<unsigned int,Point<DH::space_dimension> >  x_support_points;
+          map_dofs_to_support_points(mapping, dof_handler, x_support_points);
+
+          // now convert from the map to the linear vector. make sure every
+          // entry really appeared in the map
+          for (unsigned int i=0; i<dof_handler.n_dofs(); ++i)
+            {
+              Assert (x_support_points.find(i) != x_support_points.end(),
+                      ExcInternalError());
+              support_points[i] = x_support_points[i];
+            }
+        }
     }
   }
 
@@ -5810,32 +5832,76 @@ namespace DoFTools
       std::vector<Point<spacedim> > &support_points)
   {
     AssertDimension(support_points.size(), dof_handler.n_dofs());
+    Assert ((dynamic_cast<const parallel::distributed::Triangulation<dim,spacedim>*>
+             (&dof_handler.get_tria())
+             ==
+             0),
+            ExcMessage ("This function can not be used with distributed triangulations."
+                        "See the documentation for more information."));
 
     //Let the internal function do all the work,
     //just make sure that it gets a MappingCollection
     const hp::MappingCollection<dim, spacedim> mapping_collection(mapping);
 
-    internal::map_dofs_to_support_points_internal(
-        mapping_collection,
-        dof_handler, support_points);
+    internal::map_dofs_to_support_points (mapping_collection,
+                                          dof_handler,
+                                          support_points);
   }
 
 
   template<int dim, int spacedim>
   void
-  map_dofs_to_support_points(
-      const hp::MappingCollection<dim, spacedim> &mapping,
-      const hp::DoFHandler<dim, spacedim> &dof_handler,
-      std::vector<Point<spacedim> > &support_points)
+  map_dofs_to_support_points(const hp::MappingCollection<dim, spacedim> &mapping,
+                             const hp::DoFHandler<dim, spacedim> &dof_handler,
+                             std::vector<Point<spacedim> > &support_points)
   {
     AssertDimension(support_points.size(), dof_handler.n_dofs());
+    Assert ((dynamic_cast<const parallel::distributed::Triangulation<dim,spacedim>*>
+             (&dof_handler.get_tria())
+             ==
+             0),
+            ExcMessage ("This function can not be used with distributed triangulations."
+                        "See the documentation for more information."));
+
+    //Let the internal function do all the work,
+    //just make sure that it gets a MappingCollection
+    internal::map_dofs_to_support_points (mapping,
+                                          dof_handler,
+                                          support_points);
+  }
+
+
+  template <int dim, int spacedim>
+  void
+  map_dofs_to_support_points (const Mapping<dim,spacedim>       &mapping,
+                              const DoFHandler<dim,spacedim>    &dof_handler,
+                              std::map<unsigned int, Point<spacedim> > &support_points)
+  {
+    support_points.clear();
+
+    //Let the internal function do all the work,
+    //just make sure that it gets a MappingCollection
+    const hp::MappingCollection<dim, spacedim> mapping_collection(mapping);
+
+    internal::map_dofs_to_support_points (mapping_collection,
+                                          dof_handler,
+                                          support_points);
+  }
+
+
+  template<int dim, int spacedim>
+  void
+  map_dofs_to_support_points(const hp::MappingCollection<dim, spacedim> &mapping,
+                             const hp::DoFHandler<dim, spacedim> &dof_handler,
+                             std::map<unsigned int, Point<spacedim> > &support_points)
+  {
+    support_points.clear();
 
     //Let the internal function do all the work,
     //just make sure that it gets a MappingCollection
-    internal::map_dofs_to_support_points_internal(
-        mapping,
-        dof_handler,
-        support_points);
+    internal::map_dofs_to_support_points (mapping,
+                                          dof_handler,
+                                          support_points);
   }
 
 
index fe1259fe0389cf3b14f0238525c1d364fa462f3e..b4c14494df7b3b5151ffe533a868e87edc419dbd 100644 (file)
@@ -812,6 +812,22 @@ DoFTools::map_dofs_to_support_points<deal_II_dimension>
  const hp::DoFHandler<deal_II_dimension>&,
  std::vector<Point<deal_II_dimension> >&);
 
+
+template
+void
+DoFTools::map_dofs_to_support_points<deal_II_dimension>
+(const Mapping<deal_II_dimension,deal_II_dimension>&,
+ const DoFHandler<deal_II_dimension>&,
+ std::map<unsigned int, Point<deal_II_dimension> >&);
+
+
+template
+void
+DoFTools::map_dofs_to_support_points<deal_II_dimension>
+(const hp::MappingCollection<deal_II_dimension,deal_II_dimension>&,
+ const hp::DoFHandler<deal_II_dimension>&,
+ std::map<unsigned int, Point<deal_II_dimension> >&);
+
 #if deal_II_dimension < 3
 
 template
@@ -828,6 +844,21 @@ DoFTools::map_dofs_to_support_points<deal_II_dimension,deal_II_dimension+1>
  const hp::DoFHandler<deal_II_dimension, deal_II_dimension+1>&,
  std::vector<Point<deal_II_dimension+1> >&);
 
+template
+void
+DoFTools::map_dofs_to_support_points<deal_II_dimension,deal_II_dimension+1>
+(const Mapping<deal_II_dimension,deal_II_dimension+1>&,
+ const DoFHandler<deal_II_dimension, deal_II_dimension+1>&,
+ std::map<unsigned int, Point<deal_II_dimension+1> >&);
+
+template
+void
+DoFTools::map_dofs_to_support_points<deal_II_dimension,deal_II_dimension+1>
+(const hp::MappingCollection<deal_II_dimension,deal_II_dimension+1>&,
+ const hp::DoFHandler<deal_II_dimension, deal_II_dimension+1>&,
+ std::map<unsigned int, Point<deal_II_dimension+1> >&);
+
+
 template
 void
 DoFTools::count_dofs_per_block<DoFHandler<deal_II_dimension,deal_II_dimension+1> > (
diff --git a/tests/mpi/map_dofs_to_support_points.cc b/tests/mpi/map_dofs_to_support_points.cc
new file mode 100644 (file)
index 0000000..ea51ad6
--- /dev/null
@@ -0,0 +1,115 @@
+//---------------------------------------------------------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2009, 2010, 2011, 2012 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.
+//
+//---------------------------------------------------------------------------
+
+
+// Test DoFTools::map_dofs_to_support_points for parallel DoFHandlers
+
+#include "../tests.h"
+#include "coarse_grid_common.h"
+#include <deal.II/base/logstream.h>
+#include <deal.II/distributed/tria.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/dofs/dof_handler.h>
+#include <deal.II/dofs/dof_tools.h>
+#include <deal.II/fe/fe_q.h>
+
+#include <fstream>
+
+
+template<int dim>
+void test()
+{
+  parallel::distributed::Triangulation<dim> tr(MPI_COMM_WORLD);
+
+  GridGenerator::hyper_ball(tr);
+  tr.refine_global (1);
+
+  const FE_Q<dim> fe(1);
+  DoFHandler<dim> dofh(tr);
+  dofh.distribute_dofs (fe);
+
+  std::map<unsigned int, Point<dim> > points;
+  DoFTools::map_dofs_to_support_points (MappingQ1<dim>(),
+                                       dofh,
+                                       points);
+  if (Utilities::MPI::this_mpi_process (MPI_COMM_WORLD) == 0)
+    {
+      for (typename std::map<unsigned int, Point<dim> >::const_iterator
+            p = points.begin();
+          p != points.end();
+          ++p)
+       deallog << p->first << " -> " << p->second
+               << std::endl;
+    }
+
+                                  // the result of the call above is
+                                  // supposed to be a map that
+                                  // contains exactly the locally
+                                  // relevant dofs, so test this
+  IndexSet relevant_set;
+  DoFTools::extract_locally_relevant_dofs (dofh, relevant_set);
+
+  for (unsigned int i=0; i<dofh.n_dofs(); ++i)
+    if (relevant_set.is_element(i))
+      Assert (points.find(i) != points.end(),
+             ExcInternalError())
+      else
+      Assert (points.find(i) == points.end(),
+             ExcInternalError());
+}
+
+
+int main(int argc, char *argv[])
+{
+#ifdef DEAL_II_COMPILER_SUPPORTS_MPI
+  MPI_Init (&argc,&argv);
+#else
+  (void)argc;
+  (void)argv;
+#endif
+
+  unsigned int myid = Utilities::MPI::this_mpi_process (MPI_COMM_WORLD);
+
+
+  deallog.push(Utilities::int_to_string(myid));
+
+  if (myid == 0)
+    {
+      std::ofstream logfile(output_file_for_mpi("map_dofs_to_support_points").c_str());
+      deallog.attach(logfile);
+      deallog.depth_console(0);
+      deallog.threshold_double(1.e-10);
+
+      deallog.push("2d");
+      test<2>();
+      deallog.pop();
+
+      deallog.push("3d");
+      test<3>();
+      deallog.pop();
+    }
+  else
+    {
+      deallog.push("2d");
+      test<2>();
+      deallog.pop();
+
+      deallog.push("3d");
+      test<3>();
+      deallog.pop();
+    }
+
+#ifdef DEAL_II_COMPILER_SUPPORTS_MPI
+  MPI_Finalize();
+#endif
+}
diff --git a/tests/mpi/map_dofs_to_support_points/ncpu_10/cmp/generic b/tests/mpi/map_dofs_to_support_points/ncpu_10/cmp/generic
new file mode 100644 (file)
index 0000000..072424f
--- /dev/null
@@ -0,0 +1,69 @@
+
+DEAL:0:2d::0 -> -0.707107 -0.707107
+DEAL:0:2d::1 -> 0.00000 -0.707107
+DEAL:0:2d::2 -> -0.500000 -0.500000
+DEAL:0:2d::3 -> 0.00000 -0.500000
+DEAL:0:2d::4 -> 0.707107 -0.707107
+DEAL:0:2d::5 -> 0.500000 -0.500000
+DEAL:0:2d::6 -> -0.292893 -0.292893
+DEAL:0:2d::7 -> 0.00000 -0.292893
+DEAL:0:2d::8 -> 0.292893 -0.292893
+DEAL:0:2d::9 -> -0.707107 0.00000
+DEAL:0:2d::10 -> -0.500000 0.00000
+DEAL:0:2d::11 -> -0.292893 0.00000
+DEAL:0:2d::15 -> 0.707107 0.00000
+DEAL:0:2d::16 -> 0.500000 0.00000
+DEAL:0:2d::19 -> 0.292893 0.00000
+DEAL:0:2d::21 -> 0.00000 0.00000
+DEAL:0:3d::0 -> -0.577350 -0.577350 -0.577350
+DEAL:0:3d::1 -> 0.00000 -0.577350 -0.577350
+DEAL:0:3d::2 -> -0.577350 0.00000 -0.577350
+DEAL:0:3d::3 -> 0.00000 0.00000 -0.577350
+DEAL:0:3d::4 -> -0.394338 -0.394338 -0.394338
+DEAL:0:3d::5 -> 1.73472e-18 -0.394338 -0.394338
+DEAL:0:3d::6 -> -0.394338 0.00000 -0.394338
+DEAL:0:3d::7 -> 2.89121e-19 0.00000 -0.394338
+DEAL:0:3d::8 -> 0.577350 -0.577350 -0.577350
+DEAL:0:3d::9 -> 0.577350 0.00000 -0.577350
+DEAL:0:3d::10 -> 0.394338 -0.394338 -0.394338
+DEAL:0:3d::11 -> 0.394338 0.00000 -0.394338
+DEAL:0:3d::12 -> -0.577350 0.577350 -0.577350
+DEAL:0:3d::13 -> 0.00000 0.577350 -0.577350
+DEAL:0:3d::14 -> -0.394338 0.394338 -0.394338
+DEAL:0:3d::15 -> 1.73472e-18 0.394338 -0.394338
+DEAL:0:3d::16 -> 0.577350 0.577350 -0.577350
+DEAL:0:3d::17 -> 0.394338 0.394338 -0.394338
+DEAL:0:3d::18 -> -0.211325 -0.211325 -0.211325
+DEAL:0:3d::19 -> 0.00000 -0.211325 -0.211325
+DEAL:0:3d::20 -> -0.211325 0.00000 -0.211325
+DEAL:0:3d::21 -> 0.00000 0.00000 -0.211325
+DEAL:0:3d::22 -> 0.211325 -0.211325 -0.211325
+DEAL:0:3d::23 -> 0.211325 0.00000 -0.211325
+DEAL:0:3d::24 -> -0.211325 0.211325 -0.211325
+DEAL:0:3d::25 -> 0.00000 0.211325 -0.211325
+DEAL:0:3d::26 -> 0.211325 0.211325 -0.211325
+DEAL:0:3d::27 -> 0.577350 -0.577350 0.00000
+DEAL:0:3d::28 -> 0.577350 0.00000 0.00000
+DEAL:0:3d::29 -> 0.394338 -0.394338 1.73472e-18
+DEAL:0:3d::30 -> 0.394338 0.00000 6.93889e-18
+DEAL:0:3d::31 -> 0.577350 0.577350 0.00000
+DEAL:0:3d::32 -> 0.394338 0.394338 1.73472e-18
+DEAL:0:3d::33 -> 0.211325 -0.211325 0.00000
+DEAL:0:3d::34 -> 0.211325 0.00000 0.00000
+DEAL:0:3d::35 -> 0.211325 0.211325 0.00000
+DEAL:0:3d::45 -> -0.577350 -0.577350 0.00000
+DEAL:0:3d::46 -> -0.394338 -0.394338 0.00000
+DEAL:0:3d::47 -> -0.577350 0.00000 0.00000
+DEAL:0:3d::48 -> -0.394338 1.44560e-19 -6.93889e-18
+DEAL:0:3d::49 -> -0.211325 -0.211325 0.00000
+DEAL:0:3d::50 -> -0.211325 0.00000 0.00000
+DEAL:0:3d::51 -> -0.577350 0.577350 0.00000
+DEAL:0:3d::52 -> -0.394338 0.394338 0.00000
+DEAL:0:3d::53 -> -0.211325 0.211325 0.00000
+DEAL:0:3d::63 -> 0.00000 -0.577350 0.00000
+DEAL:0:3d::64 -> 1.44560e-19 -0.394338 6.93889e-18
+DEAL:0:3d::65 -> 0.00000 -0.211325 0.00000
+DEAL:0:3d::69 -> 0.00000 0.577350 0.00000
+DEAL:0:3d::70 -> 1.44560e-19 0.394338 0.00000
+DEAL:0:3d::71 -> 0.00000 0.211325 0.00000
+DEAL:0:3d::75 -> 0.00000 0.00000 0.00000
diff --git a/tests/mpi/map_dofs_to_support_points/ncpu_4/cmp/generic b/tests/mpi/map_dofs_to_support_points/ncpu_4/cmp/generic
new file mode 100644 (file)
index 0000000..bde99de
--- /dev/null
@@ -0,0 +1,87 @@
+
+DEAL:0:2d::0 -> -0.707107 -0.707107
+DEAL:0:2d::1 -> 0.00000 -0.707107
+DEAL:0:2d::2 -> -0.500000 -0.500000
+DEAL:0:2d::3 -> 0.00000 -0.500000
+DEAL:0:2d::4 -> 0.707107 -0.707107
+DEAL:0:2d::5 -> 0.500000 -0.500000
+DEAL:0:2d::6 -> -0.292893 -0.292893
+DEAL:0:2d::7 -> 0.00000 -0.292893
+DEAL:0:2d::8 -> 0.292893 -0.292893
+DEAL:0:2d::9 -> -0.707107 0.00000
+DEAL:0:2d::10 -> -0.500000 0.00000
+DEAL:0:2d::11 -> -0.292893 0.00000
+DEAL:0:2d::15 -> 0.707107 0.00000
+DEAL:0:2d::16 -> 0.500000 0.00000
+DEAL:0:2d::19 -> 0.292893 0.00000
+DEAL:0:2d::21 -> 0.00000 0.00000
+DEAL:0:3d::0 -> -0.577350 -0.577350 -0.577350
+DEAL:0:3d::1 -> 0.00000 -0.577350 -0.577350
+DEAL:0:3d::2 -> -0.577350 0.00000 -0.577350
+DEAL:0:3d::3 -> 0.00000 0.00000 -0.577350
+DEAL:0:3d::4 -> -0.394338 -0.394338 -0.394338
+DEAL:0:3d::5 -> 1.73472e-18 -0.394338 -0.394338
+DEAL:0:3d::6 -> -0.394338 0.00000 -0.394338
+DEAL:0:3d::7 -> 2.89121e-19 0.00000 -0.394338
+DEAL:0:3d::8 -> 0.577350 -0.577350 -0.577350
+DEAL:0:3d::9 -> 0.577350 0.00000 -0.577350
+DEAL:0:3d::10 -> 0.394338 -0.394338 -0.394338
+DEAL:0:3d::11 -> 0.394338 0.00000 -0.394338
+DEAL:0:3d::12 -> -0.577350 0.577350 -0.577350
+DEAL:0:3d::13 -> 0.00000 0.577350 -0.577350
+DEAL:0:3d::14 -> -0.394338 0.394338 -0.394338
+DEAL:0:3d::15 -> 1.73472e-18 0.394338 -0.394338
+DEAL:0:3d::16 -> 0.577350 0.577350 -0.577350
+DEAL:0:3d::17 -> 0.394338 0.394338 -0.394338
+DEAL:0:3d::18 -> -0.211325 -0.211325 -0.211325
+DEAL:0:3d::19 -> 0.00000 -0.211325 -0.211325
+DEAL:0:3d::20 -> -0.211325 0.00000 -0.211325
+DEAL:0:3d::21 -> 0.00000 0.00000 -0.211325
+DEAL:0:3d::22 -> 0.211325 -0.211325 -0.211325
+DEAL:0:3d::23 -> 0.211325 0.00000 -0.211325
+DEAL:0:3d::24 -> -0.211325 0.211325 -0.211325
+DEAL:0:3d::25 -> 0.00000 0.211325 -0.211325
+DEAL:0:3d::26 -> 0.211325 0.211325 -0.211325
+DEAL:0:3d::27 -> 0.577350 -0.577350 0.00000
+DEAL:0:3d::28 -> 0.577350 0.00000 0.00000
+DEAL:0:3d::29 -> 0.394338 -0.394338 1.73472e-18
+DEAL:0:3d::30 -> 0.394338 0.00000 6.93889e-18
+DEAL:0:3d::31 -> 0.577350 0.577350 0.00000
+DEAL:0:3d::32 -> 0.394338 0.394338 1.73472e-18
+DEAL:0:3d::33 -> 0.211325 -0.211325 0.00000
+DEAL:0:3d::34 -> 0.211325 0.00000 0.00000
+DEAL:0:3d::35 -> 0.211325 0.211325 0.00000
+DEAL:0:3d::36 -> 0.577350 -0.577350 0.577350
+DEAL:0:3d::37 -> 0.577350 0.00000 0.577350
+DEAL:0:3d::38 -> 0.394338 -0.394338 0.394338
+DEAL:0:3d::39 -> 0.394338 0.00000 0.394338
+DEAL:0:3d::40 -> 0.577350 0.577350 0.577350
+DEAL:0:3d::41 -> 0.394338 0.394338 0.394338
+DEAL:0:3d::42 -> 0.211325 -0.211325 0.211325
+DEAL:0:3d::43 -> 0.211325 0.00000 0.211325
+DEAL:0:3d::44 -> 0.211325 0.211325 0.211325
+DEAL:0:3d::45 -> -0.577350 -0.577350 0.00000
+DEAL:0:3d::46 -> -0.394338 -0.394338 0.00000
+DEAL:0:3d::47 -> -0.577350 0.00000 0.00000
+DEAL:0:3d::48 -> -0.394338 1.44560e-19 -6.93889e-18
+DEAL:0:3d::49 -> -0.211325 -0.211325 0.00000
+DEAL:0:3d::50 -> -0.211325 0.00000 0.00000
+DEAL:0:3d::51 -> -0.577350 0.577350 0.00000
+DEAL:0:3d::52 -> -0.394338 0.394338 0.00000
+DEAL:0:3d::53 -> -0.211325 0.211325 0.00000
+DEAL:0:3d::63 -> 0.00000 -0.577350 0.00000
+DEAL:0:3d::64 -> 1.44560e-19 -0.394338 6.93889e-18
+DEAL:0:3d::65 -> 0.00000 -0.211325 0.00000
+DEAL:0:3d::66 -> 0.00000 -0.577350 0.577350
+DEAL:0:3d::67 -> 0.00000 -0.394338 0.394338
+DEAL:0:3d::68 -> 0.00000 -0.211325 0.211325
+DEAL:0:3d::69 -> 0.00000 0.577350 0.00000
+DEAL:0:3d::70 -> 1.44560e-19 0.394338 0.00000
+DEAL:0:3d::71 -> 0.00000 0.211325 0.00000
+DEAL:0:3d::72 -> 0.00000 0.577350 0.577350
+DEAL:0:3d::73 -> 0.00000 0.394338 0.394338
+DEAL:0:3d::74 -> 0.00000 0.211325 0.211325
+DEAL:0:3d::75 -> 0.00000 0.00000 0.00000
+DEAL:0:3d::76 -> 0.00000 0.00000 0.211325
+DEAL:0:3d::77 -> 0.00000 0.00000 0.577350
+DEAL:0:3d::78 -> 0.00000 6.93889e-18 0.394338

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.