From 44e94b7dee0d83da825df2e5e9dfb171d92099d6 Mon Sep 17 00:00:00 2001 From: Marc Fehling Date: Sat, 20 Aug 2022 21:03:55 -0600 Subject: [PATCH] Added get/set functions for future FE indices. --- doc/news/changes/minor/20220820Fehling | 6 ++ include/deal.II/dofs/dof_handler.h | 16 ++++ source/dofs/dof_handler.cc | 41 ++++++++ tests/mpi/hp_getset_future_fe_indices_01.cc | 95 +++++++++++++++++++ ...indices_01.with_p4est=true.mpirun=2.output | 7 ++ 5 files changed, 165 insertions(+) create mode 100644 doc/news/changes/minor/20220820Fehling create mode 100644 tests/mpi/hp_getset_future_fe_indices_01.cc create mode 100644 tests/mpi/hp_getset_future_fe_indices_01.with_p4est=true.mpirun=2.output diff --git a/doc/news/changes/minor/20220820Fehling b/doc/news/changes/minor/20220820Fehling new file mode 100644 index 0000000000..d181b8dadc --- /dev/null +++ b/doc/news/changes/minor/20220820Fehling @@ -0,0 +1,6 @@ +New: Functions DoFHandler::get_future_fe_indices() and +DoFHandler::set_future_fe_indices() that store and recover future FE +indices on a DoFHandler. They can also be used to transfer future FE +indices between different DoFHandler objects. +
+(Marc Fehling, 2022/08/20) diff --git a/include/deal.II/dofs/dof_handler.h b/include/deal.II/dofs/dof_handler.h index 2c9417b7be..a38e256b49 100644 --- a/include/deal.II/dofs/dof_handler.h +++ b/include/deal.II/dofs/dof_handler.h @@ -586,6 +586,22 @@ public: void get_active_fe_indices(std::vector &active_fe_indices) const; + /** + * Go through the triangulation and set the future FE indices of all + * locally owned cells to the values given in @p future_fe_indices. + * Cells corresponding to invalid_active_fe_index will be skipped. + */ + void + set_future_fe_indices(const std::vector &future_fe_indices); + + /** + * Go through the triangulation and return a vector of future FE indices of + * all locally owned cells. If no future FE index has been set on a cell, + * its value will be invalid_active_fe_index. + */ + std::vector + get_future_fe_indices() const; + /** * Assign a Triangulation to the DoFHandler. * diff --git a/source/dofs/dof_handler.cc b/source/dofs/dof_handler.cc index bd2b7bc949..9f15f7a2b3 100644 --- a/source/dofs/dof_handler.cc +++ b/source/dofs/dof_handler.cc @@ -2661,6 +2661,47 @@ DoFHandler::get_active_fe_indices( +template +void +DoFHandler::set_future_fe_indices( + const std::vector &future_fe_indices) +{ + Assert(future_fe_indices.size() == this->get_triangulation().n_active_cells(), + ExcDimensionMismatch(future_fe_indices.size(), + this->get_triangulation().n_active_cells())); + + this->create_active_fe_table(); + // we could set the values directly, since they are stored as + // protected data of this object, but for simplicity we use the + // cell-wise access. this way we also have to pass some debug-mode + // tests which we would have to duplicate ourselves otherwise + for (const auto &cell : this->active_cell_iterators()) + if (cell->is_locally_owned() && + future_fe_indices[cell->active_cell_index()] != invalid_active_fe_index) + cell->set_future_fe_index(future_fe_indices[cell->active_cell_index()]); +} + + + +template +std::vector +DoFHandler::get_future_fe_indices() const +{ + std::vector future_fe_indices( + this->get_triangulation().n_active_cells(), invalid_active_fe_index); + + // we could try to extract the values directly, since they are + // stored as protected data of this object, but for simplicity we + // use the cell-wise access. + for (const auto &cell : this->active_cell_iterators()) + if (cell->is_locally_owned() && cell->future_fe_index_set()) + future_fe_indices[cell->active_cell_index()] = cell->future_fe_index(); + + return future_fe_indices; +} + + + template void DoFHandler::connect_to_triangulation_signals() diff --git a/tests/mpi/hp_getset_future_fe_indices_01.cc b/tests/mpi/hp_getset_future_fe_indices_01.cc new file mode 100644 index 0000000000..6ac3b61517 --- /dev/null +++ b/tests/mpi/hp_getset_future_fe_indices_01.cc @@ -0,0 +1,95 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2022 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. +// +// --------------------------------------------------------------------- + + + +// Transfer future FE indices between DoFHandler objects. + + +#include + +#include + +#include + +#include + +#include "../tests.h" + +#include "../test_grids.h" + + +template +void +test() +{ + MPI_Comm mpi_comm = MPI_COMM_WORLD; + + const unsigned int n_procs = Utilities::MPI::n_mpi_processes(mpi_comm); + + parallel::distributed::Triangulation tria(mpi_comm); + TestGrids::hyper_line(tria, n_procs); + tria.refine_global(1); + + DoFHandler dofh_src(tria); + DoFHandler dofh_dst(tria); + + // set future fe index on first cell of every process + for (const auto &cell : dofh_src.active_cell_iterators()) + if (cell->is_locally_owned()) + { + cell->set_future_fe_index(1); + break; + } + + // transfer future FE indices + dofh_dst.set_future_fe_indices(dofh_src.get_future_fe_indices()); + + // count copied future FE indices + const unsigned int n_indices_copied = + std::count_if(dofh_dst.active_cell_iterators().begin(), + dofh_dst.active_cell_iterators().end(), + [](const auto &cell) { + return (cell->is_locally_owned() && + cell->future_fe_index_set()); + }); + deallog << "indices copied: " << n_indices_copied << std::endl; + + // verify that future FE indices match + for (const auto &cell : tria.active_cell_iterators()) + if (cell->is_locally_owned()) + { + const auto cell_src = cell->as_dof_handler_iterator(dofh_src); + const auto cell_dst = cell->as_dof_handler_iterator(dofh_dst); + + AssertThrow(cell_src->future_fe_index_set() == + cell_dst->future_fe_index_set(), + ExcInternalError()); + AssertThrow(cell_src->future_fe_index() == cell_dst->future_fe_index(), + ExcInternalError()); + } + + deallog << "OK" << std::endl; +} + + +int +main(int argc, char *argv[]) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + MPILogInitAll log; + + test<2, 2>(); +} diff --git a/tests/mpi/hp_getset_future_fe_indices_01.with_p4est=true.mpirun=2.output b/tests/mpi/hp_getset_future_fe_indices_01.with_p4est=true.mpirun=2.output new file mode 100644 index 0000000000..a6f21ad876 --- /dev/null +++ b/tests/mpi/hp_getset_future_fe_indices_01.with_p4est=true.mpirun=2.output @@ -0,0 +1,7 @@ + +DEAL:0::indices copied: 1 +DEAL:0::OK + +DEAL:1::indices copied: 1 +DEAL:1::OK + -- 2.39.5