]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Make DoFTools::extract_hanging_node_dofs for p::d/s::Triangualations
authorDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Sat, 9 Dec 2017 01:38:00 +0000 (02:38 +0100)
committerDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Sat, 9 Dec 2017 01:39:24 +0000 (02:39 +0100)
doc/news/changes/minor/20171209DanielArndt [new file with mode: 0644]
include/deal.II/dofs/dof_tools.h
source/dofs/dof_tools.cc
tests/dofs/dof_tools_04a.cc [new file with mode: 0644]
tests/dofs/dof_tools_04a.with_p4est=true.mpirun=1.output [new file with mode: 0644]
tests/dofs/dof_tools_04a.with_p4est=true.mpirun=3.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20171209DanielArndt b/doc/news/changes/minor/20171209DanielArndt
new file mode 100644 (file)
index 0000000..17850e1
--- /dev/null
@@ -0,0 +1,5 @@
+Improved: DoFTools::extract_hanging_node_dofs works for 
+parallel::shared::Triangualtion and parallel::distributed::Triangulation
+and reports DoFs belonging to locally owned cells.
+<br>
+(Daniel Arndt, 2017/12/09)
index 05d5c29632dde97756b2878e36a4ab707786c488..3c88742a294a5ad4f2ef3941cc92ab9f7d34fdac 100644 (file)
@@ -1248,6 +1248,10 @@ namespace DoFTools
    *
    * The size of @p selected_dofs shall equal <tt>dof_handler.n_dofs()</tt>.
    * Previous contents of this array or overwritten.
+   *
+   * In case of a parallel::shared::Triangualtion or a
+   * parallel::distributed::Triangulation only dofs belonging to locally owned
+   * cells are reported.
    */
   template <int dim, int spacedim>
   void
index c115da2d8038d1a321ff51330bd2b066ae9d53ef..f26c89dacd3fc6b495292ed6f2f7b05b0d50d1d9 100644 (file)
@@ -898,23 +898,23 @@ namespace DoFTools
 
         // this function is similar to the make_sparsity_pattern function,
         // see there for more information
-        typename dealii::DoFHandler<dim,spacedim>::active_cell_iterator
-        cell = dof_handler.begin_active(),
-        endc = dof_handler.end();
-        for (; cell!=endc; ++cell)
-          for (unsigned int face=0; face<GeometryInfo<dim>::faces_per_cell; ++face)
-            if (cell->face(face)->has_children())
-              {
-                const typename dealii::DoFHandler<dim,spacedim>::line_iterator
-                line = cell->face(face);
+        for (auto cell : dof_handler.active_cell_iterators())
+          if (cell->is_locally_owned())
+            {
+              for (unsigned int face=0; face<GeometryInfo<dim>::faces_per_cell; ++face)
+                if (cell->face(face)->has_children())
+                  {
+                    const typename dealii::DoFHandler<dim,spacedim>::line_iterator
+                    line = cell->face(face);
 
-                for (unsigned int dof=0; dof!=fe.dofs_per_vertex; ++dof)
-                  selected_dofs[line->child(0)->vertex_dof_index(1,dof)] = true;
+                    for (unsigned int dof=0; dof!=fe.dofs_per_vertex; ++dof)
+                      selected_dofs[line->child(0)->vertex_dof_index(1,dof)] = true;
 
-                for (unsigned int child=0; child<2; ++child)
-                  for (unsigned int dof=0; dof!=fe.dofs_per_line; ++dof)
-                    selected_dofs[line->child(child)->dof_index(dof)] = true;
-              }
+                    for (unsigned int child=0; child<2; ++child)
+                      for (unsigned int dof=0; dof!=fe.dofs_per_line; ++dof)
+                        selected_dofs[line->child(child)->dof_index(dof)] = true;
+                  }
+            }
       }
 
 
@@ -933,49 +933,46 @@ namespace DoFTools
 
         // this function is similar to the make_sparsity_pattern function,
         // see there for more information
