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
--- /dev/null
+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)
+
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.
*
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;
&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;
#include <algorithm>
#include <set>
+#include <unordered_set>
DEAL_II_NAMESPACE_OPEN
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();
}
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();
}
#include <algorithm>
#include <functional>
#include <set>
+#include <unordered_set>
DEAL_II_NAMESPACE_OPEN
{
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
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();
}
// 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();
}
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// 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();
+}
--- /dev/null
+
+DEAL:2d::64
+DEAL:3d::1538
--- /dev/null
+
+DEAL:2d::64
+DEAL:3d::1538
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// 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();
+}
--- /dev/null
+
+DEAL:2d::64
+DEAL:3d::1538
--- /dev/null
+
+DEAL:2d::64
+DEAL:3d::1538