]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
prepared the MGTransferSelect for adaptive refinement, cleaned up in mg_tools.h
authorjanssen <janssen@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 5 Feb 2010 13:50:34 +0000 (13:50 +0000)
committerjanssen <janssen@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 5 Feb 2010 13:50:34 +0000 (13:50 +0000)
git-svn-id: https://svn.dealii.org/trunk@20505 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/multigrid/mg_tools.h
deal.II/deal.II/include/multigrid/mg_transfer.h
deal.II/deal.II/include/multigrid/mg_transfer.templates.h
deal.II/deal.II/include/multigrid/mg_transfer_component.h
deal.II/deal.II/include/multigrid/mg_transfer_component.templates.h
deal.II/deal.II/include/multigrid/multigrid.templates.h

index d421dffcdfd6128d192e92ccce47c60ff375771f..6967a36a6d9208b40bae2242442be9fb8db635c6 100644 (file)
@@ -243,6 +243,7 @@ class MGTools
       std::vector<std::set<unsigned int> >& boundary_indices,
       const std::vector<bool>& component_mask = std::vector<bool>());
                                       /**
+                                       * Maybe no longer needed.
                                        */
 
     template <typename number>
@@ -282,6 +283,16 @@ class MGTools
     extract_inner_interface_dofs (const MGDoFHandler<dim,spacedim> &mg_dof_handler,
                                  std::vector<std::vector<bool> >  &interface_dofs,
                                  std::vector<std::vector<bool> >  &boundary_interface_dofs);
+
+                                    /**
+                                      * Does the same as the function above, 
+                                      * but fills only the interface_dofs.
+                                     */
+    template <int dim, int spacedim>
+    static
+    void
+    extract_inner_interface_dofs (const MGDoFHandler<dim,spacedim> &mg_dof_handler,
+                                 std::vector<std::vector<bool> >  &interface_dofs);
 };
 
 /*@}*/
index 670ed5f0b3bccbc4b47a1fd7315c385fe8097ff4..0a7e9b618c49425f0d6c79cf957917a687136b2a 100644 (file)
@@ -86,7 +86,9 @@ class MGTransferPrebuilt : public MGTransferBase<VECTOR>
                                      * matrices for each level.
                                      */
     template <int dim, int spacedim>
