]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Allow for FE_DGQLegendre in MGTransferMatrixFree 3773/head
authorMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Wed, 11 Jan 2017 21:58:42 +0000 (22:58 +0100)
committerMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Wed, 11 Jan 2017 21:58:42 +0000 (22:58 +0100)
include/deal.II/multigrid/mg_transfer_internal.h
include/deal.II/multigrid/mg_transfer_matrix_free.h
source/multigrid/mg_transfer_internal.cc
source/multigrid/mg_transfer_matrix_free.cc
tests/multigrid/transfer_matrix_free_07.cc [new file with mode: 0644]
tests/multigrid/transfer_matrix_free_07.with_mpi=true.with_p4est=true.with_trilinos=true.mpirun=5.output [new file with mode: 0644]

index 80f9d291bcd54ea0086c7c89252828fbe9434d2e..0f5f6c97ebd3c2d96ae63395676b2c0ba955f2d1 100644 (file)
@@ -20,7 +20,6 @@
 #include <deal.II/base/mg_level_object.h>
 #include <deal.II/dofs/dof_handler.h>
 #include <deal.II/lac/la_parallel_vector.h>
-#include <deal.II/matrix_free/shape_info.h>
 #include <deal.II/multigrid/mg_constrained_dofs.h>
 
 DEAL_II_NAMESPACE_OPEN
@@ -87,11 +86,18 @@ namespace internal
        */
       unsigned int n_child_cell_dofs;
 
+      /**
+       * Holds the numbering between the numbering of degrees of freedom in
+       * the finite element and the lexicographic numbering needed for the
+       * tensor product application.
+       */
+      std::vector<unsigned int> lexicographic_numbering;
+
       /**
        * Holds the one-dimensional embedding (prolongation) matrix from mother
-       * element to the children.
+       * element to all the children.
        */
-      internal::MatrixFreeFunctions::ShapeInfo<Number> shape_info;
+      AlignedVector<VectorizedArray<Number> > prolongation_matrix_1d;
 
     };
 
index 96ac877e9e09d35bb5687b7e539b3c81b487d02e..26ce0935ed63c7e4ca672d3e403d45d79714d778 100644 (file)
@@ -23,7 +23,6 @@
 #include <deal.II/multigrid/mg_constrained_dofs.h>
 #include <deal.II/base/mg_level_object.h>
 #include <deal.II/multigrid/mg_transfer.h>
-#include <deal.II/matrix_free/shape_info.h>
 
 #include <deal.II/dofs/dof_handler.h>
 
@@ -188,9 +187,9 @@ private:
 
   /**
    * Holds the one-dimensional embedding (prolongation) matrix from mother
-   * element to the children.
+   * element to all the children.
    */
