]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Fix DoFHandler::n_boundary_dofs() in issue #8897 8901/head
authorDoug Shi-Dong <doug.shidong@gmail.com>
Fri, 11 Oct 2019 19:19:10 +0000 (15:19 -0400)
committerDoug Shi-Dong <doug.shidong@gmail.com>
Mon, 14 Oct 2019 15:12:23 +0000 (11:12 -0400)
Only count locally owned dofs in the case of a parallel::distributed::Triangulation.

Use suggested syntax on cell iterator.

Co-Authored-By: Daniel Arndt <arndtd@ornl.gov>
Add test case for n_boundary_dofs in parallel

Fix n_boundary_dofs() for hp version.

Also implement suggested change

    Start 3212: dofs/n_boundary_dofs_01.mpirun=1.debug
1/8 Test # 3212: dofs/n_boundary_dofs_01.mpirun=1.debug .....   Passed    0.50 sec
    Start 3213: dofs/n_boundary_dofs_01.mpirun=1.release
2/8 Test # 3213: dofs/n_boundary_dofs_01.mpirun=1.release ...   Passed    0.49 sec
    Start 3214: dofs/n_boundary_dofs_01.mpirun=2.debug
3/8 Test # 3214: dofs/n_boundary_dofs_01.mpirun=2.debug .....   Passed    0.50 sec
    Start 3215: dofs/n_boundary_dofs_01.mpirun=2.release
4/8 Test # 3215: dofs/n_boundary_dofs_01.mpirun=2.release ...   Passed    0.51 sec
    Start 3216: dofs/n_boundary_dofs_02.debug
5/8 Test # 3216: dofs/n_boundary_dofs_02.debug ..............   Passed    0.51 sec
    Start 3217: dofs/n_boundary_dofs_02.release
6/8 Test # 3217: dofs/n_boundary_dofs_02.release ............   Passed    0.52 sec
    Start 3218: dofs/n_boundary_dofs_03.debug
7/8 Test # 3218: dofs/n_boundary_dofs_03.debug ..............   Passed    0.52 sec
    Start 3219: dofs/n_boundary_dofs_03.release
8/8 Test # 3219: dofs/n_boundary_dofs_03.release ............   Passed    0.54 sec

100% tests passed, 0 tests failed out of 8

Add hp n_boundary_dofs test.

Use unordered_set instead.
Fix a || to a && from switching the if-continue statement.

      Start 3212: dofs/n_boundary_dofs_01.mpirun=1.debug
 1/18 Test # 3212: dofs/n_boundary_dofs_01.mpirun=1.debug .....   Passed   51.27 sec
      Start 3213: dofs/n_boundary_dofs_01.mpirun=1.release
 2/18 Test # 3213: dofs/n_boundary_dofs_01.mpirun=1.release ...   Passed   17.03 sec
      Start 3214: dofs/n_boundary_dofs_01.mpirun=2.debug
 3/18 Test # 3214: dofs/n_boundary_dofs_01.mpirun=2.debug .....   Passed    1.90 sec
      Start 3215: dofs/n_boundary_dofs_01.mpirun=2.release
 4/18 Test # 3215: dofs/n_boundary_dofs_01.mpirun=2.release ...   Passed    1.45 sec
      Start 3216: dofs/n_boundary_dofs_02.debug
 5/18 Test # 3216: dofs/n_boundary_dofs_02.debug ..............   Passed   50.17 sec
      Start 3217: dofs/n_boundary_dofs_02.release
 6/18 Test # 3217: dofs/n_boundary_dofs_02.release ............   Passed   16.30 sec
      Start 3218: dofs/n_boundary_dofs_03.debug
 7/18 Test # 3218: dofs/n_boundary_dofs_03.debug ..............   Passed   50.67 sec
      Start 3219: dofs/n_boundary_dofs_03.release
 8/18 Test # 3219: dofs/n_boundary_dofs_03.release ............   Passed   16.71 sec
      Start 5498: hp/n_boundary_dofs.debug
 9/18 Test # 5498: hp/n_boundary_dofs.debug ...................   Passed   44.68 sec
      Start 5499: hp/n_boundary_dofs.release
10/18 Test # 5499: hp/n_boundary_dofs.release .................   Passed   11.83 sec
      Start 5500: hp/n_boundary_dofs_01.mpirun=1.debug