-    void build_matrices (const MGDoFHandler<dim,spacedim> &mg_dof);
+    void build_matrices (const MGDoFHandler<dim,spacedim> &mg_dof,
+    const std::vector<std::set<unsigned int> >&boundary_indices
+        );
 
     virtual void prolongate (const unsigned int    to_level,
                             VECTOR       &dst,
@@ -220,7 +222,7 @@ class MGTransferPrebuilt : public MGTransferBase<VECTOR>
                                      * The data is first the global
                                      * index, then the level index.
                                     */
-    std::vector<std::vector<std::pair<unsigned int, unsigned int> > >
+    std::vector<std::map<unsigned int, unsigned int> >
       copy_indices;
 
                                     /**
@@ -231,36 +233,12 @@ class MGTransferPrebuilt : public MGTransferBase<VECTOR>
                                      */
     std::vector<unsigned int> component_to_block_map;
 
-                                    /**
-                                     * Fill the two boolean vectors
-                                     * #dofs_on_refinement_edge and
-                                     * #dofs_on_refinement_boundary.
-                                     */
-    template <int dim, int spacedim>
-    void find_dofs_on_refinement_edges (
-        const MGDoFHandler<dim,spacedim>& mg_dof);
-
                                     /**
                                      * Degrees of freedom on the
                                      * refinement edge excluding
                                      * those on the boundary.
-                                     *
-                                     * @todo Clean up names
-                                     */
-    std::vector<std::vector<bool> > dofs_on_refinement_edge;
-                                    /**
-                                     * The degrees of freedom on the
-                                     * the refinement edges. For each
-                                     * level (outer vector) and each
-                                     * dof index (inner vector), this
-                                     * bool is true if the level
-                                     * degree of freedom is on the
-                                     * refinement edge towards the
-                                     * lower level.
-                                     *
-                                     * @todo Clean up names
                                      */
-    std::vector<std::vector<bool> > dofs_on_refinement_boundary;
+    std::vector<std::vector<bool> > interface_dofs;
                                     /**
                                      * The constraints of the global
                                      * system.
index 96cdc1c4f74ab3933850f1fb1af26d59bf66a532..ac3e9308d8a284629b6dc208a5312e6b2ae2f605 100644 (file)
@@ -28,7 +28,6 @@ DEAL_II_NAMESPACE_OPEN
 
 /* --------------------- MGTransferPrebuilt -------------- */
 
-typedef std::vector<std::pair<unsigned int, unsigned int> >::const_iterator IT;
 
 
 
@@ -50,6 +49,7 @@ MGTransferPrebuilt<VECTOR>::copy_from_mg(
   dst = 0;
   for (unsigned int level=0;level<mg_dof_handler.get_tria().n_levels();++level)
   {
+    typedef std::map<unsigned int, unsigned int>::const_iterator IT;
     for (IT i= copy_indices[level].begin();
         i != copy_indices[level].end();++i)
       dst(i->first) = src[level](i->second);
@@ -76,9 +76,12 @@ MGTransferPrebuilt<VECTOR>::copy_from_mg_add (
                                       // have fine level basis
                                       // functions
   for (unsigned int level=0;level<mg_dof_handler.get_tria().n_levels();++level)
+  {
+    typedef std::map<unsigned int, unsigned int>::const_iterator IT;
     for (IT i= copy_indices[level].begin();
         i != copy_indices[level].end();++i)
       dst(i->first) += src[level](i->second);
+  }
 }
 
 
@@ -91,70 +94,6 @@ set_component_to_block_map (const std::vector<unsigned int> &map)
   component_to_block_map = map;
 }
 
-
-
-template <class VECTOR>
-template <int dim, int spacedim>
-void MGTransferPrebuilt<VECTOR>::find_dofs_on_refinement_edges (
-    const MGDoFHandler<dim,spacedim>& mg_dof)
-{
-  for (unsigned int level = 0; level<mg_dof.get_tria().n_levels(); ++level)
-  {
-    std::vector<bool> tmp (mg_dof.n_dofs(level));
-    dofs_on_refinement_edge.push_back (tmp);
-    dofs_on_refinement_boundary.push_back (tmp);
-  }
-
-  const unsigned int   dofs_per_cell   = mg_dof.get_fe().dofs_per_cell;
-  const unsigned int   dofs_per_face   = mg_dof.get_fe().dofs_per_face;
-
-  std::vector<unsigned int> local_dof_indices (dofs_per_cell);
-  std::vector<unsigned int> face_dof_indices (dofs_per_face);
-
-  typename MGDoFHandler<dim>::cell_iterator cell = mg_dof.begin(),
-           endc = mg_dof.end();
-
-  for (; cell!=endc; ++cell)
-  {
-    std::vector<bool> cell_dofs(dofs_per_cell);
-    std::vector<bool> boundary_cell_dofs(dofs_per_cell);
-    const unsigned int level = cell->level();
-    cell->get_mg_dof_indices (local_dof_indices);
-    for (unsigned int face_nr=0;
-        face_nr<GeometryInfo<dim>::faces_per_cell; ++face_nr)
-    {
-      typename DoFHandler<dim>::face_iterator face = cell->face(face_nr);
-      if(!cell->at_boundary(face_nr))
-      {
-        //interior face
-        typename MGDoFHandler<dim>::cell_iterator neighbor
-          = cell->neighbor(face_nr);
-        // Do refinement face
-        // from the coarse side
-        if (neighbor->level() < cell->level())
-        {
-          for(unsigned int j=0; j<dofs_per_face; ++j)
-            cell_dofs[mg_dof.get_fe().face_to_cell_index(j,face_nr)] = true;
-        }
-      }
-      //boundary face
-      else
-        for(unsigned int j=0; j<dofs_per_face; ++j)
-          boundary_cell_dofs[mg_dof.get_fe().face_to_cell_index(j,face_nr)] = true;
-    }//faces
-    for(unsigned int i=0; i<dofs_per_cell; ++i)
-    {
-      if(cell_dofs[i])
-      {
-        dofs_on_refinement_edge[level][local_dof_indices[i]] = true;
-        dofs_on_refinement_boundary[level][local_dof_indices[i]] = true;
-      }
-      if(boundary_cell_dofs[i])
-        dofs_on_refinement_boundary[level][local_dof_indices[i]] = true;
-    }
-  }
-}
-
 template <class VECTOR>
 unsigned int
 MGTransferPrebuilt<VECTOR>::memory_consumption () const
index f1d24840d1ce5dcbfea557bf2f1f210f342947a9..ea98f0b9577b98ed4f4cab46cff0ff19b3e8ef60 100644 (file)
@@ -16,6 +16,7 @@
 #include <base/config.h>
 
 #include <lac/block_vector.h>
+#include <lac/constraint_matrix.h>
 #ifdef DEAL_PREFER_MATRIX_EZ
 #  include <lac/sparse_matrix_ez.h>
 #  include <lac/block_sparse_matrix_ez.h>
@@ -135,13 +136,13 @@ class MGTransferComponentBase
                                    * Start index of each component.
                                    */
     std::vector<unsigned int> component_start;
-  
+
                                   /**
                                    * Start index of each component on
                                    * all levels.
                                    */
     std::vector<std::vector<unsigned int> > mg_component_start;
-    
+
                                     /**
                                      * Call build_matrices()
                                      * function first.
@@ -164,12 +165,21 @@ class MGTransferComponentBase
     std::vector<std_cxx1x::shared_ptr<BlockSparseMatrix<double> > > prolongation_matrices;
     
                                     /**
-                                     * Unused now, but intended to
-                                     * hold the mapping for the
+                                     * Holds the mapping for the
                                      * <tt>copy_to/from_mg</tt>-functions.
+                                     * The data is first the global
+                                     * index, then the level index.
                                      */
     std::vector<std::map<unsigned int, unsigned int> >
     copy_to_and_from_indices;
+
+                                    /**
+                                      * Store the boundary_indices. 
+                                      * These are needed for the 
+                                      * boundary values in the 
+                                      * restriction matrix.
+                                     */
+    std::vector<std::set<unsigned int> > boundary_indices;
 };
 
 //TODO:[GK] Update documentation for copy_* functions
@@ -191,6 +201,20 @@ class MGTransferSelect : public MGTransferBase<Vector<number> >,
                         private MGTransferComponentBase
 {
   public:
+                                    /**
+                                     * Constructor without constraint
+                                     * matrices. Use this constructor
+                                     * only with discontinuous finite
+                                     * elements or with no local
+                                     * refinement.
+                                     */
+    MGTransferSelect ();
+
+                                    /**
+                                     * Constructor with constraint matrices.
+                                     */
+    MGTransferSelect (const ConstraintMatrix& constraints);
+
                                     /**
                                      * Destructor.
                                      */
@@ -233,6 +257,9 @@ class MGTransferSelect : public MGTransferBase<Vector<number> >,
                                      * DoFRenumbering::component_wise). It
                                      * also affects the behavior of
                                      * the <tt>selected</tt> argument
+                                      *
+                                      * @arg boundary_indices: holds the
+                                      * boundary indices on each level.
                                      */
     template <int dim, int spacedim>
     void build_matrices (const DoFHandler<dim,spacedim> &dof,
@@ -242,7 +269,10 @@ class MGTransferSelect : public MGTransferBase<Vector<number> >,
                         const std::vector<unsigned int>& target_component
                         = std::vector<unsigned int>(),
                         const std::vector<unsigned int>& mg_target_component
-                        = std::vector<unsigned int>());
+                        = std::vector<unsigned int>(),
+                        const std::vector<std::set<unsigned int> >& boundary_indices
+                        = std::vector<std::set<unsigned int> >()
+                         );
 
                                     /**
                                      * Change selected
@@ -352,8 +382,7 @@ class MGTransferSelect : public MGTransferBase<Vector<number> >,
     void
     do_copy_from_mg (const MGDoFHandler<dim,spacedim>              &mg_dof,
                     OutVector                            &dst,
-                    const MGLevelObject<Vector<number> > &src,
-                    const unsigned int offset) const;
+                    const MGLevelObject<Vector<number> > &src) const;
 
                                     /**
                                      * Implementation of the public
@@ -363,8 +392,7 @@ class MGTransferSelect : public MGTransferBase<Vector<number> >,
     void
     do_copy_from_mg_add (const MGDoFHandler<dim,spacedim>              &mg_dof,
                         OutVector                            &dst,
-                        const MGLevelObject<Vector<number> > &src,
-                        const unsigned int offset) const;
+                        const MGLevelObject<Vector<number> > &src) const;
 
                                     /**
                                      * Actual implementation of
@@ -374,8 +402,7 @@ class MGTransferSelect : public MGTransferBase<Vector<number> >,
     void
     do_copy_to_mg (const MGDoFHandler<dim,spacedim>        &mg_dof,
                   MGLevelObject<Vector<number> > &dst,
-                  const InVector                 &src,
-                  const unsigned int              offset) const;
+                  const InVector                 &src) const;
                                      /**
                                       * Selected component of global vector.
                                       */
