]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add Joshua White's contribution of a MappingQEulerian class.
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Tue, 5 Feb 2008 21:07:09 +0000 (21:07 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Tue, 5 Feb 2008 21:07:09 +0000 (21:07 +0000)
git-svn-id: https://svn.dealii.org/trunk@15712 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/fe/mapping_q_eulerian.h [new file with mode: 0644]
deal.II/deal.II/source/fe/mapping_q_eulerian.cc [new file with mode: 0644]
deal.II/doc/authors.html
deal.II/doc/news/changes.h
tests/fe/Makefile
tests/fe/mapping_q_eulerian.cc [new file with mode: 0644]
tests/fe/mapping_q_eulerian/cmp/generic [new file with mode: 0644]

diff --git a/deal.II/deal.II/include/fe/mapping_q_eulerian.h b/deal.II/deal.II/include/fe/mapping_q_eulerian.h
new file mode 100644 (file)
index 0000000..17f7e61
--- /dev/null
@@ -0,0 +1,216 @@
+//---------------------------------------------------------------------------
+//    $Id: mapping_q.h 15711 2008-02-05 20:43:14Z bangerth $
+//    Version: $Name$
+//
+//    Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008 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.
+//
+//---------------------------------------------------------------------------
+
+#ifndef __deal2__mapping_q_eulerian_h
+#define __deal2__mapping_q_eulerian_h
+
+#include <base/smartpointer.h>
+#include <grid/tria_iterator.h>
+#include <dofs/dof_handler.h>
+#include <dofs/dof_accessor.h>
+#include <fe/fe.h>
+#include <fe/fe_values.h>
+#include <fe/mapping_q.h> 
+
+
+DEAL_II_NAMESPACE_OPEN
+
+
+/*!@addtogroup mapping */
+/*@{*/
+
+/**
+ * This class is an extension of the MappingQ1Eulerian
+ * class to higher order Qp mappings.  It is useful
+ * when one wants to calculate shape function information on
+ * a domain that is deforming as the computation proceeds.
+ *
+ * <h3>Usage</h3> 
+ *
+ * The constructor of this class takes three arguments: the polynomial
+ * degree of the desire Qp mapping, a reference to
+ * the vector that defines the mapping from the initial
+ * configuration to the current configuration, and a reference to the
+ * DoFHandler. The most common case is to use the solution 
+ * vector for the problem under consideration as the shift vector.  
+ * The key reqirement is that the number of components 
+ * of the given vector field be equal to (or possibly greater than) the 
+ * number of space dimensions. If there are more components than space 
+ * dimensions (for example, if one is working with a coupled problem
+ * where there are additional solution variables), the
+ * first <tt>dim</tt> components are assumed to represent the displacement 
+ * field, and the remaining components are ignored.  If this assumption
+ * does not hold one may need to set up a separate DoFHandler on
+ * the triangulation and associate the desired shift vector to it.
+ *
+ * Typically, the DoFHandler operates on a finite element that
+ * is constructed as a system element (FESystem) from continuous FE_Q()
+ * objects. An example is shown below:
+ * @verbatim
+ *    FESystem<dim> fe(FE_Q<dim>(2), dim, FE_Q<dim>(1), 1);
+ *    DoFHandler<dim> dof_handler(triangulation);
+ *    dof_handler.distribute_dofs(fe);
+ *    Vector<double> soln_vector(dof_handler.n_dofs());
+ *    MappingQEulerian<dim> q2_mapping(2,soln_vector,dof_handler);
+ * @endverbatim
+ *
+ * In this example, our element consists of <tt>(dim+1)</tt> components.
+ * Only the first <tt>dim</tt> components will be used, however, to define
+ * the Q2 mapping.  The remaining components are ignored.
+ *
+ * Note that it is essential to call the distribute_dofs(...) function
+ * before constructing a mapping object.
+ *
+ * Also note that since the vector of shift values and the dof handler are
+ * only associated to this object at construction time, you have to
+ * make sure that whenever you use this object, the given objects
+ * still represent valid data.
+ *
+ * To enable the use of the MappingQ1Eulerian class also in the context
+ * of parallel codes using the PETSc wrapper classes, the type of
+ * the vector can be specified as template parameter <tt>EulerVectorType</tt>
+ * Not specifying this template argument in applications using the PETSc
+ * vector classes leads to the construction of a copy of the vector
+ * which is not acccessible afterwards!
+ *
+ * @author Joshua White, 2008
+ */
+template <int dim, class EulerVectorType = Vector<double> >
+class MappingQEulerian : public MappingQ<dim>
+{
+  public:
+                                     /**
+                                      * Constructor. The first argument is
+                                      * the polynomical degree of the desired
+                                      * Qp mapping.  It then takes a
+                                      * <tt>Vector<double> &</tt> to specify the
+                                      * transformation of the domain
+                                      * from the reference to
+                                      * the current configuration.
+                                      * The organization of the
+                                      * elements in the @p Vector
+                                      * must follow the concept how
+                                      * deal.II stores solutions that
+                                      * are associated to a
+                                      * triangulation.  This is
+                                      * automatically the case if the
+                                      * @p Vector represents the
+                                      * solution of the previous step
+                                      * of a nonlinear problem.
+                                      * Alternatively, the @p Vector
+                                      * can be initialized by
+                                      * <tt>DoFObjectAccessor::set_dof_values()</tt>.
+                                      */
+
+    MappingQEulerian (const unsigned int     degree,
+                      const EulerVectorType  &euler_vector,
+                      const DoFHandler<dim>  &euler_dof_handler);
+
+                                     /**
+                                      * Exception
+                                      */
+
+    DeclException0 (ExcWrongNoOfComponents);
+
+                                    /**
+                                      * Exception
+                                      */
+
+    DeclException0 (ExcInactiveCell);
+    
+                                    /**
+                                      * Exception
+                                      */
+
+    DeclException2 (ExcWrongVectorSize, int, int,
+                   << "Vector has wrong size " << arg1
+                   << "-- expected size " << arg2);
+
+  protected:
+                                     /**
+                                      * Reference to the vector of
+                                      * shifts.
+                                      */
+
+    const EulerVectorType &euler_vector;
+    
+                                     /**
+                                      * Pointer to the DoFHandler to
+                                      * which the mapping vector is
+                                      * associated.
+                                      */
+
+    const SmartPointer<const DoFHandler<dim> > euler_dof_handler;
+
+
+  private:   
+                                     /**
+                                      * Special quadrature rule used
+                                      * to define the support points
+                                      * in the reference configuration.
+                                      */
+
+    class SupportQuadrature : public Quadrature<dim>
+    { 
+      public:
+                                        /**
+                                         * Constructor, with an argument
+                                         * defining the desired polynomial
+                                         * degree.
+                                         */
+
+        SupportQuadrature (const unsigned int map_degree); 
+
+    };
+
+                                    /**
+                                     * A member variable holding the
+                                     * quadrature points in the right
+                                     * order.
+                                     */
+    const SupportQuadrature support_quadrature;
+
+                                    /**
+                                      * FEValues object used to query the
+                                      * the given finite element field
+                                      * at the support points in the
+                                      * reference configuration.
+                                     *
+                                     * The variable is marked as
+                                     * mutable since we have to call
+                                     * FEValues::reinit from
+                                     * compute_mapping_support_points,
+                                     * a function that is 'const'.
+                                      */
+    mutable FEValues<dim> fe_values;
+
+                                     /**
+                                      * Compute the positions of the
+                                      * support points in the current
+                                      * configuration
+                                      */
+    virtual void compute_mapping_support_points(
+      const typename Triangulation<dim>::cell_iterator &cell,
+      std::vector<Point<dim> > &a) const;
+    
+};
+
+/*@}*/
+
+
+DEAL_II_NAMESPACE_CLOSE
+
+
+#endif // __deal2__mapping_q_eulerian_h
+
diff --git a/deal.II/deal.II/source/fe/mapping_q_eulerian.cc b/deal.II/deal.II/source/fe/mapping_q_eulerian.cc
new file mode 100644 (file)
index 0000000..2177ab5
--- /dev/null
@@ -0,0 +1,175 @@
+//---------------------------------------------------------------------------
+//    $Id: mapping_q_eulerian.cc 14038 2006-10-23 02:46:34Z bangerth $
+//    Version: $Name$
+//
+//    Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006, 2008 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.
+//
+//---------------------------------------------------------------------------
+
+
+#include <base/utilities.h>
+#include <base/quadrature_lib.h>
+#include <lac/vector.h>
+#include <lac/petsc_vector.h>
+#include <grid/tria_iterator.h>
+#include <dofs/dof_handler.h>
+#include <dofs/dof_accessor.h>
+#include <fe/fe.h>
+#include <fe/fe_tools.h>
+#include <fe/mapping_q_eulerian.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+
+
+// .... MAPPING Q EULERIAN CONSTRUCTOR
+
+template <int dim, class EulerVectorType>
+MappingQEulerian<dim, EulerVectorType>::
+MappingQEulerian (const unsigned int degree,
+                 const EulerVectorType &euler_vector,
+                 const DoFHandler<dim> &euler_dof_handler)
+               :
+               MappingQ<dim>(degree, true),
+               euler_vector(euler_vector),
+               euler_dof_handler(&euler_dof_handler),
+               support_quadrature(degree),
+               fe_values(euler_dof_handler.get_fe(),
+                         support_quadrature,
+                         update_values | update_q_points)
+{ }
+
+
+// .... SUPPORT QUADRATURE CONSTRUCTOR
+
+template <int dim, class EulerVectorType>
+MappingQEulerian<dim,EulerVectorType>::
+SupportQuadrature::
+SupportQuadrature (const unsigned int map_degree)
+               :
+               Quadrature<dim>(Utilities::fixed_power<dim>(map_degree+1))
+{
+                                  // first we determine the support points
+                                  // on the unit cell in lexicographic order.
+                                  // for this purpose we can use an interated
+                                  // trapezoidal quadrature rule.
+
+  QIterated<dim> q_iterated(QTrapez<1>(),map_degree);
+  const unsigned int n_q_points = q_iterated.n_quadrature_points;
+
+                                  // we then need to define a renumbering 
+                                  // vector that allows us to go from a 
+                                  // lexicographic numbering scheme to a hierarchic
+                                  // one.  this fragment is taking almost verbatim
+                                  // from the MappingQ class.
+
+  std::vector<unsigned int> renumber(n_q_points);
+  std::vector<unsigned int> dpo(dim+1, 1U);
+  for (unsigned int i=1; i<dpo.size(); ++i)
+    dpo[i]=dpo[i-1]*(map_degree-1);
+
+  FETools::lexicographic_to_hierarchic_numbering (
+    FiniteElementData<dim> (dpo, 1, map_degree), renumber);
+
+                                  // finally we assign the quadrature points in the
+                                  // required order.
+
+  for(unsigned int q=0; q<n_q_points; ++q)
+    this->quadrature_points[renumber[q]] = q_iterated.point(q);
+}
+
+
+
+// .... COMPUTE MAPPING SUPPORT POINTS
+
+template <int dim, class EulerVectorType>
+void
+MappingQEulerian<dim, EulerVectorType>::
+compute_mapping_support_points
+(const typename Triangulation<dim>::cell_iterator &cell,
+ std::vector<Point<dim> > &a) const
+{
+
+                                  // first, basic assertion
+                                  // with respect to vector size,
+
+  const unsigned int n_dofs      = euler_dof_handler->n_dofs();
+  const unsigned int vector_size = euler_vector.size();
+
+  Assert (n_dofs == vector_size,ExcWrongVectorSize(vector_size,n_dofs));
+
+                                  // we then transform our tria iterator
+                                  // into a dof iterator so we can 
+                                  // access data not associated with
+                                  // triangulations
+
+  typename DoFHandler<dim>::cell_iterator dof_cell 
+    (const_cast<Triangulation<dim> *> (&(cell->get_triangulation())),
+     cell->level(),
+     cell->index(),
+     euler_dof_handler);
+
+  Assert (dof_cell->active() == true, ExcInactiveCell());
+
+                                  // our quadrature rule is chosen
+                                  // so that each quadrature point
+                                  // corresponds to a support point
+                                  // in the undeformed configuration.
+                                  // we can then query the given
+                                  // displacement field at these points
+                                  // to determine the shift vector that
+                                  // maps the support points to the 
+                                  // deformed configuration.
+
+                                  // we assume that the given field contains
+                                  // dim displacement components, but
+                                  // that there may be other solution
+                                  // components as well (e.g. pressures).
+                                  // this class therefore assumes that the
+                                  // first dim components represent the
+                                  // actual shift vector we need, and simply
+                                  // ignores any components after that.  
+                                  // this implies that the user should order
+                                  // components appropriately, or create a
+                                  // separate dof handler for the displacements.
+
+  const unsigned int n_support_pts = support_quadrature.n_quadrature_points;
+  const unsigned int n_components  = euler_dof_handler->get_fe().n_components();
+
+  Assert (n_components >= dim, ExcWrongNoOfComponents() );
+
+  std::vector<Vector<double> > shift_vector(n_support_pts,Vector<double>(n_components));
+
+                                  // fill shift vector for each support
+                                  // point using an fe_values object.
+
+  fe_values.reinit(dof_cell);
+  fe_values.get_function_values(euler_vector,shift_vector);
+
+                                  // and finally compute the positions of the 
+                                  // support points in the deformed
+                                  // configuration.
+
+  a.resize(n_support_pts);
+  for(unsigned int q=0; q<n_support_pts; ++q)
+    {
+      a[q] = fe_values.quadrature_point(q);
+      for(unsigned int d=0; d<dim; ++d) 
+       a[q](d) += shift_vector[q](d);
+    }
+}
+
+
+
+// explicit instantiation
+template class MappingQEulerian<deal_II_dimension, Vector<double> >;
+#ifdef DEAL_II_USE_PETSC
+template class MappingQEulerian<deal_II_dimension, PETScWrappers::Vector>;
+#endif
+
+DEAL_II_NAMESPACE_CLOSE
index 8cde514ae213f2f33aad414ee455472e4a5e4c58..1e797eedfe5890e94919ebb827bd687447fabeba 100644 (file)
@@ -134,8 +134,13 @@ alphabetical order):
 
 <li><em>Yaqi Wang:</em>
     Some grid generation functions.