11/18 Test # 5500: hp/n_boundary_dofs_01.mpirun=1.debug .......   Passed   34.48 sec
      Start 5501: hp/n_boundary_dofs_01.mpirun=1.release
12/18 Test # 5501: hp/n_boundary_dofs_01.mpirun=1.release .....   Passed   12.46 sec
      Start 5502: hp/n_boundary_dofs_01.mpirun=2.debug
13/18 Test # 5502: hp/n_boundary_dofs_01.mpirun=2.debug .......   Passed    1.75 sec
      Start 5503: hp/n_boundary_dofs_01.mpirun=2.release
14/18 Test # 5503: hp/n_boundary_dofs_01.mpirun=2.release .....   Passed    1.36 sec
      Start 5504: hp/n_boundary_dofs_02.debug
15/18 Test # 5504: hp/n_boundary_dofs_02.debug ................   Passed   34.33 sec
      Start 5505: hp/n_boundary_dofs_02.release
16/18 Test # 5505: hp/n_boundary_dofs_02.release ..............   Passed   12.25 sec
      Start 5506: hp/n_boundary_dofs_03.debug
17/18 Test # 5506: hp/n_boundary_dofs_03.debug ................   Passed   41.88 sec
      Start 5507: hp/n_boundary_dofs_03.release
18/18 Test # 5507: hp/n_boundary_dofs_03.release ..............   Passed   15.45 sec

100% tests passed, 0 tests failed out of 18

doc/news/changes/minor/20191011DougShiDong [new file with mode: 0644]
include/deal.II/dofs/dof_handler.h
include/deal.II/hp/dof_handler.h
source/dofs/dof_handler.cc
source/hp/dof_handler.cc
tests/dofs/n_boundary_dofs_01.cc [new file with mode: 0644]
tests/dofs/n_boundary_dofs_01.with_p4est=true.mpirun=1.output [new file with mode: 0644]
tests/dofs/n_boundary_dofs_01.with_p4est=true.mpirun=2.output [new file with mode: 0644]
tests/hp/n_boundary_dofs_01.cc [new file with mode: 0644]
tests/hp/n_boundary_dofs_01.with_p4est=true.mpirun=1.output [new file with mode: 0644]
tests/hp/n_boundary_dofs_01.with_p4est=true.mpirun=2.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20191011DougShiDong b/doc/news/changes/minor/20191011DougShiDong
new file mode 100644 (file)
index 0000000..d7c9956
--- /dev/null
@@ -0,0 +1,7 @@
+Fixed: Take in account p::d::Triangulation when using
+DoFHandler<dim, spacedim>::n_boundary_dofs() and its
+hp version. It now returns locally owned dofs instead
+of attempting a global count.
+<br>
+(Doug Shi-Dong, 2019/10/11)
+
index a29432a77d031a9c60bc7906bf4c14f28a58a63e..fc0a5763d0d08048550bd0449a660e6568184bd3 100644 (file)
@@ -954,15 +954,17 @@ public:
   n_dofs(const unsigned int level) const;
 
   /**
-   * Return the number of degrees of freedom located on the boundary.
+   * Return the number of locally owned degrees of freedom located on the
+   * boundary.
    */
   types::global_dof_index
   n_boundary_dofs() const;
 
   /**
-   * Return the number of degrees of freedom located on those parts of the
-   * boundary which have a boundary indicator listed in the given set. The
-   * reason that a @p map rather than a @p set is used is the same as
+   * Return the number of locally owned degrees of freedom located on those
+   * parts of the boundary which have a boundary indicator listed in the given
+   * set.
+   * The reason that a @p map rather than a @p set is used is the same as
    * described in the documentation of that variant of
    * DoFTools::make_boundary_sparsity_pattern() that takes a map.
    *
index e4850d9fd1b28d2fed08adddbc98c622ecd77687..badc4bdf1a44a845d3d9bef045c054a56f67e9bf 100644 (file)
@@ -778,7 +778,8 @@ namespace hp
     n_dofs(const unsigned int level) const;
 
     /**
-     * Return the number of degrees of freedom located on the boundary.
+     * Return the number of locally owned degrees of freedom located on the
+     * boundary.
      */
     types::global_dof_index
     n_boundary_dofs() const;