@@ -384,6 +411,25 @@ class MGTransferSelect : public MGTransferBase<Vector<number> >,
                                       * Selected component inside multigrid.
                                       */
     unsigned int mg_selected_component;
+
+                                    /**
+                                     * The degrees of freedom on the
+                                     * the refinement edges. For each
+                                     * level (outer vector) and each
+                                     * dof index (inner vector), this
+                                     * bool is true if the level
+                                     * degree of freedom is on the
+                                     * refinement edge towards the
+                                     * lower level excluding boundary dofs.
+                                     */
+    std::vector<std::vector<bool> > interface_dofs;
+
+                                    /**
+                                     * The constraints of the global
+                                     * system.
+                                     */
+  public:
+    SmartPointer<const ConstraintMatrix> constraints;
 };
 
 /*@}*/
index 4765f323c9df1d977d9bc0a963ef363d502a4c49..301c987de5b8cb58472583ef4d49b4e1484eb6eb 100644 (file)
 #include <multigrid/mg_dof_accessor.h>
 #include <multigrid/mg_tools.h>
 #include <multigrid/mg_transfer_component.h>
+#include <numerics/data_out.h>
 
 #include <algorithm>
+#include <sstream>
+#include <fstream>
 
 DEAL_II_NAMESPACE_OPEN
 