+        for (auto cell : dof_handler.active_cell_iterators())
+          if (cell->is_locally_owned())
+            for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
+              if (cell->face(f)->has_children())
+                {
+                  const typename dealii::DoFHandler<dim,spacedim>::face_iterator
+                  face = cell->face(f);
 
-        typename dealii::DoFHandler<dim,spacedim>::active_cell_iterator
-        cell = dof_handler.begin_active(),
-        endc = dof_handler.end();
-        for (; cell!=endc; ++cell)
-          for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
-            if (cell->face(f)->has_children())
-              {
-                const typename dealii::DoFHandler<dim,spacedim>::face_iterator
-                face = cell->face(f);
-
-                for (unsigned int dof=0; dof!=fe.dofs_per_vertex; ++dof)
-                  selected_dofs[face->child(0)->vertex_dof_index(2,dof)] = true;
-
-                // dof numbers on the centers of the lines bounding this
-                // face
-                for (unsigned int line=0; line<4; ++line)
                   for (unsigned int dof=0; dof!=fe.dofs_per_vertex; ++dof)
-                    selected_dofs[face->line(line)->child(0)->vertex_dof_index(1,dof)] = true;
-
-                // next the dofs on the lines interior to the face; the
-                // order of these lines is laid down in the FiniteElement
-                // class documentation
-                for (unsigned int dof=0; dof<fe.dofs_per_line; ++dof)
-                  selected_dofs[face->child(0)->line(1)->dof_index(dof)] = true;
-                for (unsigned int dof=0; dof<fe.dofs_per_line; ++dof)
-                  selected_dofs[face->child(1)->line(2)->dof_index(dof)] = true;
-                for (unsigned int dof=0; dof<fe.dofs_per_line; ++dof)
-                  selected_dofs[face->child(2)->line(3)->dof_index(dof)] = true;
-                for (unsigned int dof=0; dof<fe.dofs_per_line; ++dof)
-                  selected_dofs[face->child(3)->line(0)->dof_index(dof)] = true;
-
-                // dofs on the bordering lines
-                for (unsigned int line=0; line<4; ++line)
-                  for (unsigned int child=0; child<2; ++child)
-                    for (unsigned int dof=0; dof!=fe.dofs_per_line; ++dof)
-                      selected_dofs[face->line(line)->child(child)->dof_index(dof)] = true;
-
-                // finally, for the dofs interior to the four child faces
-                for (unsigned int child=0; child<4; ++child)
-                  for (unsigned int dof=0; dof!=fe.dofs_per_quad; ++dof)
-                    selected_dofs[face->child(child)->dof_index(dof)] = true;
-              }
+                    selected_dofs[face->child(0)->vertex_dof_index(2,dof)] = true;
+
+                  // dof numbers on the centers of the lines bounding this
+                  // face
+                  for (unsigned int line=0; line<4; ++line)
+                    for (unsigned int dof=0; dof!=fe.dofs_per_vertex; ++dof)
+                      selected_dofs[face->line(line)->child(0)->vertex_dof_index(1,dof)] = true;
+
+                  // next the dofs on the lines interior to the face; the
+                  // order of these lines is laid down in the FiniteElement
+                  // class documentation
+                  for (unsigned int dof=0; dof<fe.dofs_per_line; ++dof)
+                    selected_dofs[face->child(0)->line(1)->dof_index(dof)] = true;
+                  for (unsigned int dof=0; dof<fe.dofs_per_line; ++dof)
+                    selected_dofs[face->child(1)->line(2)->dof_index(dof)] = true;
+                  for (unsigned int dof=0; dof<fe.dofs_per_line; ++dof)
+                    selected_dofs[face->child(2)->line(3)->dof_index(dof)] = true;
+                  for (unsigned int dof=0; dof<fe.dofs_per_line; ++dof)
+                    selected_dofs[face->child(3)->line(0)->dof_index(dof)] = true;
+
+                  // dofs on the bordering lines
+                  for (unsigned int line=0; line<4; ++line)
+                    for (unsigned int child=0; child<2; ++child)
+                      for (unsigned int dof=0; dof!=fe.dofs_per_line; ++dof)
+                        selected_dofs[face->line(line)->child(child)->dof_index(dof)] = true;
+
+                  // finally, for the dofs interior to the four child faces
+                  for (unsigned int child=0; child<4; ++child)
+                    for (unsigned int dof=0; dof!=fe.dofs_per_quad; ++dof)
+                      selected_dofs[face->child(child)->dof_index(dof)] = true;
+                }
       }
     }
   }