@@ -800,8 +801,9 @@ namespace hp
         &boundary_ids) const;
 
     /**
-     * Return the number of degrees of freedom located on those parts of the
-     * boundary which have a boundary indicator listed in the given set. The
+     * Return the number of locally owned degrees of freedom located on those
+     * parts of the boundary which have a boundary indicator listed in the given
+     * set.
      */
     types::global_dof_index
     n_boundary_dofs(const std::set<types::boundary_id> &boundary_ids) const;
index 21bd935fbf0f9313e1657c52bf9be6333fe1fed5..2145bd5b90cef36fa906e781b64ca8c1a531db5a 100644 (file)
@@ -35,6 +35,7 @@
 
 #include <algorithm>
 #include <set>
+#include <unordered_set>
 
 DEAL_II_NAMESPACE_OPEN
 
@@ -1096,36 +1097,51 @@ template <int dim, int spacedim>
 types::global_dof_index
 DoFHandler<dim, spacedim>::n_boundary_dofs() const
 {
-  std::set<int> boundary_dofs;
-
-  const unsigned int                   dofs_per_face = get_fe().dofs_per_face;
-  std::vector<types::global_dof_index> dofs_on_face(dofs_per_face);
-
-  // loop over all faces of all cells
-  // and see whether they are at a
-  // boundary. note (i) that we visit
-  // interior faces twice (which we
-  // don't care about) but exterior
-  // faces only once as is
-  // appropriate, and (ii) that we
-  // need not take special care of
-  // single lines (using
-  // @p{cell->has_boundary_lines}),
-  // since we do not support
-  // boundaries of dimension dim-2,
-  // and so every boundary line is
+  const FiniteElement<dim, spacedim> &fe            = get_fe();
+  const unsigned int                  dofs_per_face = fe.n_dofs_per_face();
+  std::vector<dealii::types::global_dof_index> dofs_on_face(dofs_per_face);
+
+  const IndexSet &owned_dofs = locally_owned_dofs();
+
+  std::unordered_set<unsigned int> boundary_dof_indices;
+  // Loop over all faces of all cells and see whether they are at a boundary.
+  // note (i) that we visit interior faces twice (which we don't care about) but
+  // exterior faces only once as is appropriate, and (ii) that we need not take
+  // special care of single lines (using @p{cell->has_boundary_lines}), since we
+  // do not support boundaries of dimension dim-2, and so every boundary line is
   // also part of a boundary face.
-  active_cell_iterator cell = begin_active(), endc = end();
-  for (; cell != endc; ++cell)
-    for (unsigned int f = 0; f < GeometryInfo<dim>::faces_per_cell; ++f)
-      if (cell->at_boundary(f))
+  for (const auto &cell : this->active_cell_iterators())
+    {
+      if (cell->is_locally_owned() && cell->at_boundary())
         {
-          cell->face(f)->get_dof_indices(dofs_on_face);
-          for (unsigned int i = 0; i < dofs_per_face; ++i)
-            boundary_dofs.insert(dofs_on_face[i]);
+          for (unsigned int iface = 0;
+               iface < dealii::GeometryInfo<dim>::faces_per_cell;
+               ++iface)
+            {
+              const auto face = cell->face(iface);
+              if (face->at_boundary())
+                {
+                  cell->face(iface)->get_dof_indices(dofs_on_face);
+                  for (unsigned int idof_face = 0; idof_face < dofs_per_face;
+                       ++idof_face)
+                    {
+                      const unsigned int global_idof_index =
+                        dofs_on_face[idof_face];
+
+                      // If dof is not already in our hashmap, then insert it
+                      // into our set of surface indices. However, do not add if
+                      // it is not locally owned, it will get picked up by
+                      // another processor
+                      if (owned_dofs.is_element(global_idof_index))
+                        {
+                          boundary_dof_indices.insert(global_idof_index);
+                        }
+                    }
+                }
+            }
         }
-
-  return boundary_dofs.size();
+    }
+  return boundary_dof_indices.size();
 }
 
 
@@ -1139,28 +1155,56 @@ DoFHandler<dim, spacedim>::n_boundary_dofs(
            boundary_ids.end(),
          ExcInvalidBoundaryIndicator());
 