+
+<li><em>Joshua White:</em>
+    Arbitrary order Eulerian mappings.
 </ul>
 
+
+
 <h2>deal.II mailing list members</h2>
 
 <p>
index 0fb045ad8d31f14cbc7010d75a11c425423c31f0..3ebaeb842b759ef1d204fd8ffa65d213945fcd90 100644 (file)
@@ -233,6 +233,13 @@ constraints individually.
 
 <ol> 
 
+  <li> <p>New: The MappingQEulerian class provides an arbitrary order Eulerian
+  mapping where the nodes of a mesh are displaced by a previously computed
+  vector field.
+  <br>
+  (Joshua White 2008/02/05)
+  </p></li>
+
   <li> <p>New: The MappingQ class constructor now takes a parameter that determines
   whether a higher order mapping should also be used for interior cells. By default,
   the higher order mapping is only used for cells on the boundary; cells in the
index 5067740a58af86b1f7d59524ee6dd099d9839473..912cf3279a39038f221ee8677e38a8a12fe84af8 100644 (file)
@@ -1,6 +1,6 @@
 ############################################################
 # Makefile,v 1.15 2002/06/13 12:51:13 hartmann Exp
-# Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007 by the deal.II authors
+# Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008 by the deal.II authors
 ############################################################
 
 ############################################################