-  internal::MatrixFreeFunctions::ShapeInfo<Number> shape_info;
+  AlignedVector<VectorizedArray<Number> > prolongation_matrix_1d;
 
   /**
    * Holds the temporary values for the tensor evaluation
index c560bb5e595f1e9681811e652e209e5c11c4ded8..7639c8d5953da24b90a3674e96d2e5f1caf0a1cd 100644 (file)
@@ -17,6 +17,7 @@
 #include <deal.II/distributed/tria.h>
 #include <deal.II/dofs/dof_tools.h>
 #include <deal.II/fe/fe_tools.h>
+#include <deal.II/matrix_free/shape_info.h>
 #include <deal.II/multigrid/mg_transfer_internal.h>
 
 DEAL_II_NAMESPACE_OPEN
@@ -414,35 +415,30 @@ namespace internal
           renumbering[fe.dofs_per_cell-fe.dofs_per_vertex] = fe.dofs_per_vertex;
       }
 
-      // step 1.3: create a 1D quadrature formula from the finite element that
-      // collects the support points of the basis functions on the two children.
+      // step 1.3: create a dummy 1D quadrature formula to extract the
+      // lexicographic numbering for the elements
       std::vector<Point<1> > basic_support_points = fe.get_unit_support_points();
       Assert(fe.dofs_per_vertex == 0 || fe.dofs_per_vertex == 1,
              ExcNotImplemented());
-      std::vector<Point<1> > points_refined(fe.dofs_per_vertex > 0 ?
+      const unsigned int shift = fe.dofs_per_cell - fe.dofs_per_vertex;
+      const unsigned int n_child_dofs_1d = (fe.dofs_per_vertex > 0 ?
                                             (2 * fe.dofs_per_cell - 1) :
                                             (2 * fe.dofs_per_cell));
-      const unsigned int shift = fe.dofs_per_cell - fe.dofs_per_vertex;
-      for (unsigned int c=0; c<GeometryInfo<1>::max_children_per_cell; ++c)
-        for (unsigned int j=0; j<basic_support_points.size(); ++j)
-          points_refined[shift*c+j][0] =
-            c*0.5 + 0.5 * basic_support_points[renumbering[j]][0];
-
-      const unsigned int n_child_dofs_1d = points_refined.size();
 
       elem_info.n_child_cell_dofs = elem_info.n_components*Utilities::fixed_power<dim>(n_child_dofs_1d);
+      const Quadrature<1> dummy_quadrature(std::vector<Point<1> >(1, Point<1>()));
+      internal::MatrixFreeFunctions::ShapeInfo<Number> shape_info;
+      shape_info.reinit(dummy_quadrature, mg_dof.get_fe(), 0);
+      elem_info.lexicographic_numbering = shape_info.lexicographic_numbering;
 
-      // step 1.4: evaluate the polynomials and store the data in ShapeInfo
-      const Quadrature<1> quadrature(points_refined);
-      elem_info.shape_info.reinit(quadrature, mg_dof.get_fe(), 0);
+      // step 1.4: get the 1d prolongation matrix and combine from both children
+      elem_info.prolongation_matrix_1d.resize(fe.dofs_per_cell*n_child_dofs_1d);
 
       for (unsigned int c=0; c<GeometryInfo<1>::max_children_per_cell; ++c)
         for (unsigned int i=0; i<fe.dofs_per_cell; ++i)
           for (unsigned int j=0; j<fe.dofs_per_cell; ++j)
-            Assert(std::abs(elem_info.shape_info.shape_values[i*n_child_dofs_1d+j+c*shift][0] -
-                            fe.get_prolongation_matrix(c)(renumbering[j],renumbering[i]))
-                   < std::max(2.*(double)std::numeric_limits<Number>::epsilon(),1e-12),
-                   ExcInternalError());
+            elem_info.prolongation_matrix_1d[i*n_child_dofs_1d+j+c*shift] =
+              fe.get_prolongation_matrix(c)(renumbering[j],renumbering[i]);
     }
 
 
@@ -567,7 +563,7 @@ namespace internal
                       ghosted_level_dofs.push_back(local_dof_indices[i]);
 
                   add_child_indices<dim>(c, fe.dofs_per_cell - fe.dofs_per_vertex,
-                                         fe.degree, elem_info.shape_info.lexicographic_numbering,
+                                         fe.degree, elem_info.lexicographic_numbering,
                                          local_dof_indices,
                                          &next_indices[start_index]);
 
@@ -598,7 +594,7 @@ namespace internal
                       if (mg_constrained_dofs != 0)
                         for (unsigned int i=0; i<mg_dof.get_fe().dofs_per_cell; ++i)
                           if (mg_constrained_dofs->is_boundary_index(level,
-                                                                     local_dof_indices[elem_info.shape_info.lexicographic_numbering[i]]))
+                                                                     local_dof_indices[elem_info.lexicographic_numbering[i]]))
                             dirichlet_indices[level][child_index].push_back(i);
                     }
                 }
@@ -626,14 +622,14 @@ namespace internal
                   global_level_dof_indices_l0.resize(start_index+elem_info.n_child_cell_dofs,
                                                      numbers::invalid_dof_index);
                   add_child_indices<dim>(0, fe.dofs_per_cell - fe.dofs_per_vertex,
-                                         fe.degree, elem_info.shape_info.lexicographic_numbering,
+                                         fe.degree, elem_info.lexicographic_numbering,
                                          local_dof_indices,
                                          &global_level_dof_indices_l0[start_index]);
 
                   dirichlet_indices[0].push_back(std::vector<unsigned short>());
                   if (mg_constrained_dofs != 0)
                     for (unsigned int i=0; i<mg_dof.get_fe().dofs_per_cell; ++i)
-                      if (mg_constrained_dofs->is_boundary_index(0, local_dof_indices[elem_info.shape_info.lexicographic_numbering[i]]))
+                      if (mg_constrained_dofs->is_boundary_index(0, local_dof_indices[elem_info.lexicographic_numbering[i]]))
                         dirichlet_indices[0].back().push_back(i);
                 }
             }
index fefe0ec73b864a607c5ce5814f4f924e0ecb3a17..7bf82ed31b7715ae358c42eebd871b510966f6b2 100644 (file)
@@ -27,7 +27,6 @@
 #include <deal.II/multigrid/mg_transfer_matrix_free.h>
 #include <deal.II/multigrid/mg_transfer_internal.h>
 
-#include <deal.II/matrix_free/shape_info.h>
 #include <deal.II/matrix_free/fe_evaluation.h>
 
 #include <algorithm>
@@ -85,7 +84,7 @@ void MGTransferMatrixFree<dim,Number>::clear ()
   level_dof_indices.clear();
   parent_child_connect.clear();
   n_owned_level_cells.clear();
-  shape_info = internal::MatrixFreeFunctions::ShapeInfo<Number>();
+  prolongation_matrix_1d.clear();
   evaluation_data.clear();
   weights_on_refined.clear();
 }
@@ -116,7 +115,7 @@ void MGTransferMatrixFree<dim,Number>::build
   element_is_continuous    = elem_info.element_is_continuous;
   n_components             = elem_info.n_components;
   n_child_cell_dofs        = elem_info.n_child_cell_dofs;
-  shape_info               = elem_info.shape_info;
+  prolongation_matrix_1d   = elem_info.prolongation_matrix_1d;
 
 
   // reshuffle into aligned vector of vectorized arrays
@@ -387,16 +386,14 @@ void MGTransferMatrixFree<dim,Number>
         }
 
       // perform tensorized operation
-      Assert(shape_info.element_type ==
-             internal::MatrixFreeFunctions::tensor_symmetric, ExcNotImplemented());
       if (element_is_continuous)
         {
-          AssertDimension(shape_info.shape_val_evenodd.size(),
-                          (degree+1)*(degree+1));
-          typedef internal::EvaluatorTensorProduct<internal::evaluate_evenodd,dim,degree,2*degree+1,VectorizedArray<Number> > Evaluator;
-          Evaluator evaluator(shape_info.shape_val_evenodd,
-                              shape_info.shape_val_evenodd,
-                              shape_info.shape_val_evenodd);
+          AssertDimension(prolongation_matrix_1d.size(),
+                          (2*degree+1)*(degree+1));
+          typedef internal::EvaluatorTensorProduct<internal::evaluate_general,dim,degree,2*degree+1,VectorizedArray<Number> > Evaluator;
+          Evaluator evaluator(prolongation_matrix_1d,
+                              prolongation_matrix_1d,
+                              prolongation_matrix_1d);
           perform_tensorized_op<dim,Evaluator,Number,true>(evaluator,
                                                            n_child_cell_dofs,
                                                            n_components,
@@ -407,12 +404,12 @@ void MGTransferMatrixFree<dim,Number>
         }
       else
         {
-          AssertDimension(shape_info.shape_val_evenodd.size(),
-                          (degree+1)*(degree+1));
-          typedef internal::EvaluatorTensorProduct<internal::evaluate_evenodd,dim,degree,2*degree+2,VectorizedArray<Number> > Evaluator;
-          Evaluator evaluator(shape_info.shape_val_evenodd,
-                              shape_info.shape_val_evenodd,
-                              shape_info.shape_val_evenodd);
+          AssertDimension(prolongation_matrix_1d.size(),
+                          2*(degree+1)*(degree+1));
+          typedef internal::EvaluatorTensorProduct<internal::evaluate_general,dim,degree,2*(degree+1),VectorizedArray<Number> > Evaluator;
+          Evaluator evaluator(prolongation_matrix_1d,
+                              prolongation_matrix_1d,
+                              prolongation_matrix_1d);
           perform_tensorized_op<dim,Evaluator,Number,true>(evaluator,
                                                            n_child_cell_dofs,
                                                            n_components,
@@ -464,16 +461,14 @@ void MGTransferMatrixFree<dim,Number>
       }
 
       // perform tensorized operation
-      Assert(shape_info.element_type ==
-             internal::MatrixFreeFunctions::tensor_symmetric, ExcNotImplemented());
       if (element_is_continuous)
         {
-          AssertDimension(shape_info.shape_val_evenodd.size(),
-                          (degree+1)*(degree+1));
-          typedef internal::EvaluatorTensorProduct<internal::evaluate_evenodd,dim,degree,2*degree+1,VectorizedArray<Number> > Evaluator;
-          Evaluator evaluator(shape_info.shape_val_evenodd,
-                              shape_info.shape_val_evenodd,
-                              shape_info.shape_val_evenodd);
+          AssertDimension(prolongation_matrix_1d.size(),
+                          (2*degree+1)*(degree+1));
+          typedef internal::EvaluatorTensorProduct<internal::evaluate_general,dim,degree,2*degree+1,VectorizedArray<Number> > Evaluator;
+          Evaluator evaluator(prolongation_matrix_1d,
+                              prolongation_matrix_1d,
+                              prolongation_matrix_1d);
           weight_dofs_on_child<dim,degree,Number>(&weights_on_refined[from_level-1][(cell/vec_size)*three_to_dim],
                                                   n_components,
                                                   &evaluation_data[0]);
@@ -484,12 +479,12 @@ void MGTransferMatrixFree<dim,Number>
         }
       else
         {
-          AssertDimension(shape_info.shape_val_evenodd.size(),
-                          (degree+1)*(degree+1));
-          typedef internal::EvaluatorTensorProduct<internal::evaluate_evenodd,dim,degree,2*degree+2,VectorizedArray<Number> > Evaluator;
-          Evaluator evaluator(shape_info.shape_val_evenodd,
-                              shape_info.shape_val_evenodd,
-                              shape_info.shape_val_evenodd);
+          AssertDimension(prolongation_matrix_1d.size(),
+                          2*(degree+1)*(degree+1));
+          typedef internal::EvaluatorTensorProduct<internal::evaluate_general,dim,degree,2*(degree+1),VectorizedArray<Number> > Evaluator;
+          Evaluator evaluator(prolongation_matrix_1d,
+                              prolongation_matrix_1d,
+                              prolongation_matrix_1d);
           perform_tensorized_op<dim,Evaluator,Number,false>(evaluator,
                                                             n_child_cell_dofs,
                                                             n_components,
@@ -534,7 +529,7 @@ MGTransferMatrixFree<dim,Number>::memory_consumption() const
   memory += MemoryConsumption::memory_consumption(level_dof_indices);
   memory += MemoryConsumption::memory_consumption(parent_child_connect);
   memory += MemoryConsumption::memory_consumption(n_owned_level_cells);
-  memory += shape_info.memory_consumption();
+  memory += MemoryConsumption::memory_consumption(prolongation_matrix_1d);
   memory += MemoryConsumption::memory_consumption(evaluation_data);
   memory += MemoryConsumption::memory_consumption(weights_on_refined);
   memory += MemoryConsumption::memory_consumption(dirichlet_indices);
diff --git a/tests/multigrid/transfer_matrix_free_07.cc b/tests/multigrid/transfer_matrix_free_07.cc
new file mode 100644 (file)
index 0000000..b8fed53
--- /dev/null
@@ -0,0 +1,153 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2017 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 MGTransferMatrixFree by comparison with MGTransferPrebuilt on a
+// series of meshes with adaptive meshes for FE_DGQLegendre (except for the
+// different element the same test as transfer_matrix_free_05)
+
+#include "../tests.h"
+#include <deal.II/base/logstream.h>
+#include <deal.II/lac/la_parallel_vector.h>
+#include <deal.II/distributed/tria.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/fe/fe_dgq.h>
+#include <deal.II/multigrid/mg_transfer.h>
+#include <deal.II/multigrid/mg_transfer_matrix_free.h>
+
+
+template <int dim, typename Number>
+void check(const unsigned int fe_degree)
+{
+  deallog.threshold_double(std::max(5e2*(double)std::numeric_limits<Number>::epsilon(),
+                                    1e-11));
+  FE_DGQLegendre<dim> fe(fe_degree);
+  deallog << "FE: " << fe.get_name() << std::endl;
+
+  // run a few different sizes...
+  unsigned int sizes [] = {1, 2, 3};
+  for (unsigned int cycle=0; cycle<sizeof(sizes)/sizeof(unsigned int); ++cycle)
+    {
+      unsigned int n_refinements = 0;
+      unsigned int n_subdiv = sizes[cycle];
+      if (n_subdiv > 1)
+        while (n_subdiv%2 == 0)
+          {
+            n_refinements += 1;
+            n_subdiv /= 2;
+          }
+      n_refinements += 3-dim;
+      if (fe_degree < 3)
+        n_refinements += 1;
+
+      parallel::distributed::Triangulation<dim>
+      tr(MPI_COMM_WORLD,
+         Triangulation<dim>::limit_level_difference_at_vertices,
+         parallel::distributed::Triangulation<dim>::construct_multigrid_hierarchy);
+      GridGenerator::subdivided_hyper_cube(tr, n_subdiv);
+      tr.refine_global(n_refinements);
+
+      // adaptive refinement into a circle
+      for (typename Triangulation<dim>::active_cell_iterator cell=tr.begin_active(); cell != tr.end(); ++cell)
+        if (cell->is_locally_owned() &&
+            cell->center().norm() < 0.5)
+          cell->set_refine_flag();
+      tr.execute_coarsening_and_refinement();
+      for (typename Triangulation<dim>::active_cell_iterator cell=tr.begin_active(); cell != tr.end(); ++cell)
+        if (cell->is_locally_owned() &&
+            cell->center().norm() > 0.3 && cell->center().norm() < 0.4)
+          cell->set_refine_flag();
+      tr.execute_coarsening_and_refinement();
+      for (typename Triangulation<dim>::active_cell_iterator cell=tr.begin_active(); cell != tr.end(); ++cell)
+        if (cell->is_locally_owned() &&
+            cell->center().norm() > 0.33 && cell->center().norm() < 0.37)
+          cell->set_refine_flag();
+      tr.execute_coarsening_and_refinement();
+
+      deallog << "no. cells: " << tr.n_global_active_cells() << std::endl;
+
+      DoFHandler<dim> mgdof(tr);
+      mgdof.distribute_dofs(fe);
+      mgdof.distribute_mg_dofs(fe);
+
+      // build reference
+      MGTransferPrebuilt<LinearAlgebra::distributed::Vector<double> > transfer_ref;
+      transfer_ref.build_matrices(mgdof);
+
+      // build matrix-free transfer
+      MGTransferMatrixFree<dim, Number> transfer;
+      transfer.build(mgdof);
+
+      // check prolongation for all levels using random vector
+      for (unsigned int level=1; level<mgdof.get_triangulation().n_global_levels(); ++level)
+        {
+          LinearAlgebra::distributed::Vector<Number> v1, v2;
+          LinearAlgebra::distributed::Vector<double> v1_cpy, v2_cpy, v3;
+          v1.reinit(mgdof.locally_owned_mg_dofs(level-1), MPI_COMM_WORLD);
+          v2.reinit(mgdof.locally_owned_mg_dofs(level), MPI_COMM_WORLD);
+          v3.reinit(mgdof.locally_owned_mg_dofs(level), MPI_COMM_WORLD);
+          for (unsigned int i=0; i<v1.local_size(); ++i)
+            v1.local_element(i) = (double)Testing::rand()/RAND_MAX;
+          v1_cpy = v1;
+          transfer.prolongate(level, v2, v1);
+          transfer_ref.prolongate(level, v3, v1_cpy);
+          v2_cpy = v2;
+          v3 -= v2_cpy;
+          deallog << "Diff prolongate   l" << level << ": " << v3.l2_norm() << std::endl;
+        }
+
+      // check restriction for all levels using random vector
+      for (unsigned int level=1; level<mgdof.get_triangulation().n_global_levels(); ++level)
+        {
+          LinearAlgebra::distributed::Vector<Number> v1, v2;
+          LinearAlgebra::distributed::Vector<double> v1_cpy, v2_cpy, v3;
+          v1.reinit(mgdof.locally_owned_mg_dofs(level), MPI_COMM_WORLD);
+          v2.reinit(mgdof.locally_owned_mg_dofs(level-1), MPI_COMM_WORLD);
+          v3.reinit(mgdof.locally_owned_mg_dofs(level-1), MPI_COMM_WORLD);
+          for (unsigned int i=0; i<v1.local_size(); ++i)
+            v1.local_element(i) = (double)Testing::rand()/RAND_MAX;
+          v1_cpy = v1;
+          transfer.restrict_and_add(level, v2, v1);
+          transfer_ref.restrict_and_add(level, v3, v1_cpy);
+          v2_cpy = v2;
+          v3 -= v2_cpy;
+          deallog << "Diff restrict     l" << level << ": " << v3.l2_norm() << std::endl;
+
+          v2 = 1.;
+          v3 = 1.;
+          transfer.restrict_and_add(level, v2, v1);
+          transfer_ref.restrict_and_add(level, v3, v1_cpy);
+          v2_cpy = v2;
+          v3 -= v2_cpy;
+          deallog << "Diff restrict add l" << level << ": " << v3.l2_norm() << std::endl;
+        }
+      deallog << std::endl;
+    }
+}
+
+
+int main(int argc, char **argv)
+{
+  // no threading in this test...
+  Utilities::MPI::MPI_InitFinalize mpi(argc, argv, 1);
+  mpi_initlog();
+
+  check<2,double>(1);
+  check<2,double>(3);
+  check<3,double>(1);
+  check<3,double>(3);
+  check<2,float> (2);
+  check<3,float> (2);
+}
diff --git a/tests/multigrid/transfer_matrix_free_07.with_mpi=true.with_p4est=true.with_trilinos=true.mpirun=5.output b/tests/multigrid/transfer_matrix_free_07.with_mpi=true.with_p4est=true.with_trilinos=true.mpirun=5.output
new file mode 100644 (file)
index 0000000..a01f58a
--- /dev/null
@@ -0,0 +1,259 @@
+
+DEAL::FE: FE_DGQLegendre<2>(1)
+DEAL::no. cells: 88
+DEAL::Diff prolongate   l1: 0
+DEAL::Diff prolongate   l2: 0
+DEAL::Diff prolongate   l3: 0
+DEAL::Diff prolongate   l4: 0
+DEAL::Diff prolongate   l5: 0
+DEAL::Diff restrict     l1: 0
+DEAL::Diff restrict add l1: 0
+DEAL::Diff restrict     l2: 0
+DEAL::Diff restrict add l2: 0
+DEAL::Diff restrict     l3: 0
+DEAL::Diff restrict add l3: 0
+DEAL::Diff restrict     l4: 0
+DEAL::Diff restrict add l4: 0
+DEAL::Diff restrict     l5: 0
+DEAL::Diff restrict add l5: 0
+DEAL::
+DEAL::no. cells: 250
+DEAL::Diff prolongate   l1: 0
+DEAL::Diff prolongate   l2: 0
+DEAL::Diff prolongate   l3: 0
+DEAL::Diff prolongate   l4: 0
+DEAL::Diff prolongate   l5: 0
+DEAL::Diff prolongate   l6: 0
+DEAL::Diff restrict     l1: 0
+DEAL::Diff restrict add l1: 0
+DEAL::Diff restrict     l2: 0
+DEAL::Diff restrict add l2: 0
+DEAL::Diff restrict     l3: 0
+DEAL::Diff restrict add l3: 0
+DEAL::Diff restrict     l4: 0
+DEAL::Diff restrict add l4: 0
+DEAL::Diff restrict     l5: 0
+DEAL::Diff restrict add l5: 0
+DEAL::Diff restrict     l6: 0
+DEAL::Diff restrict add l6: 0
+DEAL::
+DEAL::no. cells: 510
+DEAL::Diff prolongate   l1: 0
+DEAL::Diff prolongate   l2: 0
+DEAL::Diff prolongate   l3: 0
+DEAL::Diff prolongate   l4: 0
+DEAL::Diff prolongate   l5: 0
+DEAL::Diff restrict     l1: 0
+DEAL::Diff restrict add l1: 0
+DEAL::Diff restrict     l2: 0
+DEAL::Diff restrict add l2: 0
+DEAL::Diff restrict     l3: 0
+DEAL::Diff restrict add l3: 0
+DEAL::Diff restrict     l4: 0
+DEAL::Diff restrict add l4: 0
+DEAL::Diff restrict     l5: 0
+DEAL::Diff restrict add l5: 0
+DEAL::
+DEAL::FE: FE_DGQLegendre<2>(3)
+DEAL::no. cells: 34
+DEAL::Diff prolongate   l1: 0
+DEAL::Diff prolongate   l2: 0
+DEAL::Diff prolongate   l3: 0
+DEAL::Diff prolongate   l4: 0
+DEAL::Diff restrict     l1: 0
+DEAL::Diff restrict add l1: 0
+DEAL::Diff restrict     l2: 0
+DEAL::Diff restrict add l2: 0
+DEAL::Diff restrict     l3: 0
+DEAL::Diff restrict add l3: 0
+DEAL::Diff restrict     l4: 0
+DEAL::Diff restrict add l4: 0
+DEAL::
+DEAL::no. cells: 88
+DEAL::Diff prolongate   l1: 0
+DEAL::Diff prolongate   l2: 0
+DEAL::Diff prolongate   l3: 0
+DEAL::Diff prolongate   l4: 0
+DEAL::Diff prolongate   l5: 0
+DEAL::Diff restrict     l1: 0
+DEAL::Diff restrict add l1: 0
+DEAL::Diff restrict     l2: 0
+DEAL::Diff restrict add l2: 0
+DEAL::Diff restrict     l3: 0
+DEAL::Diff restrict add l3: 0
+DEAL::Diff restrict     l4: 0
+DEAL::Diff restrict add l4: 0
+DEAL::Diff restrict     l5: 0
+DEAL::Diff restrict add l5: 0
+DEAL::
+DEAL::no. cells: 141
+DEAL::Diff prolongate   l1: 0
+DEAL::Diff prolongate   l2: 0
+DEAL::Diff prolongate   l3: 0
+DEAL::Diff prolongate   l4: 0
+DEAL::Diff restrict     l1: 0
+DEAL::Diff restrict add l1: 0
+DEAL::Diff restrict     l2: 0
+DEAL::Diff restrict add l2: 0
+DEAL::Diff restrict     l3: 0
+DEAL::Diff restrict add l3: 0
+DEAL::Diff restrict     l4: 0
+DEAL::Diff restrict add l4: 0
+DEAL::
+DEAL::FE: FE_DGQLegendre<3>(1)
+DEAL::no. cells: 15
+DEAL::Diff prolongate   l1: 0
+DEAL::Diff prolongate   l2: 0
+DEAL::Diff restrict     l1: 0
+DEAL::Diff restrict add l1: 0
+DEAL::Diff restrict     l2: 0
+DEAL::Diff restrict add l2: 0
+DEAL::
+DEAL::no. cells: 694
+DEAL::Diff prolongate   l1: 0
+DEAL::Diff prolongate   l2: 0
+DEAL::Diff prolongate   l3: 0
+DEAL::Diff prolongate   l4: 0
+DEAL::Diff prolongate   l5: 0
+DEAL::Diff restrict     l1: 0
+DEAL::Diff restrict add l1: 0
+DEAL::Diff restrict     l2: 0
+DEAL::Diff restrict add l2: 0
+DEAL::Diff restrict     l3: 0
+DEAL::Diff restrict add l3: 0
+DEAL::Diff restrict     l4: 0
+DEAL::Diff restrict add l4: 0
+DEAL::Diff restrict     l5: 0
+DEAL::Diff restrict add l5: 0
+DEAL::
+DEAL::no. cells: 1791
+DEAL::Diff prolongate   l1: 0
+DEAL::Diff prolongate   l2: 0
+DEAL::Diff prolongate   l3: 0
+DEAL::Diff prolongate   l4: 0
+DEAL::Diff restrict     l1: 0
+DEAL::Diff restrict add l1: 0
+DEAL::Diff restrict     l2: 0
+DEAL::Diff restrict add l2: 0
+DEAL::Diff restrict     l3: 0
+DEAL::Diff restrict add l3: 0
+DEAL::Diff restrict     l4: 0
+DEAL::Diff restrict add l4: 0
+DEAL::
+DEAL::FE: FE_DGQLegendre<3>(3)
+DEAL::no. cells: 1
+DEAL::
+DEAL::no. cells: 15
+DEAL::Diff prolongate   l1: 0
+DEAL::Diff prolongate   l2: 0
+DEAL::Diff restrict     l1: 0
+DEAL::Diff restrict add l1: 0
+DEAL::Diff restrict     l2: 0
+DEAL::Diff restrict add l2: 0
+DEAL::
+DEAL::no. cells: 223
+DEAL::Diff prolongate   l1: 0
+DEAL::Diff prolongate   l2: 0
+DEAL::Diff prolongate   l3: 0
+DEAL::Diff restrict     l1: 0
+DEAL::Diff restrict add l1: 0
+DEAL::Diff restrict     l2: 0
+DEAL::Diff restrict add l2: 0
+DEAL::Diff restrict     l3: 0
+DEAL::Diff restrict add l3: 0
+DEAL::
+DEAL::FE: FE_DGQLegendre<2>(2)
+DEAL::no. cells: 88
+DEAL::Diff prolongate   l1: 0
+DEAL::Diff prolongate   l2: 0
+DEAL::Diff prolongate   l3: 0
+DEAL::Diff prolongate   l4: 0
+DEAL::Diff prolongate   l5: 0
+DEAL::Diff restrict     l1: 0
+DEAL::Diff restrict add l1: 0
+DEAL::Diff restrict     l2: 0
+DEAL::Diff restrict add l2: 0
+DEAL::Diff restrict     l3: 0
+DEAL::Diff restrict add l3: 0
+DEAL::Diff restrict     l4: 0
+DEAL::Diff restrict add l4: 0
+DEAL::Diff restrict     l5: 0
+DEAL::Diff restrict add l5: 0
+DEAL::
+DEAL::no. cells: 250
+DEAL::Diff prolongate   l1: 0
+DEAL::Diff prolongate   l2: 0
+DEAL::Diff prolongate   l3: 0
+DEAL::Diff prolongate   l4: 0
+DEAL::Diff prolongate   l5: 0
+DEAL::Diff prolongate   l6: 0
+DEAL::Diff restrict     l1: 0
+DEAL::Diff restrict add l1: 0
+DEAL::Diff restrict     l2: 0
+DEAL::Diff restrict add l2: 0
+DEAL::Diff restrict     l3: 0
+DEAL::Diff restrict add l3: 0
+DEAL::Diff restrict     l4: 0
+DEAL::Diff restrict add l4: 0
+DEAL::Diff restrict     l5: 0
+DEAL::Diff restrict add l5: 0
+DEAL::Diff restrict     l6: 0
+DEAL::Diff restrict add l6: 0
+DEAL::
+DEAL::no. cells: 510
+DEAL::Diff prolongate   l1: 0
+DEAL::Diff prolongate   l2: 0
+DEAL::Diff prolongate   l3: 0
+DEAL::Diff prolongate   l4: 0
+DEAL::Diff prolongate   l5: 0
+DEAL::Diff restrict     l1: 0
+DEAL::Diff restrict add l1: 0
+DEAL::Diff restrict     l2: 0
+DEAL::Diff restrict add l2: 0
+DEAL::Diff restrict     l3: 0
+DEAL::Diff restrict add l3: 0
+DEAL::Diff restrict     l4: 0
+DEAL::Diff restrict add l4: 0
+DEAL::Diff restrict     l5: 0
+DEAL::Diff restrict add l5: 0
+DEAL::
+DEAL::FE: FE_DGQLegendre<3>(2)
+DEAL::no. cells: 15
+DEAL::Diff prolongate   l1: 0
+DEAL::Diff prolongate   l2: 0
+DEAL::Diff restrict     l1: 0
+DEAL::Diff restrict add l1: 0
+DEAL::Diff restrict     l2: 0
+DEAL::Diff restrict add l2: 0
+DEAL::
+DEAL::no. cells: 694
+DEAL::Diff prolongate   l1: 0
+DEAL::Diff prolongate   l2: 0
+DEAL::Diff prolongate   l3: 0
+DEAL::Diff prolongate   l4: 0
+DEAL::Diff prolongate   l5: 0
+DEAL::Diff restrict     l1: 0
+DEAL::Diff restrict add l1: 0
+DEAL::Diff restrict     l2: 0
+DEAL::Diff restrict add l2: 0
+DEAL::Diff restrict     l3: 0
+DEAL::Diff restrict add l3: 0
+DEAL::Diff restrict     l4: 0
+DEAL::Diff restrict add l4: 0
+DEAL::Diff restrict     l5: 0
+DEAL::Diff restrict add l5: 0
+DEAL::
+DEAL::no. cells: 1791
+DEAL::Diff prolongate   l1: 0
+DEAL::Diff prolongate   l2: 0
+DEAL::Diff prolongate   l3: 0
+DEAL::Diff prolongate   l4: 0
+DEAL::Diff restrict     l1: 0
+DEAL::Diff restrict add l1: 0
+DEAL::Diff restrict     l2: 0
+DEAL::Diff restrict add l2: 0
+DEAL::Diff restrict     l3: 0
+DEAL::Diff restrict add l3: 0
+DEAL::Diff restrict     l4: 0
+DEAL::Diff restrict add l4: 0
+DEAL::

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.