-  std::set<types::global_dof_index> boundary_dofs;
-
-  const unsigned int                   dofs_per_face = get_fe().dofs_per_face;
-  std::vector<types::global_dof_index> dofs_on_face(dofs_per_face);
-
-  // same as in the previous
-  // function, but with a different
-  // check for the boundary indicator
-  active_cell_iterator cell = begin_active(), endc = end();
-  for (; cell != endc; ++cell)
-    for (unsigned int f = 0; f < GeometryInfo<dim>::faces_per_cell; ++f)
-      if (cell->at_boundary(f) &&
-          (std::find(boundary_ids.begin(),
-                     boundary_ids.end(),
-                     cell->face(f)->boundary_id()) != boundary_ids.end()))
+  const FiniteElement<dim, spacedim> &fe            = get_fe();
+  const unsigned int                  dofs_per_face = fe.n_dofs_per_face();
+  std::vector<dealii::types::global_dof_index> dofs_on_face(dofs_per_face);
+
+  const IndexSet &owned_dofs = locally_owned_dofs();
+
+  std::unordered_set<unsigned int> boundary_dof_indices;
+  // Loop over all faces of all cells and see whether they are at a boundary.
+  // note (i) that we visit interior faces twice (which we don't care about) but
+  // exterior faces only once as is appropriate, and (ii) that we need not take
+  // special care of single lines (using @p{cell->has_boundary_lines}), since we
+  // do not support boundaries of dimension dim-2, and so every boundary line is
+  // also part of a boundary face.
+  for (const auto &cell : this->active_cell_iterators())
+    {
+      if (cell->is_locally_owned() && cell->at_boundary())
         {
-          cell->face(f)->get_dof_indices(dofs_on_face);
-          for (unsigned int i = 0; i < dofs_per_face; ++i)
-            boundary_dofs.insert(dofs_on_face[i]);
+          for (unsigned int iface = 0;
+               iface < dealii::GeometryInfo<dim>::faces_per_cell;
+               ++iface)
+            {
+              const auto         face        = cell->face(iface);
+              const unsigned int boundary_id = face->boundary_id();
+
+              if (face->at_boundary() &&
+                  std::find(boundary_ids.begin(),
+                            boundary_ids.end(),
+                            boundary_id) != boundary_ids.end())
+                {
+                  cell->face(iface)->get_dof_indices(dofs_on_face);
+                  for (unsigned int idof_face = 0; idof_face < dofs_per_face;
+                       ++idof_face)
+                    {
+                      const unsigned int global_idof_index =
+                        dofs_on_face[idof_face];
+
+                      // If dof is not already in our hashmap, then insert it
+                      // into our set of surface indices. However, do not add if
+                      // it is not locally owned, it will get picked up by
+                      // another processor
+                      if (owned_dofs.is_element(global_idof_index))
+                        {
+                          boundary_dof_indices.insert(global_idof_index);
+                        }
+                    }
+                }
+            }
         }
-
-  return boundary_dofs.size();
+    }
+  return boundary_dof_indices.size();
 }
 
 
index 27a0948ab2c8f9de7b26c42930f886f01db4f594..c5e1678176b17c59454d4b0a0cd1ab417468e70f 100644 (file)
@@ -42,6 +42,7 @@
 #include <algorithm>
 #include <functional>
 #include <set>
+#include <unordered_set>
 
 DEAL_II_NAMESPACE_OPEN
 