@@ -22,7 +22,7 @@ default: run-tests
 ############################################################
 
 tests_x=fe_data_test traits fe_tools fe_tools_* mapping \
-        mapping_c1 shapes derivatives function numbering mapping_q1_eulerian \
+        mapping_c1 shapes derivatives function numbering mapping_*_eulerian \
         transfer internals derivatives_face \
        dgq_1 \
        q_* \
diff --git a/tests/fe/mapping_q_eulerian.cc b/tests/fe/mapping_q_eulerian.cc
new file mode 100644 (file)
index 0000000..cba0658
--- /dev/null
@@ -0,0 +1,257 @@
+//----------------------------  mapping_q_eulerian.cc  ---------------------------
+//    mapping_q_eulerian.cc,v 1.3 2003/06/09 16:00:38 wolf Exp
+//    Version: 
+//
+//    Copyright (C) 2008 by the deal.II authors and Joshua White
+//
+//    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.
+//
+//----------------------------  mapping_q_eulerian.cc  ---------------------------
+
+// compute some convergence results from computing pi on a mesh that
+// is deformed to represent a quarter of a ring
+
+#include "../tests.h"
+#include <base/quadrature_lib.h>
+#include <base/function.h>
+#include <base/numbers.h>
+#include <base/convergence_table.h>
+#include <lac/vector.h>
+#include <grid/tria.h>
+#include <grid/grid_generator.h>
+#include <grid/grid_refinement.h>
+#include <grid/tria_accessor.h>
+#include <grid/tria_iterator.h>
+#include <grid/grid_out.h>
+#include <dofs/dof_handler.h>
+#include <dofs/dof_accessor.h>
+#include <dofs/dof_tools.h>
+#include <fe/fe_values.h>
+#include <fe/mapping.h>
+#include <fe/mapping_q1.h>
+#include <numerics/vectors.h>
+#include <numerics/matrices.h>
+#include <numerics/data_out.h>
+#include <fe/fe_system.h>
+#include <fe/fe_q.h>
+#include <fstream>
+#include <iostream>
+
+#include <fe/mapping_q_eulerian.h>
+
+
+// .... IMPOSED DISPLACEMENT
+
+template <int dim>
+class ImposedDisplacement : public Function<dim>
+{
+  public:
+    ImposedDisplacement() : Function<dim> (dim) { }
+    virtual void vector_value(const Point<dim> &p,
+                              Vector<double> &value) const;
+};
+
+template <>
+void ImposedDisplacement<2>::vector_value(const Point<2> &p,
+                                          Vector<double> &value) const
+{
+  double radius = 1 + (sqrt(5)-1)*p(0);
+  double angle  = 0.5*numbers::PI*(1-p(1));
+  value(0) = radius*sin(angle)-p(0);
+  value(1) = radius*cos(angle)-p(1);
+}
+
+
+// .... MAPPING TEST CLASS
+
+template <int dim>
+class MappingTest
+{
+  public:
+    MappingTest (unsigned int degree);
+    ~MappingTest ();
+
+    void run_test();
+    void graphical_output();
+    
+  private:
+    double compute_area();
+    void explicitly_move_mesh();
+    void write_tria_to_eps(std::string id);
+
+    Triangulation<dim>     triangulation;
+    DoFHandler<dim>        dof_handler;
+    FESystem<dim>          fe;
+
+    unsigned int           degree;
+
+    ImposedDisplacement<dim> imposed_displacement;
+    Vector<double>           displacements;
+};
+
+
+// .... CONSTRUCTOR
+
+template <int dim>
+MappingTest<dim>::MappingTest (unsigned int degree)
+               :
+               dof_handler (triangulation),
+               fe (FE_Q<dim>(degree),dim),
+                degree(degree)
+{ }
+
+
+// .... DESTRUCTOR
+
+template <int dim>
+MappingTest<dim>::~MappingTest () 
+{
+  dof_handler.clear ();
+}
+
+
+// .... COMPUTE AREA
+
+template <int dim>
+double MappingTest<dim>::compute_area () 
+{  
+  QGauss<dim>  quadrature_formula(degree+1);
+
+  MappingQEulerian<dim> mapping(degree,displacements,dof_handler);
+
+  FEValues<dim> fe_values (mapping, fe, quadrature_formula, 
+                          update_JxW_values);
+
+  const unsigned int   n_q_points = quadrature_formula.n_quadrature_points;
+
+  long double area = 0.;
+
+  typename DoFHandler<dim>::active_cell_iterator 
+                              cell = dof_handler.begin_active(),
+                              endc = dof_handler.end();
+
+  for (; cell!=endc; ++cell) {
+    fe_values.reinit (cell);
+    for(unsigned int q=0; q<n_q_points; ++q) area += fe_values.JxW(q);
+  }      
+
+  return area;
+}
+
+
+// .... RUN TEST
+
+template <int dim>
+void MappingTest<dim>::run_test () 
+{
+  GridGenerator::hyper_cube (triangulation,0, 1);
+
+  ConvergenceTable table;
+
+  for(unsigned int ref_level = 0; 
+                   ref_level < 5; 
+                 ++ref_level, triangulation.refine_global(1)){
+
+    dof_handler.distribute_dofs (fe);
+    displacements.reinit (dof_handler.n_dofs());
+
+    VectorTools::interpolate(MappingQ1<dim>(),dof_handler,
+                             imposed_displacement,displacements);
+
+
+    table.add_value("cells",triangulation.n_active_cells());
+    table.add_value("dofs",dof_handler.n_dofs());
+
+    long double area  = compute_area();
+    long double error = std::fabs(numbers::PI-area)/numbers::PI;
+
+    table.add_value("area",  static_cast<double> (area));
+    table.add_value("error", static_cast<double> (error));
+  }    
+
+  table.set_precision("area", 8);
+    table.set_precision("error", 4);
+  table.set_scientific("error", true);
+  table.evaluate_convergence_rates("error",
+    ConvergenceTable::reduction_rate_log2); 
+  table.write_text(deallog.get_file_stream());
+  deallog << std::endl;
+
+}
+
+
+// .... EXPLICITLY MOVE MESH
+
+template <int dim>
+void MappingTest<dim>::explicitly_move_mesh () 
+{  
+  std::vector<bool> moved (triangulation.n_vertices(),false);
+  unsigned int vpc = GeometryInfo<dim>::vertices_per_cell;
+
+  typename DoFHandler<dim>::active_cell_iterator
+    cell = dof_handler.begin_active (),
+    endc = dof_handler.end();
+
+  for (; cell != endc; cell++) {
+  for (unsigned int v=0; v < vpc; v++) {
+  if (moved[cell->vertex_index(v)] == false){
+      moved[cell->vertex_index(v)] =  true;
+      Point<dim> vertex_disp;
+      for (unsigned int d=0; d<dim; d++) {
+        vertex_disp[d] = displacements(cell->vertex_dof_index(v,d));
+      }
+      cell->vertex(v) += vertex_disp;
+  }}}
+}
+
+
+
+// .... GRAPHICAL OUTPUT
+
+template <int dim>
+void MappingTest<dim>::graphical_output () 
+{
+  GridGenerator::hyper_cube (triangulation,0, 1);
+  triangulation.refine_global(4);
+
+  dof_handler.distribute_dofs (fe);
+  displacements.reinit (dof_handler.n_dofs());
+
+  VectorTools::interpolate(MappingQ1<dim>(),dof_handler,
+                           imposed_displacement,displacements);
+
+  explicitly_move_mesh();
+}
+
+
+// .... MAIN
+
+int main () 
+{
+  std::ofstream logfile ("mapping_q_eulerian/output");
+  logfile.precision (2);
+  logfile.setf(std::ios::fixed);  
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-10);
+
+                                  // convergence studies
+
+  for(unsigned int degree = 1; degree <=4; ++degree)
+  {
+    deallog << ".... Q" << degree << " Mapping ...." << std::endl;
+    MappingTest<2> test_one(degree);
+    test_one.run_test();
+  }
+
+                       // graphical output
+
+  MappingTest<2> test_two(1);
+  test_two.graphical_output();
+
+  return 0;
+}
+
diff --git a/tests/fe/mapping_q_eulerian/cmp/generic b/tests/fe/mapping_q_eulerian/cmp/generic
new file mode 100644 (file)
index 0000000..55c69da
--- /dev/null
@@ -0,0 +1,33 @@
+
+DEAL::.... Q1 Mapping ....
+cells dofs    area      error    
+1     8    2.00000000 3.6338e-01 -    
+4     18   2.82842712 9.9684e-02 1.87 
+16    50   3.06146746 2.5505e-02 1.97 
+64    162  3.12144515 6.4131e-03 1.99 
+256   578  3.13654849 1.6056e-03 2.00 
+DEAL::
+DEAL::.... Q2 Mapping ....
+cells dofs    area      error    
+1     18   3.10456950 1.1785e-02 -    
+4     50   3.13914757 7.7829e-04 3.92 
+16    162  3.14143772 4.9318e-05 3.98 
+64    578  3.14158294 3.0930e-06 4.00 
+256   2178 3.14159205 1.9348e-07 4.00 
+DEAL::
+DEAL::.... Q3 Mapping ....
+cells dofs    area      error    
+1     32   3.14653903 1.5745e-03 -    
+4     98   3.14194613 1.1251e-04 3.81 
+16    338  3.14161547 7.2623e-06 3.95 
+64    1250 3.14159409 4.5753e-07 3.99 
+256   4802 3.14159274 2.8653e-08 4.00 
+DEAL::
+DEAL::.... Q4 Mapping ....
+cells dofs    area      error    
+1     50   3.14181857 7.1913e-05 -    
+4     162  3.14159639 1.1900e-06 5.92 
+16    578  3.14159271 1.8860e-08 5.98 
+64    2178 3.14159265 2.9590e-10 5.99 
+256   8450 3.14159265 4.9533e-12 5.90 
+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.