From d5279b85a61c18d3de54de9f181fa682a0ec6d54 Mon Sep 17 00:00:00 2001 From: Timo Heister Date: Sun, 31 Jul 2016 15:27:21 +0200 Subject: [PATCH] Fix locally_owned_dofs_per_subdomain DoFTools::locally_owned_dofs_per_subdomain would wrongly trigger an Assertion if running with a dealii::Triangulation (not a shared one) and with more than one MPI rank. This is now fixed. --- source/dofs/dof_tools.cc | 3 +- .../locally_owned_dofs_per_subdomain_02.cc | 60 +++++++++++++++++++ ...wned_dofs_per_subdomain_02.mpirun=2.output | 11 ++++ 3 files changed, 73 insertions(+), 1 deletion(-) create mode 100644 tests/mpi/locally_owned_dofs_per_subdomain_02.cc create mode 100644 tests/mpi/locally_owned_dofs_per_subdomain_02.mpirun=2.output diff --git a/source/dofs/dof_tools.cc b/source/dofs/dof_tools.cc index 251e0d641e..738da85a2d 100644 --- a/source/dofs/dof_tools.cc +++ b/source/dofs/dof_tools.cc @@ -1149,7 +1149,8 @@ namespace DoFTools = (dynamic_cast *> (&dof_handler.get_triangulation()) == 0 ? - 1 + (1+*std::max_element (subdomain_association.begin (), + subdomain_association.end ())) : Utilities::MPI::n_mpi_processes (dynamic_cast *> diff --git a/tests/mpi/locally_owned_dofs_per_subdomain_02.cc b/tests/mpi/locally_owned_dofs_per_subdomain_02.cc new file mode 100644 index 0000000000..b02a062f93 --- /dev/null +++ b/tests/mpi/locally_owned_dofs_per_subdomain_02.cc @@ -0,0 +1,60 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2016 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 at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + + +// Test DoFTools::locally_owned_dofs_per_subdomain (for a standard Triangulation) + +#include "../tests.h" +#include +#include +#include +#include +#include +#include + + +template +void test () +{ + Triangulation triangulation; + GridGenerator::hyper_cube(triangulation, -1.0, 1.0); + triangulation.refine_global(1); + + const unsigned int nproc = Utilities::MPI::n_mpi_processes(MPI_COMM_WORLD); + + GridTools::partition_triangulation (nproc, + triangulation); + + FE_Q fe(1); + DoFHandler dof_handler (triangulation); + dof_handler.distribute_dofs(fe); + + std::vector locally_owned_dofs_per_proc = DoFTools::locally_owned_dofs_per_subdomain(dof_handler); + + for (unsigned int p=0; p(); +} diff --git a/tests/mpi/locally_owned_dofs_per_subdomain_02.mpirun=2.output b/tests/mpi/locally_owned_dofs_per_subdomain_02.mpirun=2.output new file mode 100644 index 0000000000..ffe0b4dc22 --- /dev/null +++ b/tests/mpi/locally_owned_dofs_per_subdomain_02.mpirun=2.output @@ -0,0 +1,11 @@ + +DEAL:0::proc 0: +{[0,1], [4,5]} +DEAL:0::proc 1: +{[2,3], [6,8]} + +DEAL:1::proc 0: +{[0,1], [4,5]} +DEAL:1::proc 1: +{[2,3], [6,8]} + -- 2.39.5