@@ -1431,10 +1432,12 @@ namespace hp
   {
     Assert(fe_collection.size() > 0, ExcNoFESelected());
 
-    std::set<types::global_dof_index>    boundary_dofs;
-    std::vector<types::global_dof_index> dofs_on_face;
+    std::unordered_set<types::global_dof_index> boundary_dofs;
+    std::vector<types::global_dof_index>        dofs_on_face;
     dofs_on_face.reserve(this->get_fe_collection().max_dofs_per_face());
 
+    const IndexSet &owned_dofs = locally_owned_dofs();
+
     // loop over all faces to check whether they are at a
     // boundary. note that we need not take special care of single
     // lines in 3d (using @p{cell->has_boundary_lines}), since we do
@@ -1444,17 +1447,26 @@ namespace hp
       cell = this->begin_active(),
       endc = this->end();
     for (; cell != endc; ++cell)
-      for (unsigned int f = 0; f < GeometryInfo<dim>::faces_per_cell; ++f)
-        if (cell->at_boundary(f))
-          {
-            const unsigned int dofs_per_face = cell->get_fe().dofs_per_face;
-            dofs_on_face.resize(dofs_per_face);
+      if (cell->is_locally_owned() && cell->at_boundary())
+        {
+          for (unsigned int f = 0; f < GeometryInfo<dim>::faces_per_cell; ++f)
+            if (cell->at_boundary(f))
+              {
+                const unsigned int dofs_per_face = cell->get_fe().dofs_per_face;
+                dofs_on_face.resize(dofs_per_face);
 
-            cell->face(f)->get_dof_indices(dofs_on_face,
-                                           cell->active_fe_index());
-            for (unsigned int i = 0; i < dofs_per_face; ++i)
-              boundary_dofs.insert(dofs_on_face[i]);
-          }
+                cell->face(f)->get_dof_indices(dofs_on_face,
+                                               cell->active_fe_index());
+                for (unsigned int i = 0; i < dofs_per_face; ++i)
+                  {
+                    const unsigned int global_idof_index = dofs_on_face[i];
+                    if (owned_dofs.is_element(global_idof_index))
+                      {
+                        boundary_dofs.insert(global_idof_index);
+                      }
+                  }
+              }
+        }
     return boundary_dofs.size();
   }
 
@@ -1472,27 +1484,38 @@ namespace hp
 
     // same as above, but with additional checks for set of boundary
     // indicators
-    std::set<types::global_dof_index>    boundary_dofs;
-    std::vector<types::global_dof_index> dofs_on_face;
+    std::unordered_set<types::global_dof_index> boundary_dofs;
+    std::vector<types::global_dof_index>        dofs_on_face;
     dofs_on_face.reserve(this->get_fe_collection().max_dofs_per_face());
 
+    const IndexSet &owned_dofs = locally_owned_dofs();
+
     typename HpDoFHandler<dim, spacedim>::active_cell_iterator
       cell = this->begin_active(),
       endc = this->end();
     for (; cell != endc; ++cell)
-      for (unsigned int f = 0; f < GeometryInfo<dim>::faces_per_cell; ++f)
-        if (cell->at_boundary(f) &&
-            (boundary_ids.find(cell->face(f)->boundary_id()) !=
-             boundary_ids.end()))
-          {
-            const unsigned int dofs_per_face = cell->get_fe().dofs_per_face;
-            dofs_on_face.resize(dofs_per_face);
+      if (cell->is_locally_owned() && cell->at_boundary())
+        {
+          for (unsigned int f = 0; f < GeometryInfo<dim>::faces_per_cell; ++f)
+            if (cell->at_boundary(f) &&
+                (boundary_ids.find(cell->face(f)->boundary_id()) !=
+                 boundary_ids.end()))
+              {
+                const unsigned int dofs_per_face = cell->get_fe().dofs_per_face;
+                dofs_on_face.resize(dofs_per_face);
 
-            cell->face(f)->get_dof_indices(dofs_on_face,
-                                           cell->active_fe_index());
-            for (unsigned int i = 0; i < dofs_per_face; ++i)
-              boundary_dofs.insert(dofs_on_face[i]);
-          }
+                cell->face(f)->get_dof_indices(dofs_on_face,
+                                               cell->active_fe_index());
+                for (unsigned int i = 0; i < dofs_per_face; ++i)
+                  {
+                    const unsigned int global_idof_index = dofs_on_face[i];
+                    if (owned_dofs.is_element(global_idof_index))
+                      {
+                        boundary_dofs.insert(global_idof_index);
+                      }
+                  }
+              }
+        }
     return boundary_dofs.size();
   }
 
