From: Peter Munch Date: Wed, 24 Jun 2020 05:23:54 +0000 (+0200) Subject: Fix parallel::distributed::SolutionTransfer for FENothing X-Git-Tag: v9.3.0-rc1~1383^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=refs%2Fpull%2F10592%2Fhead;p=dealii.git Fix parallel::distributed::SolutionTransfer for FENothing --- diff --git a/doc/news/changes/minor/20200624SoldnerMunch b/doc/news/changes/minor/20200624SoldnerMunch new file mode 100644 index 0000000000..3354fd86ec --- /dev/null +++ b/doc/news/changes/minor/20200624SoldnerMunch @@ -0,0 +1,4 @@ +Fixed: The class parallel::distributed::SolutionTransfer can now +also handle FE_Nothing. +
+(Dominic Soldner, Peter Munch, 2020/06/24) diff --git a/source/distributed/solution_transfer.cc b/source/distributed/solution_transfer.cc index 52ca51e57f..d4b7672851 100644 --- a/source/distributed/solution_transfer.cc +++ b/source/distributed/solution_transfer.cc @@ -365,6 +365,9 @@ namespace parallel const unsigned int dofs_per_cell = dof_handler->get_fe(fe_index).dofs_per_cell; + if (dofs_per_cell == 0) + return std::vector(); // nothing to do for FE_Nothing + auto it_input = input_vectors.cbegin(); auto it_output = dof_values.begin(); for (; it_input != input_vectors.cend(); ++it_input, ++it_output) @@ -437,6 +440,9 @@ namespace parallel const unsigned int dofs_per_cell = dof_handler->get_fe(fe_index).dofs_per_cell; + if (dofs_per_cell == 0) + return; // nothing to do for FE_Nothing + const std::vector<::dealii::Vector> dof_values = unpack_dof_values(data_range, diff --git a/tests/hp/distributed_solution_transfer_01.cc b/tests/hp/distributed_solution_transfer_01.cc new file mode 100644 index 0000000000..bda8020149 --- /dev/null +++ b/tests/hp/distributed_solution_transfer_01.cc @@ -0,0 +1,116 @@ +// --------------------------------------------------------------------- +// +// 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. +// +// --------------------------------------------------------------------- + + + +// Test parallel::distributed::SolutionTransfer for FE_Nothing +// (see also https://github.com/dealii/dealii/issues/10570). + +#include + +#include +#include + +#include +#include +#include + +#include +#include + +#include +#include +#include +#include +#include + +#include + +#include +#include + +#include "../tests.h" + + +template +void +transfer(const MPI_Comm &comm) +{ + AssertDimension(Utilities::MPI::n_mpi_processes(comm), 1); + + parallel::distributed::Triangulation tria(comm); + GridGenerator::hyper_cube(tria); + tria.refine_global(1); + + hp::FECollection fe; + fe.push_back(FE_Q(1)); + fe.push_back(FE_Nothing()); + + // create a DoFHandler on which we have both cells with FE_Q as well as + // FE_Nothing + DoFHandler dof_handler(tria, true); + dof_handler.begin(0)->child(0)->set_active_fe_index(1); + + IndexSet locally_relevant_dofs; + DoFTools::extract_locally_relevant_dofs(dof_handler, locally_relevant_dofs); + + LinearAlgebra::distributed::Vector solution( + dof_handler.locally_owned_dofs(), locally_relevant_dofs, comm); + AffineConstraints cm; + cm.close(); + + dof_handler.distribute_dofs(fe); + solution.reinit(dof_handler.n_dofs()); + + for (unsigned int i = 0; i < solution.size(); ++i) + solution(i) = i; + + parallel::distributed:: + SolutionTransfer> + soltrans(dof_handler); + + for (const auto &cell : tria.active_cell_iterators()) + cell->set_refine_flag(); + + LinearAlgebra::distributed::Vector old_solution = solution; + + tria.prepare_coarsening_and_refinement(); + soltrans.prepare_for_coarsening_and_refinement(old_solution); + tria.execute_coarsening_and_refinement(); + dof_handler.distribute_dofs(fe); + solution.reinit(dof_handler.n_dofs()); + + soltrans.interpolate(solution); + + deallog << "OK" << std::endl; +} + + +int +main(int argc, char *argv[]) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + + const MPI_Comm comm = MPI_COMM_WORLD; + initlog(); + + deallog.push("2D solution transfer"); + transfer<2>(comm); + deallog.pop(); + + deallog.push("3D solution transfer"); + transfer<3>(comm); + deallog.pop(); +} diff --git a/tests/hp/distributed_solution_transfer_01.mpirun=1.output b/tests/hp/distributed_solution_transfer_01.mpirun=1.output new file mode 100644 index 0000000000..f368bb0bb4 --- /dev/null +++ b/tests/hp/distributed_solution_transfer_01.mpirun=1.output @@ -0,0 +1,3 @@ + +DEAL:2D solution transfer::OK +DEAL:3D solution transfer::OK