@@ -39,7 +42,7 @@ MGTransferSelect<number>::copy_to_mg (
   MGLevelObject<Vector<number> > &dst,
   const BlockVector<number2>     &src) const
 {
-  do_copy_to_mg (mg_dof_handler, dst, src, 0);
+  do_copy_to_mg (mg_dof_handler, dst, src);
 }
 
 
@@ -52,7 +55,7 @@ MGTransferSelect<number>::copy_to_mg (
   MGLevelObject<Vector<number> > &dst,
   const Vector<number2>          &src) const
 {
-  do_copy_to_mg (mg_dof_handler, dst, src, component_start[selected_component]);
+  do_copy_to_mg (mg_dof_handler, dst, src);
 }
 
 
@@ -65,7 +68,9 @@ MGTransferSelect<number>::copy_from_mg (
   BlockVector<number2>&                 dst,
   const MGLevelObject<Vector<number> >& src) const
 {
-  do_copy_from_mg (mg_dof_handler, dst, src, 0);
+  do_copy_from_mg (mg_dof_handler, dst, src);
+  if (constraints != 0)
+    constraints->condense(dst);
 }
 
 
@@ -78,8 +83,30 @@ MGTransferSelect<number>::copy_from_mg (
   Vector<number2>&                      dst,
   const MGLevelObject<Vector<number> >& src) const
 {
-  do_copy_from_mg (mg_dof_handler, dst, src,
-                  component_start[selected_component]);
+  do_copy_from_mg (mg_dof_handler, dst, src);
+  if (constraints != 0)
+  {
+    //If we were given constraints 
+    //apply them to the dst that goes 
+    //back now to the linear solver.
+    //Since constraints are globally
+    //defined create a global vector here
+    //and copy dst to the right component, 
+    //apply the constraints then and copy 
+    //the block back to dst.
+    const unsigned int n_blocks = 
+      *std::max_element(target_component.begin(), target_component.end()) + 1;
+    std::vector<unsigned int> dofs_per_block (n_blocks);
+    DoFTools::count_dofs_per_block (mg_dof_handler, dofs_per_block, target_component);
+    BlockVector<number> tmp;
+    tmp.reinit(n_blocks);
+    for(unsigned int b=0; b<n_blocks; ++b)
+      tmp.block(b).reinit(dofs_per_block[b]);
+    tmp.collect_sizes ();
+    tmp.block(target_component[selected_component]) = dst;
+    constraints->condense(tmp);
+    dst = tmp.block(target_component[selected_component]);
+  }
 }
 
 