diff --git a/tests/dofs/n_boundary_dofs_01.cc b/tests/dofs/n_boundary_dofs_01.cc
new file mode 100644 (file)
index 0000000..80b12ed
--- /dev/null
@@ -0,0 +1,65 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2008 - 2018 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.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+
+// test that DoFHandler::n_boundary_dofs() yields the correct
+// results for parallel::distributed::Triangulations
+
+#include <deal.II/distributed/tria.h>
+
+#include <deal.II/dofs/dof_handler.h>
+
+#include <deal.II/fe/fe_q.h>
+
+#include <deal.II/grid/grid_generator.h>
+
+#include "../tests.h"
+
+
+template <int dim>
+void
+test()
+{
+  parallel::distributed::Triangulation<dim> triangulation(MPI_COMM_WORLD);
+  GridGenerator::hyper_cube(triangulation);
+  triangulation.refine_global(3);
+
+  FE_Q<dim>       fe(2);
+  DoFHandler<dim> dof_handler(triangulation);
+
+  dof_handler.distribute_dofs(fe);
+
+  const unsigned int local_boundary_dofs = dof_handler.n_boundary_dofs();
+  const unsigned int global_boundary_dofs =
+    Utilities::MPI::sum(local_boundary_dofs, MPI_COMM_WORLD);
+  deallog << global_boundary_dofs << std::endl;
+}
+
+
+int
+main(int argc, char *argv[])
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
+  mpi_initlog();
+
+  deallog.push("2d");
+  test<2>();
+  deallog.pop();
+
+  deallog.push("3d");
+  test<3>();
+  deallog.pop();
+}
diff --git a/tests/dofs/n_boundary_dofs_01.with_p4est=true.mpirun=1.output b/tests/dofs/n_boundary_dofs_01.with_p4est=true.mpirun=1.output
new file mode 100644 (file)
index 0000000..1108787
--- /dev/null
@@ -0,0 +1,3 @@
+
+DEAL:2d::64
+DEAL:3d::1538
diff --git a/tests/dofs/n_boundary_dofs_01.with_p4est=true.mpirun=2.output b/tests/dofs/n_boundary_dofs_01.with_p4est=true.mpirun=2.output
new file mode 100644 (file)
index 0000000..1108787
--- /dev/null
@@ -0,0 +1,3 @@
+
+DEAL:2d::64
+DEAL:3d::1538
diff --git a/tests/hp/n_boundary_dofs_01.cc b/tests/hp/n_boundary_dofs_01.cc
new file mode 100644 (file)
index 0000000..b7a7b03
--- /dev/null
@@ -0,0 +1,67 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2008 - 2018 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.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+
+// test that DoFHandler::n_boundary_dofs() yields the correct
+// results for parallel::distributed::Triangulations
+
+#include <deal.II/distributed/tria.h>
+
+#include <deal.II/fe/fe_q.h>
+
+#include <deal.II/grid/grid_generator.h>
+
+#include <deal.II/hp/dof_handler.h>
+#include <deal.II/hp/fe_collection.h>
+
+#include "../tests.h"
+
+
+template <int dim>
+void
+test()
+{
+  parallel::distributed::Triangulation<dim> triangulation(MPI_COMM_WORLD);
+  GridGenerator::hyper_cube(triangulation);
+  triangulation.refine_global(3);
+
+  FE_Q<dim>           fe(2);
+  hp::FECollection    fe_coll(fe);
+  hp::DoFHandler<dim> dof_handler(triangulation);
+
+  dof_handler.distribute_dofs(fe_coll);
+
+  const unsigned int local_boundary_dofs = dof_handler.n_boundary_dofs();
+  const unsigned int global_boundary_dofs =
+    Utilities::MPI::sum(local_boundary_dofs, MPI_COMM_WORLD);
+  deallog << global_boundary_dofs << std::endl;
+}
+
+
+int
+main(int argc, char *argv[])
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
+  mpi_initlog();
+
+  deallog.push("2d");
+  test<2>();
+  deallog.pop();
+
+  deallog.push("3d");
+  test<3>();
+  deallog.pop();
+}
diff --git a/tests/hp/n_boundary_dofs_01.with_p4est=true.mpirun=1.output b/tests/hp/n_boundary_dofs_01.with_p4est=true.mpirun=1.output
new file mode 100644 (file)
index 0000000..1108787
--- /dev/null
@@ -0,0 +1,3 @@
+
+DEAL:2d::64
+DEAL:3d::1538
diff --git a/tests/hp/n_boundary_dofs_01.with_p4est=true.mpirun=2.output b/tests/hp/n_boundary_dofs_01.with_p4est=true.mpirun=2.output
new file mode 100644 (file)
index 0000000..1108787
--- /dev/null
@@ -0,0 +1,3 @@
+
+DEAL:2d::64
+DEAL:3d::1538

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.