From: Wolfgang Bangerth Date: Sun, 20 Mar 2016 18:52:31 +0000 (-0500) Subject: Fix a segmentation fault in DoFHandler::locally_owned_dofs(). X-Git-Tag: v8.5.0-rc1~1187^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=refs%2Fpull%2F2388%2Fhead;p=dealii.git Fix a segmentation fault in DoFHandler::locally_owned_dofs(). This was caused by an incorrect size computation inside DoFTools::locally_owned_dofs_per_subdomain(). --- diff --git a/doc/news/changes.h b/doc/news/changes.h index ec5a260aee..f4384b45a2 100644 --- a/doc/news/changes.h +++ b/doc/news/changes.h @@ -104,6 +104,13 @@ inconvenience this causes.

Specific improvements

    +
  1. Fixed: DoFHandler::locally_owned_dofs() could create a segmentation + fault in cases where some processors do not own any cells. This was caused + by an incorrect computation in DoFTools::locally_owned_dofs_per_subdomain(). +
    + (Wolfgang Bangerth, 2016/03/20) +
  2. +
  3. Improved: The distribution of degrees of freedom on multigrid levels, DoFHandler::distribute_mg_dofs(), contained a few steps that scaled quadratically in the number of local cells for certain configurations. These diff --git a/source/dofs/dof_tools.cc b/source/dofs/dof_tools.cc index 1aa3e161a3..502e791cb5 100644 --- a/source/dofs/dof_tools.cc +++ b/source/dofs/dof_tools.cc @@ -1145,10 +1145,22 @@ namespace DoFTools std::vector< dealii::types::subdomain_id > subdomain_association (dof_handler.n_dofs ()); dealii::DoFTools::get_subdomain_association (dof_handler, subdomain_association); - const unsigned int n_subdomains = 1 + (*std::max_element (subdomain_association.begin (), - subdomain_association.end () )); + const unsigned int n_subdomains + = (dynamic_cast *> + (&dof_handler.get_triangulation()) == 0 + ? + 1 + : + Utilities::MPI::n_mpi_processes + (dynamic_cast *> + (&dof_handler.get_triangulation())->get_communicator())); + Assert (n_subdomains > + *std::max_element (subdomain_association.begin (), + subdomain_association.end ()), + ExcInternalError()); - std::vector index_sets (n_subdomains,dealii::IndexSet(dof_handler.n_dofs())); + std::vector index_sets (n_subdomains, + dealii::IndexSet(dof_handler.n_dofs())); // loop over subdomain_association and populate IndexSet when a // change in subdomain ID is found diff --git a/tests/mpi/locally_owned_dofs_per_subdomain.cc b/tests/mpi/locally_owned_dofs_per_subdomain.cc new file mode 100644 index 0000000000..460119e945 --- /dev/null +++ b/tests/mpi/locally_owned_dofs_per_subdomain.cc @@ -0,0 +1,59 @@ +// --------------------------------------------------------------------- +// +// 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. +// +// --------------------------------------------------------------------- + + +// DoFTools::locally_owned_dofs_per_subdomain wrongly sized the array +// it returns when we have some processors that don't own any cells + +#include "../tests.h" +#include +#include +#include +#include + + +template +void test () +{ + parallel::shared::Triangulation triangulation(MPI_COMM_WORLD); + GridGenerator::hyper_cube(triangulation, -1.0, 1.0); + + FE_Q fe(1); + DoFHandler dof_handler (triangulation); + dof_handler.distribute_dofs(fe); + + // this used to crash here: + const IndexSet locally_owned_dofs = dof_handler.locally_owned_dofs(); + + deallog << "dim=" << dim + << std::endl + << "rank=" << Utilities::MPI::this_mpi_process(MPI_COMM_WORLD) + << std::endl + << " n_dofs=" << dof_handler.n_dofs() + << std::endl + << " n_locally_owned_dofs=" + << locally_owned_dofs.n_elements() + << std::endl; +} + +int main (int argc, char **argv) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv, testing_max_num_threads()); + + MPILogInitAll all; + + test<2>(); + test<3>(); +} diff --git a/tests/mpi/locally_owned_dofs_per_subdomain.mpirun=2.output b/tests/mpi/locally_owned_dofs_per_subdomain.mpirun=2.output new file mode 100644 index 0000000000..481ba49a0d --- /dev/null +++ b/tests/mpi/locally_owned_dofs_per_subdomain.mpirun=2.output @@ -0,0 +1,19 @@ + +DEAL:0::dim=2 +DEAL:0::rank=0 +DEAL:0:: n_dofs=4 +DEAL:0:: n_locally_owned_dofs=4 +DEAL:0::dim=3 +DEAL:0::rank=0 +DEAL:0:: n_dofs=8 +DEAL:0:: n_locally_owned_dofs=8 + +DEAL:1::dim=2 +DEAL:1::rank=1 +DEAL:1:: n_dofs=4 +DEAL:1:: n_locally_owned_dofs=0 +DEAL:1::dim=3 +DEAL:1::rank=1 +DEAL:1:: n_dofs=8 +DEAL:1:: n_locally_owned_dofs=0 +