From: Marc Fehling Date: Tue, 27 Aug 2019 08:20:55 +0000 (+0200) Subject: Minor changes: tests/mpi/solution_transfer_05. X-Git-Tag: v9.2.0-rc1~1174^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=24cf4fefcc568ca3ade21d1965dfe9a7fbd71e80;p=dealii.git Minor changes: tests/mpi/solution_transfer_05. --- diff --git a/doc/news/changes/minor/20190827Fehling b/doc/news/changes/minor/20190827Fehling new file mode 100644 index 0000000000..1ac4552e20 --- /dev/null +++ b/doc/news/changes/minor/20190827Fehling @@ -0,0 +1,9 @@ +Bugfix: For parallel::distributed::Triangulation objects, p4est determines +which cells will be refined and coarsened. Thus, refinement flags are not +a good measure to predict refinement behavior for transferring active finite +element indices. We now provide a matching coarsening strategy to the +parallel::distributed::CellDataTransfer object responsible for their transfer. +This fixes assertions being triggered in +parallel::distributed::SolutionTransfer::interpolate(). +
+(Marc Fehling, 2019/08/27) diff --git a/tests/mpi/solution_transfer_05.cc b/tests/mpi/solution_transfer_05.cc index cdc4a8aaff..0c4d7d4c72 100644 --- a/tests/mpi/solution_transfer_05.cc +++ b/tests/mpi/solution_transfer_05.cc @@ -1,6 +1,6 @@ // --------------------------------------------------------------------- // -// Copyright (C) 1998 - 2018 by the deal.II authors +// Copyright (C) 2019 by the deal.II authors // // This file is part of the deal.II library. // @@ -15,6 +15,13 @@ +// This test crashed at some point: We have set and sent active_fe_indices based +// on the refinement flags on the p::d::Triangulation object. However, p4est has +// the last word on deciding which cells will be refined -- and p4est makes use +// of it in the specific scenario provided as a test. A fix has been introduced +// along with this test. + + #include #include @@ -34,16 +41,14 @@ #include -#include -#include - #include "../tests.h" template void -transfer(std::ostream &out) +test() { + // setup parallel::distributed::Triangulation tria(MPI_COMM_WORLD); GridGenerator::hyper_cube(tria); tria.refine_global(2); @@ -55,32 +60,34 @@ transfer(std::ostream &out) hp::DoFHandler dgq_dof_handler(tria); - // randomly assign FE orders + // randomly assign fes for (const auto &cell : dgq_dof_handler.active_cell_iterators()) if (cell->is_locally_owned()) cell->set_active_fe_index(Testing::rand() % max_degree); dgq_dof_handler.distribute_dofs(fe_dgq); + // prepare index sets IndexSet dgq_locally_owned_dofs = dgq_dof_handler.locally_owned_dofs(); IndexSet dgq_locally_relevant_dofs; - dealii::DoFTools::extract_locally_relevant_dofs(dgq_dof_handler, - dgq_locally_relevant_dofs); + DoFTools::extract_locally_relevant_dofs(dgq_dof_handler, + dgq_locally_relevant_dofs); IndexSet dgq_ghost_dofs = dgq_locally_relevant_dofs; dgq_ghost_dofs.subtract_set(dgq_locally_owned_dofs); + // prepare dof_values LinearAlgebra::distributed::Vector dgq_solution; dgq_solution.reinit(dgq_locally_owned_dofs, dgq_ghost_dofs, MPI_COMM_WORLD); - VectorTools::interpolate(dgq_dof_handler, ZeroFunction(), dgq_solution); + dgq_solution.update_ghost_values(); parallel::distributed::SolutionTransfer< dim, LinearAlgebra::distributed::Vector, hp::DoFHandler> dgq_soltrans(dgq_dof_handler); + dgq_soltrans.prepare_for_coarsening_and_refinement(dgq_solution); - LinearAlgebra::distributed::Vector dgq_old_solution = dgq_solution; - dgq_old_solution.update_ghost_values(); + // refine and transfer { unsigned int counter = 0; for (auto cell = tria.begin_active(); cell != tria.end(); ++cell, ++counter) @@ -93,18 +100,17 @@ transfer(std::ostream &out) } } - tria.prepare_coarsening_and_refinement(); - dgq_soltrans.prepare_for_coarsening_and_refinement(dgq_old_solution); tria.execute_coarsening_and_refinement(); - dgq_dof_handler.distribute_dofs(fe_dgq); + // prepare index sets dgq_locally_owned_dofs = dgq_dof_handler.locally_owned_dofs(); - dealii::DoFTools::extract_locally_relevant_dofs(dgq_dof_handler, - dgq_locally_relevant_dofs); + DoFTools::extract_locally_relevant_dofs(dgq_dof_handler, + dgq_locally_relevant_dofs); dgq_ghost_dofs = dgq_locally_relevant_dofs; dgq_ghost_dofs.subtract_set(dgq_locally_owned_dofs); + // unpack dof_values dgq_solution.reinit(dgq_locally_owned_dofs, dgq_ghost_dofs, MPI_COMM_WORLD); dgq_soltrans.interpolate(dgq_solution); @@ -119,9 +125,9 @@ main(int argc, char *argv[]) MPILogInitAll log; deallog.push("2d"); - transfer<2>(deallog.get_file_stream()); + test<2>(); deallog.pop(); deallog.push("3d"); - transfer<3>(deallog.get_file_stream()); + test<3>(); deallog.pop(); } diff --git a/tests/mpi/solution_transfer_05.with_p4est=true.with_trilinos=true.mpirun=1.output b/tests/mpi/solution_transfer_05.with_p4est=true.mpirun=1.output similarity index 100% rename from tests/mpi/solution_transfer_05.with_p4est=true.with_trilinos=true.mpirun=1.output rename to tests/mpi/solution_transfer_05.with_p4est=true.mpirun=1.output diff --git a/tests/mpi/solution_transfer_05.with_p4est=true.with_trilinos=true.mpirun=2.output b/tests/mpi/solution_transfer_05.with_p4est=true.mpirun=2.output similarity index 100% rename from tests/mpi/solution_transfer_05.with_p4est=true.with_trilinos=true.mpirun=2.output rename to tests/mpi/solution_transfer_05.with_p4est=true.mpirun=2.output diff --git a/tests/mpi/solution_transfer_05.with_p4est=true.with_trilinos=true.mpirun=4.output b/tests/mpi/solution_transfer_05.with_p4est=true.mpirun=4.output similarity index 100% rename from tests/mpi/solution_transfer_05.with_p4est=true.with_trilinos=true.mpirun=4.output rename to tests/mpi/solution_transfer_05.with_p4est=true.mpirun=4.output