// find out which is the most dominating finite
// element of the ones that are used on this vertex
- const unsigned int most_dominating_fe_index =
+ unsigned int most_dominating_fe_index =
dof_handler.get_fe_collection().find_dominating_fe(
fe_indices,
/*codim=*/dim);
- // if we found the most dominating element, then use
- // this to eliminate some of the degrees of freedom
- // by identification. otherwise, the code that
- // computes hanging node constraints will have to
- // deal with it by computing appropriate constraints
- // along this face/edge
- if (most_dominating_fe_index !=
+ // if we haven't found a dominating finite element,
+ // choose the very first one to be dominant similar
+ // to compute_vertex_dof_identities()
+ if (most_dominating_fe_index ==
numbers::invalid_unsigned_int)
- {
- // loop over the indices of all the finite
- // elements that are not dominating, and
- // identify their dofs to the most dominating
- // one
- for (const auto &other_fe_index : fe_indices)
- if (other_fe_index != most_dominating_fe_index)
+ most_dominating_fe_index =
+ dealii::internal::DoFAccessorImplementation::
+ Implementation::nth_active_fe_index(
+ dof_handler,
+ 0,
+ vertex_index,
+ 0,
+ std::integral_constant<int, 0>());
+
+ // loop over the indices of all the finite
+ // elements that are not dominating, and
+ // identify their dofs to the most dominating
+ // one
+ for (const auto &other_fe_index : fe_indices)
+ if (other_fe_index != most_dominating_fe_index)
+ {
+ // make sure the entry in the equivalence
+ // table exists
+ ensure_existence_of_dof_identities<0>(
+ dof_handler.get_fe(most_dominating_fe_index),
+ dof_handler.get_fe(other_fe_index),
+ vertex_dof_identities[most_dominating_fe_index]
+ [other_fe_index]);
+
+ // then loop through the identities we
+ // have. first get the global numbers of the
+ // dofs we want to identify and make sure they
+ // are not yet constrained to anything else,
+ // except for to each other. use the rule that
+ // we will always constrain the dof with the
+ // higher fe index to the one with the lower,
+ // to avoid circular reasoning.
+ DoFIdentities &identities =
+ *vertex_dof_identities[most_dominating_fe_index]
+ [other_fe_index];
+ for (const auto &identity : identities)
{
- // make sure the entry in the equivalence
- // table exists
- ensure_existence_of_dof_identities<0>(
- dof_handler.get_fe(most_dominating_fe_index),
- dof_handler.get_fe(other_fe_index),
- vertex_dof_identities[most_dominating_fe_index]
- [other_fe_index]);
+ const types::global_dof_index master_dof_index =
+ dealii::internal::DoFAccessorImplementation::
+ Implementation::get_dof_index(
+ dof_handler,
+ 0,
+ vertex_index,
+ most_dominating_fe_index,
+ identity.first,
+ std::integral_constant<int, 0>());
+ const types::global_dof_index slave_dof_index =
+ dealii::internal::DoFAccessorImplementation::
+ Implementation::get_dof_index(
+ dof_handler,
+ 0,
+ vertex_index,
+ other_fe_index,
+ identity.second,
+ std::integral_constant<int, 0>());
- // then loop through the identities we
- // have. first get the global numbers of the
- // dofs we want to identify and make sure they
- // are not yet constrained to anything else,
- // except for to each other. use the rule that
- // we will always constrain the dof with the
- // higher fe index to the one with the lower,
- // to avoid circular reasoning.
- DoFIdentities &identities =
- *vertex_dof_identities[most_dominating_fe_index]
- [other_fe_index];
- for (const auto &identity : identities)
- {
- const types::global_dof_index
- master_dof_index = dealii::internal::
- DoFAccessorImplementation::
- Implementation::get_dof_index(
- dof_handler,
- 0,
- vertex_index,
- most_dominating_fe_index,
- identity.first,
- std::integral_constant<int, 0>());
- const types::global_dof_index
- slave_dof_index = dealii::internal::
- DoFAccessorImplementation::
- Implementation::get_dof_index(
- dof_handler,
- 0,
- vertex_index,
- other_fe_index,
- identity.second,
- std::integral_constant<int, 0>());
-
- // check if we are on an interface between
- // a locally owned and a ghost cell on which
- // we need to work on.
- //
- // all degrees of freedom belonging to
- // dominating fe indices or to a processor
- // with a higher rank have been set at this
- // point (either in Phase 2, or after the
- // first ghost exchange in Phase 5). thus,
- // we only have to set the indices of
- // degrees of freedom that have been
- // previously flagged invalid.
- if ((slave_dof_index ==
- numbers::invalid_dof_index) &&
- (master_dof_index !=
- numbers::invalid_dof_index))
- dealii::internal::
- DoFAccessorImplementation::
- Implementation::set_dof_index(
- dof_handler,
- 0,
- vertex_index,
- other_fe_index,
- identity.second,
- std::integral_constant<int, 0>(),
- master_dof_index);
- }
+ // check if we are on an interface between
+ // a locally owned and a ghost cell on which
+ // we need to work on.
+ //
+ // all degrees of freedom belonging to
+ // dominating fe indices or to a processor
+ // with a higher rank have been set at this
+ // point (either in Phase 2, or after the
+ // first ghost exchange in Phase 5). thus,
+ // we only have to set the indices of
+ // degrees of freedom that have been
+ // previously flagged invalid.
+ if ((slave_dof_index ==
+ numbers::invalid_dof_index) &&
+ (master_dof_index !=
+ numbers::invalid_dof_index))
+ dealii::internal::DoFAccessorImplementation::
+ Implementation::set_dof_index(
+ dof_handler,
+ 0,
+ vertex_index,
+ other_fe_index,
+ identity.second,
+ std::integral_constant<int, 0>(),
+ master_dof_index);
}
- }
+ }
}
}
}
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2020 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.
+//
+// ---------------------------------------------------------------------
+
+
+
+// have a 2x1 coarse mesh (or 2x1x1) and verify DoF indices in the hp
+// case with a FECollection that contains two finite elements that do
+// not dominate each other. Here, a (FE_Q(1) x FE_Q(2)) and a
+// (FE_Q(2) x FE_Q(1)) element on two separate subdomains face each
+// other. the hp code will unify DoF indices on boundaries between all
+// subdomains.
+
+
+#include <deal.II/distributed/tria.h>
+
+#include <deal.II/fe/fe_q.h>
+#include <deal.II/fe/fe_system.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, Triangulation<dim>::limit_level_difference_at_vertices);
+
+ std::vector<unsigned int> reps(dim, 1U);
+ reps[0] = 2;
+ Point<dim> top_right;
+ for (unsigned int d = 0; d < dim; ++d)
+ top_right[d] = (d == 0 ? 2 : 1);
+ GridGenerator::subdivided_hyper_rectangle(triangulation,
+ reps,
+ Point<dim>(),
+ top_right);
+ Assert(triangulation.n_global_active_cells() == 2, ExcInternalError());
+ Assert(triangulation.n_active_cells() == 2, ExcInternalError());
+
+ hp::FECollection<dim> fe;
+ fe.push_back(FESystem<dim>(FE_Q<dim>(1), 1, FE_Q<dim>(2), 1));
+ fe.push_back(FESystem<dim>(FE_Q<dim>(2), 1, FE_Q<dim>(1), 1));
+
+ hp::DoFHandler<dim> dof_handler(triangulation);
+ for (auto &cell : dof_handler.active_cell_iterators())
+ {
+ if (cell->is_locally_owned())
+ {
+ if (cell->id().to_string() == "0_0:")
+ cell->set_active_fe_index(0);
+ if (cell->id().to_string() == "1_0:")
+ cell->set_active_fe_index(1);
+ }
+ }
+ dof_handler.distribute_dofs(fe);
+
+ deallog << "Processor: " << Utilities::MPI::this_mpi_process(MPI_COMM_WORLD)
+ << std::endl;
+ for (auto &cell : dof_handler.active_cell_iterators())
+ {
+ deallog << " Cell: " << cell;
+ if (cell->is_locally_owned())
+ deallog << " is locally owned";
+ else if (cell->is_ghost())
+ deallog << " is ghost";
+ deallog << std::endl;
+
+ std::vector<types::global_dof_index> dof_indices(
+ cell->get_fe().dofs_per_cell);
+ cell->get_dof_indices(dof_indices);
+ deallog << " ";
+ for (auto i : dof_indices)
+ deallog << i << ' ';
+ deallog << std::endl;
+ }
+ deallog << " n_locally_owned_dofs: " << dof_handler.n_locally_owned_dofs()
+ << std::endl;
+ deallog << " n_global_dofs: " << dof_handler.n_dofs() << std::endl;
+}
+
+
+int
+main(int argc, char *argv[])
+{
+ Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
+ MPILogInitAll log;
+
+ deallog.push("2d");
+ test<2>();
+ deallog.pop();
+
+ deallog.push("3d");
+ test<3>();
+ deallog.pop();
+}
--- /dev/null
+
+DEAL:0:2d::Processor: 0
+DEAL:0:2d:: Cell: 0.0 is locally owned
+DEAL:0:2d:: 0 1 2 3 4 5 6 7 8 9 10 11 12
+DEAL:0:2d:: Cell: 0.1 is ghost
+DEAL:0:2d:: 2 3 13 14 6 7 15 16 17 18 19 20 21
+DEAL:0:2d:: n_locally_owned_dofs: 13
+DEAL:0:2d:: n_global_dofs: 22
+DEAL:0:3d::Processor: 0
+DEAL:0:3d:: Cell: 0.0 is locally owned
+DEAL:0:3d:: 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34
+DEAL:0:3d:: Cell: 0.1 is ghost
+DEAL:0:3d:: 2 3 35 36 6 7 37 38 10 11 39 40 14 15 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61
+DEAL:0:3d:: n_locally_owned_dofs: 35
+DEAL:0:3d:: n_global_dofs: 62
+
+DEAL:1:2d::Processor: 1
+DEAL:1:2d:: Cell: 0.0 is ghost
+DEAL:1:2d:: 0 1 2 3 4 5 6 7 8 9 10 11 12
+DEAL:1:2d:: Cell: 0.1 is locally owned
+DEAL:1:2d:: 2 3 13 14 6 7 15 16 17 18 19 20 21
+DEAL:1:2d:: n_locally_owned_dofs: 9
+DEAL:1:2d:: n_global_dofs: 22
+DEAL:1:3d::Processor: 1
+DEAL:1:3d:: Cell: 0.0 is ghost
+DEAL:1:3d:: 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34
+DEAL:1:3d:: Cell: 0.1 is locally owned
+DEAL:1:3d:: 2 3 35 36 6 7 37 38 10 11 39 40 14 15 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61
+DEAL:1:3d:: n_locally_owned_dofs: 27
+DEAL:1:3d:: n_global_dofs: 62
+
projects / dealii.git / commitdiff