@@ -92,7 +119,7 @@ MGTransferSelect<number>::copy_from_mg_add (
   BlockVector<number2>&                 dst,
   const MGLevelObject<Vector<number> >& src) const
 {
-  do_copy_from_mg_add (mg_dof_handler, dst, src, 0);
+  do_copy_from_mg_add (mg_dof_handler, dst, src);
 }
 
 
@@ -105,8 +132,7 @@ MGTransferSelect<number>::copy_from_mg_add (
   Vector<number2>&                      dst,
   const MGLevelObject<Vector<number> >& src) const
 {
-  do_copy_from_mg_add (mg_dof_handler, dst, src,
-                      component_start[selected_component]);
+  do_copy_from_mg_add (mg_dof_handler, dst, src);
 }
 
 
@@ -117,14 +143,13 @@ void
 MGTransferSelect<number>::do_copy_from_mg (
   const MGDoFHandler<dim,spacedim>              &mg_dof_handler,
   OutVector                            &dst,
-  const MGLevelObject<Vector<number> > &src,
-  const unsigned int offset) const
+  const MGLevelObject<Vector<number> > &src) const
 {
-  const FiniteElement<dim>& fe = mg_dof_handler.get_fe();
+//  const FiniteElement<dim>& fe = mg_dof_handler.get_fe();
 
-  const unsigned int dofs_per_cell = fe.dofs_per_cell;
-  std::vector<unsigned int> global_dof_indices (dofs_per_cell);
-  std::vector<unsigned int> level_dof_indices (dofs_per_cell);
+//  const unsigned int dofs_per_cell = fe.dofs_per_cell;
+//  std::vector<unsigned int> global_dof_indices (dofs_per_cell);
+//  std::vector<unsigned int> level_dof_indices (dofs_per_cell);
 
   typename MGDoFHandler<dim,spacedim>::active_cell_iterator
     level_cell = mg_dof_handler.begin_active();
@@ -137,44 +162,25 @@ MGTransferSelect<number>::do_copy_from_mg (
 
                                   // Note that the level is
                                   // monotonuosly increasing
+  dst = 0;
   for (; level_cell != endc; ++level_cell)
-    {
-      const unsigned int level = level_cell->level();
-
-                                      // get the dof numbers of
-                                      // this cell for the global
-                                      // and the level-wise
-                                      // numbering
-      level_cell->get_dof_indices (global_dof_indices);
-      level_cell->get_mg_dof_indices(level_dof_indices);
-
-                                      // copy level-wise data to
-                                      // global vector
-      for (unsigned int i=0; i<dofs_per_cell; ++i)
-       {
-         const unsigned int component
-           = mg_target_component[fe.system_to_component_index(i).first];
-       if (mg_selected[component])
-         {
-           const unsigned int level_start
-             = mg_component_start[level][component];
-           dst(global_dof_indices[i] - offset)
-             = src[level](level_dof_indices[i]-level_start);
-         }
-       }
-    }
+  {
+    const unsigned int level = level_cell->level();
+    typedef std::map<unsigned int, unsigned int>::const_iterator IT;
+    for (IT i=copy_to_and_from_indices[level].begin();
+        i != copy_to_and_from_indices[level].end(); ++i)
+      dst(i->first) = src[level](i->second);
+  }
 }
 
 
-
 template <typename number>
 template <int dim, class OutVector, int spacedim>
 void
 MGTransferSelect<number>::do_copy_from_mg_add (
   const MGDoFHandler<dim,spacedim>              &mg_dof_handler,
   OutVector                            &dst,
-  const MGLevelObject<Vector<number> > &src,
-  const unsigned int offset) const
+  const MGLevelObject<Vector<number> > &src) const
 {
   const FiniteElement<dim>& fe = mg_dof_handler.get_fe();
   const unsigned int dofs_per_cell = fe.dofs_per_cell;
@@ -193,30 +199,14 @@ MGTransferSelect<number>::do_copy_from_mg_add (
 
                                   // Note that the level is
                                   // monotonuosly increasing
+  dst = 0;
   for (; level_cell != endc; ++level_cell)
     {
       const unsigned int level = level_cell->level();
-
-                                      // get the dof numbers of
-                                      // this cell for the global
-                                      // and the level-wise
-                                      // numbering
-      level_cell->get_dof_indices (global_dof_indices);
-      level_cell->get_mg_dof_indices(level_dof_indices);
-                                      // copy level-wise data to
-                                      // global vector
-      for (unsigned int i=0; i<dofs_per_cell; ++i)
-       {
-         const unsigned int component
-           = mg_target_component[fe.system_to_component_index(i).first];
-         if (mg_selected[component])
-           {
-             const unsigned int level_start
-               = mg_component_start[level][component];
-             dst(global_dof_indices[i] - offset)
-               += src[level](level_dof_indices[i] - level_start);
-           }
-       }
+      typedef std::map<unsigned int, unsigned int>::const_iterator IT;
+      for (IT i=copy_to_and_from_indices[level].begin();
+        i != copy_to_and_from_indices[level].end(); ++i)
+             dst(i->first) += src[level](i->second);
     }
 }
 
index e0710af3889b733756db8cf7ae5b899551bdf2e4..34562a2e90ddc71436662785dc021bc6c1827d1c 100644 (file)
@@ -16,7 +16,7 @@
 
 #include <base/logstream.h>
 
-#include "iostream"
+#include <iostream>
 
 DEAL_II_NAMESPACE_OPEN
 
@@ -141,8 +141,11 @@ Multigrid<VECTOR>::level_v_step(const unsigned int level)
     deallog << "Smoothing on     level " << level << std::endl;
                                   // smoothing of the residual by
                                   // modifying s
+//  defect[level].print(std::cout, 2,false);
+//  std::cout<<std::endl;
   pre_smooth->smooth(level, solution[level], defect[level]);
-  
+//  solution[level].print(std::cout, 2,false);
+
   if (debug>2)
     deallog << "Solution norm          " << solution[level].l2_norm()
            << std::endl;
@@ -151,6 +154,8 @@ Multigrid<VECTOR>::level_v_step(const unsigned int level)
     deallog << "Residual on      level " << level << std::endl;
                                   // t = A*solution[level]
   matrix->vmult(level, t[level], solution[level]);
+//  std::cout<<std::endl;
+//  t[level].print(std::cout, 2,false);
   
                                   // make t rhs of lower level The
                                   // non-refined parts of the
@@ -161,13 +166,18 @@ Multigrid<VECTOR>::level_v_step(const unsigned int level)
     {
       t[l-1] = 0.;
       if (l==level && edge_out != 0)
-       edge_out->vmult_add(level, t[level], solution[level]);
-      
+      {
+        edge_out->vmult_add(level, t[level], solution[level]);
+        deallog << "Norm t[" << level << "] " << t[level].l2_norm() << std::endl;
+      }
+
       if (l==level && edge_down != 0)
        edge_down->vmult(level, t[level-1], solution[level]);
       
       transfer->restrict_and_add (l, t[l-1], t[l]);
+      deallog << "restrict t[" << l-1 << "] " << t[l-1].l2_norm() << std::endl;
       defect[l-1] -= t[l-1];
+      deallog << "defect d[" << l-1 << "] " << defect[l-1].l2_norm() << std::endl;
     }
   
                                   // do recursion
@@ -198,14 +208,19 @@ Multigrid<VECTOR>::level_v_step(const unsigned int level)
    }
   
   if (debug>2)
-    deallog << "V-cycle  Defect norm   " << defect[level].l2_norm()
+    deallog << "V-cycle  Defect norm hier  " << defect[level].l2_norm()
            << std::endl;
   
   if (debug>1)
     deallog << "Smoothing on     level " << level << std::endl;
                                   // post-smoothing
+
+//  std::cout<<std::endl;
+//  defect[level].print(std::cout, 2,false);
   post_smooth->smooth(level, solution[level], defect[level]);
-  
+//  solution[level].print(std::cout, 2,false);
+//  std::cout<<std::endl;
+
   if (debug>2)
     deallog << "Solution norm          " << solution[level].l2_norm()
            << std::endl;

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.