diff --git a/tests/dofs/dof_tools_04a.cc b/tests/dofs/dof_tools_04a.cc
new file mode 100644 (file)
index 0000000..3c99ea1
--- /dev/null
@@ -0,0 +1,163 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2003 - 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.
+//
+// ---------------------------------------------------------------------
+
+
+// The same as dof_tools_04 but for a parallel::distributed::parallel::distributed::Triangulation
+// instead of a parallel::distributed::Triangulation:
+// check
+//   DoFTools::extract_hanging_node_constraints
+// uses a slightly different refinement and less different FiniteElements
+
+#include "../tests.h"
+#include <deal.II/base/logstream.h>
+#include <deal.II/distributed/tria.h>
+#include <deal.II/grid/tria_iterator.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_out.h>
+#include <deal.II/dofs/dof_accessor.h>
+#include <deal.II/dofs/dof_handler.h>
+#include <deal.II/dofs/dof_tools.h>
+#include <deal.II/fe/fe_q.h>
+#include <deal.II/fe/mapping_q_generic.h>
+
+#include <fstream>
+#include <iomanip>
+#include <iomanip>
+#include <string>
+#include <algorithm>
+
+
+
+template <int dim>
+void
+check_this (const DoFHandler<dim> &dof_handler)
+{
+  types::global_dof_index n_dofs = dof_handler.n_dofs();
+
+  std::vector<bool> is_hanging_node_constrained (n_dofs);
+  DoFTools::extract_hanging_node_dofs (dof_handler,
+                                       is_hanging_node_constrained);
+
+  MappingQGeneric<dim> mapping(1);
+  std::map< types::global_dof_index, Point< dim > > support_point_of_dof;
+  DoFTools::map_dofs_to_support_points (mapping, dof_handler, support_point_of_dof);
+
+  std::vector<Point<dim>> constrained_points;
+
+  for (const auto &pair : support_point_of_dof)
+    if (is_hanging_node_constrained[pair.first])
+      constrained_points.push_back(pair.second);
+
+  const int root = 0;
+  const int mylen = constrained_points.size()*dim;
+  const unsigned int n_processes = Utilities::MPI::n_mpi_processes(MPI_COMM_WORLD);
+  const unsigned int my_id = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD);
+  std::vector<int> recvcounts (n_processes);
+
+  MPI_Gather(&mylen, 1, MPI_INT, recvcounts.data(), 1, MPI_INT, 0, MPI_COMM_WORLD);
+
+  int totlen = 0;
+  std::vector<int> displs;
+  std::vector<Point<dim>> all_points;
+
+  if (my_id==0)
+    {
+      displs.resize(n_processes);
+
+      displs[0] = 0;
+      totlen += recvcounts[0];
+
+      for (int i=1; i<n_processes; i++)
+        {
+          totlen += recvcounts[i];
+          displs[i] = displs[i-1] + recvcounts[i-1];
+        }
+
+      all_points.resize(totlen/dim);
+    }
+
+  MPI_Gatherv(constrained_points.data(), mylen, MPI_DOUBLE,
+              all_points.data(), recvcounts.data(), displs.data(), MPI_DOUBLE,
+              0, MPI_COMM_WORLD);
+
+  std::sort(all_points.begin(), all_points.end(),
+            [](const Point<dim> &a, const Point<dim> &b)
+  {
+    for (unsigned int i=0; i<dim; ++i)
+      {
+        if (a(i) < b(i))
+          return true;
+        if (a(i) > b(i))
+          return false;
+      }
+    return false;
+  });
+  all_points.erase(std::unique(all_points.begin(), all_points.end()), all_points.end());
+
+  deallog << all_points.size() << std::endl;
+  for (const auto &point: all_points)
+    deallog << point << std::endl;
+  deallog << std::endl;
+}
+
+
+
+template <int dim>
+void
+check (const FiniteElement<dim> &fe,
+       const std::string        &name)
+{
+  deallog << "Checking " << name
+          << " in " << dim << "d:"
+          << std::endl;
+
+  // create tria and dofhandler
+  // objects. set different boundary
+  // and sub-domain ids
+  parallel::distributed::Triangulation<dim> tria(MPI_COMM_WORLD);
+  GridGenerator::hyper_cube(tria, 0., 1.);
+  tria.refine_global (1);
+  for (unsigned int ref=0; ref<2; ++ref)
+    {
+      for (auto cell : tria.active_cell_iterators())
+        if (cell->is_locally_owned() &&
+            cell->center()(0)<.5 && cell->center()(1)<.5)
+          cell->set_refine_flag();
+      tria.execute_coarsening_and_refinement ();
+    }
+  DoFHandler<dim> dof_handler (tria);
+  dof_handler.distribute_dofs (fe);
+
+  check_this (dof_handler);
+}
+
+#define CHECK(EL,deg,dim)\
+  { FE_ ## EL<dim> EL(deg);   \
+    check(EL, #EL #deg); }
+
+#define CHECK_ALL(EL,deg)\
+  { CHECK(EL,deg,2); \
+    CHECK(EL,deg,3); \
+  }
+
+int main (int argc, char **argv)
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
+  mpi_initlog();
+  CHECK_ALL(Q,1);
+  CHECK_ALL(Q,2);
+
+  return 0;
+}
diff --git a/tests/dofs/dof_tools_04a.with_p4est=true.mpirun=1.output b/tests/dofs/dof_tools_04a.with_p4est=true.mpirun=1.output
new file mode 100644 (file)
index 0000000..6feb20d
--- /dev/null
@@ -0,0 +1,285 @@
+
+DEAL::Checking Q1 in 2d:
+DEAL::4
+DEAL::0.125000 0.500000
+DEAL::0.375000 0.500000
+DEAL::0.500000 0.125000
+DEAL::0.500000 0.375000
+DEAL::
+DEAL::Checking Q1 in 3d:
+DEAL::40
+DEAL::0.00000 0.500000 0.125000
+DEAL::0.00000 0.500000 0.375000
+DEAL::0.00000 0.500000 0.625000
+DEAL::0.00000 0.500000 0.875000
+DEAL::0.125000 0.500000 0.00000
+DEAL::0.125000 0.500000 0.250000
+DEAL::0.125000 0.500000 0.500000
+DEAL::0.125000 0.500000 0.750000
+DEAL::0.125000 0.500000 1.00000
+DEAL::0.250000 0.500000 0.125000
+DEAL::0.250000 0.500000 0.375000
+DEAL::0.250000 0.500000 0.625000
+DEAL::0.250000 0.500000 0.875000
+DEAL::0.375000 0.500000 0.00000
+DEAL::0.375000 0.500000 0.250000
+DEAL::0.375000 0.500000 0.500000
+DEAL::0.375000 0.500000 0.750000
+DEAL::0.375000 0.500000 1.00000
+DEAL::0.500000 0.00000 0.125000
+DEAL::0.500000 0.00000 0.375000
+DEAL::0.500000 0.00000 0.625000
+DEAL::0.500000 0.00000 0.875000
+DEAL::0.500000 0.125000 0.00000
+DEAL::0.500000 0.125000 0.250000
+DEAL::0.500000 0.125000 0.750000
+DEAL::0.500000 0.125000 1.00000
+DEAL::0.500000 0.125000 0.500000
+DEAL::0.500000 0.250000 0.125000
+DEAL::0.500000 0.250000 0.375000
+DEAL::0.500000 0.250000 0.625000
+DEAL::0.500000 0.250000 0.875000
+DEAL::0.500000 0.375000 0.00000
+DEAL::0.500000 0.375000 0.250000
+DEAL::0.500000 0.375000 0.750000
+DEAL::0.500000 0.375000 1.00000
+DEAL::0.500000 0.375000 0.500000
+DEAL::0.500000 0.500000 0.125000
+DEAL::0.500000 0.500000 0.375000
+DEAL::0.500000 0.500000 0.625000
+DEAL::0.500000 0.500000 0.875000
+DEAL::
+DEAL::Checking Q2 in 2d:
+DEAL::12
+DEAL::0.0625000 0.500000
+DEAL::0.125000 0.500000
+DEAL::0.187500 0.500000
+DEAL::0.312500 0.500000
+DEAL::0.375000 0.500000
+DEAL::0.437500 0.500000
+DEAL::0.500000 0.0625000
+DEAL::0.500000 0.125000
+DEAL::0.500000 0.187500
+DEAL::0.500000 0.312500
+DEAL::0.500000 0.375000
+DEAL::0.500000 0.437500
+DEAL::
+DEAL::Checking Q2 in 3d:
+DEAL::216
+DEAL::0.00000 0.500000 0.0625000
+DEAL::0.00000 0.500000 0.125000
+DEAL::0.00000 0.500000 0.187500
+DEAL::0.00000 0.500000 0.312500
+DEAL::0.00000 0.500000 0.375000
+DEAL::0.00000 0.500000 0.437500
+DEAL::0.00000 0.500000 0.562500
+DEAL::0.00000 0.500000 0.625000
+DEAL::0.00000 0.500000 0.687500
+DEAL::0.00000 0.500000 0.812500
+DEAL::0.00000 0.500000 0.875000
+DEAL::0.00000 0.500000 0.937500
+DEAL::0.0625000 0.500000 0.00000
+DEAL::0.0625000 0.500000 0.0625000
+DEAL::0.0625000 0.500000 0.125000
+DEAL::0.0625000 0.500000 0.187500
+DEAL::0.0625000 0.500000 0.250000
+DEAL::0.0625000 0.500000 0.312500
+DEAL::0.0625000 0.500000 0.375000
+DEAL::0.0625000 0.500000 0.437500
+DEAL::0.0625000 0.500000 0.500000
+DEAL::0.0625000 0.500000 0.562500
+DEAL::0.0625000 0.500000 0.625000
+DEAL::0.0625000 0.500000 0.687500
+DEAL::0.0625000 0.500000 0.750000
+DEAL::0.0625000 0.500000 0.812500
+DEAL::0.0625000 0.500000 0.875000
+DEAL::0.0625000 0.500000 0.937500
+DEAL::0.0625000 0.500000 1.00000
+DEAL::0.125000 0.500000 0.00000
+DEAL::0.125000 0.500000 0.250000
+DEAL::0.125000 0.500000 0.500000
+DEAL::0.125000 0.500000 0.750000
+DEAL::0.125000 0.500000 1.00000
+DEAL::0.187500 0.500000 0.00000
+DEAL::0.187500 0.500000 0.0625000
+DEAL::0.187500 0.500000 0.125000
+DEAL::0.187500 0.500000 0.187500
+DEAL::0.187500 0.500000 0.250000
+DEAL::0.187500 0.500000 0.312500
+DEAL::0.187500 0.500000 0.375000
+DEAL::0.187500 0.500000 0.437500
+DEAL::0.187500 0.500000 0.500000
+DEAL::0.187500 0.500000 0.562500
+DEAL::0.187500 0.500000 0.625000
+DEAL::0.187500 0.500000 0.687500
+DEAL::0.187500 0.500000 0.750000
+DEAL::0.187500 0.500000 0.812500
+DEAL::0.187500 0.500000 0.875000
+DEAL::0.187500 0.500000 0.937500
+DEAL::0.187500 0.500000 1.00000
+DEAL::0.250000 0.500000 0.0625000
+DEAL::0.250000 0.500000 0.125000
+DEAL::0.250000 0.500000 0.187500
+DEAL::0.250000 0.500000 0.312500
+DEAL::0.250000 0.500000 0.375000
+DEAL::0.250000 0.500000 0.437500
+DEAL::0.250000 0.500000 0.562500
+DEAL::0.250000 0.500000 0.625000
+DEAL::0.250000 0.500000 0.687500
+DEAL::0.250000 0.500000 0.812500
+DEAL::0.250000 0.500000 0.875000
+DEAL::0.250000 0.500000 0.937500
+DEAL::0.312500 0.500000 0.00000
+DEAL::0.312500 0.500000 0.0625000
+DEAL::0.312500 0.500000 0.125000
+DEAL::0.312500 0.500000 0.187500
+DEAL::0.312500 0.500000 0.250000
+DEAL::0.312500 0.500000 0.312500
+DEAL::0.312500 0.500000 0.375000
+DEAL::0.312500 0.500000 0.437500
+DEAL::0.312500 0.500000 0.500000
+DEAL::0.312500 0.500000 0.562500
+DEAL::0.312500 0.500000 0.625000
+DEAL::0.312500 0.500000 0.687500
+DEAL::0.312500 0.500000 0.750000
+DEAL::0.312500 0.500000 0.812500
+DEAL::0.312500 0.500000 0.875000
+DEAL::0.312500 0.500000 0.937500
+DEAL::0.312500 0.500000 1.00000
+DEAL::0.375000 0.500000 0.00000
+DEAL::0.375000 0.500000 0.250000
+DEAL::0.375000 0.500000 0.500000
+DEAL::0.375000 0.500000 0.750000
+DEAL::0.375000 0.500000 1.00000
+DEAL::0.437500 0.500000 0.00000
+DEAL::0.437500 0.500000 0.0625000
+DEAL::0.437500 0.500000 0.125000
+DEAL::0.437500 0.500000 0.187500
+DEAL::0.437500 0.500000 0.250000
+DEAL::0.437500 0.500000 0.312500
+DEAL::0.437500 0.500000 0.375000
+DEAL::0.437500 0.500000 0.437500
+DEAL::0.437500 0.500000 0.500000
+DEAL::0.437500 0.500000 0.562500
+DEAL::0.437500 0.500000 0.625000
+DEAL::0.437500 0.500000 0.687500
+DEAL::0.437500 0.500000 0.750000
+DEAL::0.437500 0.500000 0.812500
+DEAL::0.437500 0.500000 0.875000
+DEAL::0.437500 0.500000 0.937500
+DEAL::0.437500 0.500000 1.00000
+DEAL::0.500000 0.00000 0.0625000
+DEAL::0.500000 0.00000 0.125000
+DEAL::0.500000 0.00000 0.187500
+DEAL::0.500000 0.00000 0.312500
+DEAL::0.500000 0.00000 0.375000
+DEAL::0.500000 0.00000 0.437500
+DEAL::0.500000 0.00000 0.562500
+DEAL::0.500000 0.00000 0.625000
+DEAL::0.500000 0.00000 0.687500
+DEAL::0.500000 0.00000 0.812500
+DEAL::0.500000 0.00000 0.875000
+DEAL::0.500000 0.00000 0.937500
+DEAL::0.500000 0.0625000 0.00000
+DEAL::0.500000 0.0625000 0.0625000
+DEAL::0.500000 0.0625000 0.187500
+DEAL::0.500000 0.0625000 0.250000
+DEAL::0.500000 0.0625000 0.312500
+DEAL::0.500000 0.0625000 0.437500
+DEAL::0.500000 0.0625000 0.562500
+DEAL::0.500000 0.0625000 0.687500
+DEAL::0.500000 0.0625000 0.750000
+DEAL::0.500000 0.0625000 0.812500
+DEAL::0.500000 0.0625000 0.937500
+DEAL::0.500000 0.0625000 1.00000
+DEAL::0.500000 0.0625000 0.500000
+DEAL::0.500000 0.125000 0.00000
+DEAL::0.500000 0.125000 0.0625000
+DEAL::0.500000 0.125000 0.187500
+DEAL::0.500000 0.125000 0.250000
+DEAL::0.500000 0.125000 0.312500
+DEAL::0.500000 0.125000 0.437500
+DEAL::0.500000 0.125000 0.562500
+DEAL::0.500000 0.125000 0.687500
+DEAL::0.500000 0.125000 0.750000
+DEAL::0.500000 0.125000 0.812500
+DEAL::0.500000 0.125000 0.937500
+DEAL::0.500000 0.125000 1.00000
+DEAL::0.500000 0.125000 0.500000
+DEAL::0.500000 0.187500 0.00000
+DEAL::0.500000 0.187500 0.0625000
+DEAL::0.500000 0.187500 0.187500
+DEAL::0.500000 0.187500 0.250000
+DEAL::0.500000 0.187500 0.312500
+DEAL::0.500000 0.187500 0.437500
+DEAL::0.500000 0.187500 0.687500
+DEAL::0.500000 0.187500 0.750000
+DEAL::0.500000 0.187500 0.812500
+DEAL::0.500000 0.187500 0.937500
+DEAL::0.500000 0.187500 1.00000
+DEAL::0.500000 0.187500 0.562500
+DEAL::0.500000 0.187500 0.500000
+DEAL::0.500000 0.250000 0.0625000
+DEAL::0.500000 0.250000 0.125000
+DEAL::0.500000 0.250000 0.187500
+DEAL::0.500000 0.250000 0.312500
+DEAL::0.500000 0.250000 0.375000
+DEAL::0.500000 0.250000 0.437500
+DEAL::0.500000 0.250000 0.562500
+DEAL::0.500000 0.250000 0.625000
+DEAL::0.500000 0.250000 0.687500
+DEAL::0.500000 0.250000 0.812500
+DEAL::0.500000 0.250000 0.875000
+DEAL::0.500000 0.250000 0.937500
+DEAL::0.500000 0.312500 0.00000
+DEAL::0.500000 0.312500 0.0625000
+DEAL::0.500000 0.312500 0.187500
+DEAL::0.500000 0.312500 0.250000
+DEAL::0.500000 0.312500 0.312500
+DEAL::0.500000 0.312500 0.437500
+DEAL::0.500000 0.312500 0.687500
+DEAL::0.500000 0.312500 0.750000
+DEAL::0.500000 0.312500 0.812500
+DEAL::0.500000 0.312500 0.937500
+DEAL::0.500000 0.312500 1.00000
+DEAL::0.500000 0.312500 0.562500
+DEAL::0.500000 0.312500 0.500000
+DEAL::0.500000 0.375000 0.00000
+DEAL::0.500000 0.375000 0.0625000
+DEAL::0.500000 0.375000 0.187500
+DEAL::0.500000 0.375000 0.250000
+DEAL::0.500000 0.375000 0.312500
+DEAL::0.500000 0.375000 0.687500
+DEAL::0.500000 0.375000 0.750000
+DEAL::0.500000 0.375000 0.812500
+DEAL::0.500000 0.375000 0.937500
+DEAL::0.500000 0.375000 1.00000
+DEAL::0.500000 0.375000 0.437500
+DEAL::0.500000 0.375000 0.562500
+DEAL::0.500000 0.375000 0.500000
+DEAL::0.500000 0.437500 0.00000
+DEAL::0.500000 0.437500 0.0625000
+DEAL::0.500000 0.437500 0.187500
+DEAL::0.500000 0.437500 0.250000
+DEAL::0.500000 0.437500 0.312500
+DEAL::0.500000 0.437500 0.437500
+DEAL::0.500000 0.437500 0.687500
+DEAL::0.500000 0.437500 0.750000
+DEAL::0.500000 0.437500 0.812500
+DEAL::0.500000 0.437500 0.937500
+DEAL::0.500000 0.437500 1.00000
+DEAL::0.500000 0.437500 0.562500
+DEAL::0.500000 0.437500 0.500000
+DEAL::0.500000 0.500000 0.0625000
+DEAL::0.500000 0.500000 0.125000
+DEAL::0.500000 0.500000 0.187500
+DEAL::0.500000 0.500000 0.312500
+DEAL::0.500000 0.500000 0.375000
+DEAL::0.500000 0.500000 0.437500
+DEAL::0.500000 0.500000 0.562500
+DEAL::0.500000 0.500000 0.625000
+DEAL::0.500000 0.500000 0.687500
+DEAL::0.500000 0.500000 0.812500
+DEAL::0.500000 0.500000 0.875000
+DEAL::0.500000 0.500000 0.937500
+DEAL::
diff --git a/tests/dofs/dof_tools_04a.with_p4est=true.mpirun=3.output b/tests/dofs/dof_tools_04a.with_p4est=true.mpirun=3.output
new file mode 100644 (file)
index 0000000..6feb20d
--- /dev/null
@@ -0,0 +1,285 @@
+
+DEAL::Checking Q1 in 2d:
+DEAL::4
+DEAL::0.125000 0.500000
+DEAL::0.375000 0.500000
+DEAL::0.500000 0.125000
+DEAL::0.500000 0.375000
+DEAL::
+DEAL::Checking Q1 in 3d:
+DEAL::40
+DEAL::0.00000 0.500000 0.125000
+DEAL::0.00000 0.500000 0.375000
+DEAL::0.00000 0.500000 0.625000
+DEAL::0.00000 0.500000 0.875000
+DEAL::0.125000 0.500000 0.00000
+DEAL::0.125000 0.500000 0.250000
+DEAL::0.125000 0.500000 0.500000
+DEAL::0.125000 0.500000 0.750000
+DEAL::0.125000 0.500000 1.00000
+DEAL::0.250000 0.500000 0.125000
+DEAL::0.250000 0.500000 0.375000
+DEAL::0.250000 0.500000 0.625000
+DEAL::0.250000 0.500000 0.875000
+DEAL::0.375000 0.500000 0.00000
+DEAL::0.375000 0.500000 0.250000
+DEAL::0.375000 0.500000 0.500000
+DEAL::0.375000 0.500000 0.750000
+DEAL::0.375000 0.500000 1.00000
+DEAL::0.500000 0.00000 0.125000
+DEAL::0.500000 0.00000 0.375000
+DEAL::0.500000 0.00000 0.625000
+DEAL::0.500000 0.00000 0.875000
+DEAL::0.500000 0.125000 0.00000
+DEAL::0.500000 0.125000 0.250000
+DEAL::0.500000 0.125000 0.750000
+DEAL::0.500000 0.125000 1.00000
+DEAL::0.500000 0.125000 0.500000
+DEAL::0.500000 0.250000 0.125000
+DEAL::0.500000 0.250000 0.375000
+DEAL::0.500000 0.250000 0.625000
+DEAL::0.500000 0.250000 0.875000
+DEAL::0.500000 0.375000 0.00000
+DEAL::0.500000 0.375000 0.250000
+DEAL::0.500000 0.375000 0.750000
+DEAL::0.500000 0.375000 1.00000
+DEAL::0.500000 0.375000 0.500000
+DEAL::0.500000 0.500000 0.125000
+DEAL::0.500000 0.500000 0.375000
+DEAL::0.500000 0.500000 0.625000
+DEAL::0.500000 0.500000 0.875000
+DEAL::
+DEAL::Checking Q2 in 2d:
+DEAL::12
+DEAL::0.0625000 0.500000
+DEAL::0.125000 0.500000
+DEAL::0.187500 0.500000
+DEAL::0.312500 0.500000
+DEAL::0.375000 0.500000
+DEAL::0.437500 0.500000
+DEAL::0.500000 0.0625000
+DEAL::0.500000 0.125000
+DEAL::0.500000 0.187500
+DEAL::0.500000 0.312500
+DEAL::0.500000 0.375000
+DEAL::0.500000 0.437500
+DEAL::
+DEAL::Checking Q2 in 3d:
+DEAL::216
+DEAL::0.00000 0.500000 0.0625000
+DEAL::0.00000 0.500000 0.125000
+DEAL::0.00000 0.500000 0.187500
+DEAL::0.00000 0.500000 0.312500
+DEAL::0.00000 0.500000 0.375000
+DEAL::0.00000 0.500000 0.437500
+DEAL::0.00000 0.500000 0.562500
+DEAL::0.00000 0.500000 0.625000
+DEAL::0.00000 0.500000 0.687500
+DEAL::0.00000 0.500000 0.812500
+DEAL::0.00000 0.500000 0.875000
+DEAL::0.00000 0.500000 0.937500
+DEAL::0.0625000 0.500000 0.00000
+DEAL::0.0625000 0.500000 0.0625000
+DEAL::0.0625000 0.500000 0.125000
+DEAL::0.0625000 0.500000 0.187500
+DEAL::0.0625000 0.500000 0.250000
+DEAL::0.0625000 0.500000 0.312500
+DEAL::0.0625000 0.500000 0.375000
+DEAL::0.0625000 0.500000 0.437500
+DEAL::0.0625000 0.500000 0.500000
+DEAL::0.0625000 0.500000 0.562500
+DEAL::0.0625000 0.500000 0.625000
+DEAL::0.0625000 0.500000 0.687500
+DEAL::0.0625000 0.500000 0.750000
+DEAL::0.0625000 0.500000 0.812500
+DEAL::0.0625000 0.500000 0.875000
+DEAL::0.0625000 0.500000 0.937500
+DEAL::0.0625000 0.500000 1.00000
+DEAL::0.125000 0.500000 0.00000
+DEAL::0.125000 0.500000 0.250000
+DEAL::0.125000 0.500000 0.500000
+DEAL::0.125000 0.500000 0.750000
+DEAL::0.125000 0.500000 1.00000
+DEAL::0.187500 0.500000 0.00000
+DEAL::0.187500 0.500000 0.0625000
+DEAL::0.187500 0.500000 0.125000
+DEAL::0.187500 0.500000 0.187500
+DEAL::0.187500 0.500000 0.250000
+DEAL::0.187500 0.500000 0.312500
+DEAL::0.187500 0.500000 0.375000
+DEAL::0.187500 0.500000 0.437500
+DEAL::0.187500 0.500000 0.500000
+DEAL::0.187500 0.500000 0.562500
+DEAL::0.187500 0.500000 0.625000
+DEAL::0.187500 0.500000 0.687500
+DEAL::0.187500 0.500000 0.750000
+DEAL::0.187500 0.500000 0.812500
+DEAL::0.187500 0.500000 0.875000
+DEAL::0.187500 0.500000 0.937500
+DEAL::0.187500 0.500000 1.00000
+DEAL::0.250000 0.500000 0.0625000
+DEAL::0.250000 0.500000 0.125000
+DEAL::0.250000 0.500000 0.187500
+DEAL::0.250000 0.500000 0.312500
+DEAL::0.250000 0.500000 0.375000
+DEAL::0.250000 0.500000 0.437500
+DEAL::0.250000 0.500000 0.562500
+DEAL::0.250000 0.500000 0.625000
+DEAL::0.250000 0.500000 0.687500
+DEAL::0.250000 0.500000 0.812500
+DEAL::0.250000 0.500000 0.875000
+DEAL::0.250000 0.500000 0.937500
+DEAL::0.312500 0.500000 0.00000
+DEAL::0.312500 0.500000 0.0625000
+DEAL::0.312500 0.500000 0.125000
+DEAL::0.312500 0.500000 0.187500
+DEAL::0.312500 0.500000 0.250000
+DEAL::0.312500 0.500000 0.312500
+DEAL::0.312500 0.500000 0.375000
+DEAL::0.312500 0.500000 0.437500
+DEAL::0.312500 0.500000 0.500000
+DEAL::0.312500 0.500000 0.562500
+DEAL::0.312500 0.500000 0.625000
+DEAL::0.312500 0.500000 0.687500
+DEAL::0.312500 0.500000 0.750000
+DEAL::0.312500 0.500000 0.812500
+DEAL::0.312500 0.500000 0.875000
+DEAL::0.312500 0.500000 0.937500
+DEAL::0.312500 0.500000 1.00000
+DEAL::0.375000 0.500000 0.00000
+DEAL::0.375000 0.500000 0.250000
+DEAL::0.375000 0.500000 0.500000
+DEAL::0.375000 0.500000 0.750000
+DEAL::0.375000 0.500000 1.00000
+DEAL::0.437500 0.500000 0.00000
+DEAL::0.437500 0.500000 0.0625000
+DEAL::0.437500 0.500000 0.125000
+DEAL::0.437500 0.500000 0.187500
+DEAL::0.437500 0.500000 0.250000
+DEAL::0.437500 0.500000 0.312500
+DEAL::0.437500 0.500000 0.375000
+DEAL::0.437500 0.500000 0.437500
+DEAL::0.437500 0.500000 0.500000
+DEAL::0.437500 0.500000 0.562500
+DEAL::0.437500 0.500000 0.625000
+DEAL::0.437500 0.500000 0.687500
+DEAL::0.437500 0.500000 0.750000
+DEAL::0.437500 0.500000 0.812500
+DEAL::0.437500 0.500000 0.875000
+DEAL::0.437500 0.500000 0.937500
+DEAL::0.437500 0.500000 1.00000
+DEAL::0.500000 0.00000 0.0625000
+DEAL::0.500000 0.00000 0.125000
+DEAL::0.500000 0.00000 0.187500
+DEAL::0.500000 0.00000 0.312500
+DEAL::0.500000 0.00000 0.375000
+DEAL::0.500000 0.00000 0.437500
+DEAL::0.500000 0.00000 0.562500
+DEAL::0.500000 0.00000 0.625000
+DEAL::0.500000 0.00000 0.687500
+DEAL::0.500000 0.00000 0.812500
+DEAL::0.500000 0.00000 0.875000
+DEAL::0.500000 0.00000 0.937500
+DEAL::0.500000 0.0625000 0.00000
+DEAL::0.500000 0.0625000 0.0625000
+DEAL::0.500000 0.0625000 0.187500
+DEAL::0.500000 0.0625000 0.250000
+DEAL::0.500000 0.0625000 0.312500
+DEAL::0.500000 0.0625000 0.437500
+DEAL::0.500000 0.0625000 0.562500
+DEAL::0.500000 0.0625000 0.687500
+DEAL::0.500000 0.0625000 0.750000
+DEAL::0.500000 0.0625000 0.812500
+DEAL::0.500000 0.0625000 0.937500
+DEAL::0.500000 0.0625000 1.00000
+DEAL::0.500000 0.0625000 0.500000
+DEAL::0.500000 0.125000 0.00000
+DEAL::0.500000 0.125000 0.0625000
+DEAL::0.500000 0.125000 0.187500
+DEAL::0.500000 0.125000 0.250000
+DEAL::0.500000 0.125000 0.312500
+DEAL::0.500000 0.125000 0.437500
+DEAL::0.500000 0.125000 0.562500
+DEAL::0.500000 0.125000 0.687500
+DEAL::0.500000 0.125000 0.750000
+DEAL::0.500000 0.125000 0.812500
+DEAL::0.500000 0.125000 0.937500
+DEAL::0.500000 0.125000 1.00000
+DEAL::0.500000 0.125000 0.500000
+DEAL::0.500000 0.187500 0.00000
+DEAL::0.500000 0.187500 0.0625000
+DEAL::0.500000 0.187500 0.187500
+DEAL::0.500000 0.187500 0.250000
+DEAL::0.500000 0.187500 0.312500
+DEAL::0.500000 0.187500 0.437500
+DEAL::0.500000 0.187500 0.687500
+DEAL::0.500000 0.187500 0.750000
+DEAL::0.500000 0.187500 0.812500
+DEAL::0.500000 0.187500 0.937500
+DEAL::0.500000 0.187500 1.00000
+DEAL::0.500000 0.187500 0.562500
+DEAL::0.500000 0.187500 0.500000
+DEAL::0.500000 0.250000 0.0625000
+DEAL::0.500000 0.250000 0.125000
+DEAL::0.500000 0.250000 0.187500
+DEAL::0.500000 0.250000 0.312500
+DEAL::0.500000 0.250000 0.375000
+DEAL::0.500000 0.250000 0.437500
+DEAL::0.500000 0.250000 0.562500
+DEAL::0.500000 0.250000 0.625000
+DEAL::0.500000 0.250000 0.687500
+DEAL::0.500000 0.250000 0.812500
+DEAL::0.500000 0.250000 0.875000
+DEAL::0.500000 0.250000 0.937500
+DEAL::0.500000 0.312500 0.00000
+DEAL::0.500000 0.312500 0.0625000
+DEAL::0.500000 0.312500 0.187500
+DEAL::0.500000 0.312500 0.250000
+DEAL::0.500000 0.312500 0.312500
+DEAL::0.500000 0.312500 0.437500
+DEAL::0.500000 0.312500 0.687500
+DEAL::0.500000 0.312500 0.750000
+DEAL::0.500000 0.312500 0.812500
+DEAL::0.500000 0.312500 0.937500
+DEAL::0.500000 0.312500 1.00000
+DEAL::0.500000 0.312500 0.562500
+DEAL::0.500000 0.312500 0.500000
+DEAL::0.500000 0.375000 0.00000
+DEAL::0.500000 0.375000 0.0625000
+DEAL::0.500000 0.375000 0.187500
+DEAL::0.500000 0.375000 0.250000
+DEAL::0.500000 0.375000 0.312500
+DEAL::0.500000 0.375000 0.687500
+DEAL::0.500000 0.375000 0.750000
+DEAL::0.500000 0.375000 0.812500
+DEAL::0.500000 0.375000 0.937500
+DEAL::0.500000 0.375000 1.00000
+DEAL::0.500000 0.375000 0.437500
+DEAL::0.500000 0.375000 0.562500
+DEAL::0.500000 0.375000 0.500000
+DEAL::0.500000 0.437500 0.00000
+DEAL::0.500000 0.437500 0.0625000
+DEAL::0.500000 0.437500 0.187500
+DEAL::0.500000 0.437500 0.250000
+DEAL::0.500000 0.437500 0.312500
+DEAL::0.500000 0.437500 0.437500
+DEAL::0.500000 0.437500 0.687500
+DEAL::0.500000 0.437500 0.750000
+DEAL::0.500000 0.437500 0.812500
+DEAL::0.500000 0.437500 0.937500
+DEAL::0.500000 0.437500 1.00000
+DEAL::0.500000 0.437500 0.562500
+DEAL::0.500000 0.437500 0.500000
+DEAL::0.500000 0.500000 0.0625000
+DEAL::0.500000 0.500000 0.125000
+DEAL::0.500000 0.500000 0.187500
+DEAL::0.500000 0.500000 0.312500
+DEAL::0.500000 0.500000 0.375000
+DEAL::0.500000 0.500000 0.437500
+DEAL::0.500000 0.500000 0.562500
+DEAL::0.500000 0.500000 0.625000
+DEAL::0.500000 0.500000 0.687500
+DEAL::0.500000 0.500000 0.812500
+DEAL::0.500000 0.500000 0.875000
+DEAL::0.500000 0.500000 0.937500
+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.