From: Alexander Grayver Date: Mon, 11 Dec 2017 16:25:02 +0000 (+0100) Subject: Update doc, add changelog item X-Git-Tag: v9.0.0-rc1~665^2~1 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=d12834fb7493aaadb8537b36169fc0434ef8e98a;p=dealii.git Update doc, add changelog item --- diff --git a/doc/news/changes/incompatibilities/20171211AlexanderGrayver b/doc/news/changes/incompatibilities/20171211AlexanderGrayver new file mode 100644 index 0000000000..adfa5aba70 --- /dev/null +++ b/doc/news/changes/incompatibilities/20171211AlexanderGrayver @@ -0,0 +1,4 @@ +Changed: Change the logic of DoFTools::get_subdomain_association to have the degrees of freedom along a refinement edge be now all given to the processor with the smallest subdomain_id. +
+(Alexander Grayver, 2017/12/11) + diff --git a/include/deal.II/dofs/dof_tools.h b/include/deal.II/dofs/dof_tools.h index 05d5c29632..829cd6b799 100644 --- a/include/deal.II/dofs/dof_tools.h +++ b/include/deal.II/dofs/dof_tools.h @@ -1668,13 +1668,13 @@ namespace DoFTools * * Note that degrees of freedom associated with faces, edges, and vertices * may be associated with multiple subdomains if they are sitting on - * partition boundaries. In these cases, we put them into one of the - * associated partitions in an undefined way. This may sometimes lead to - * different numbers of degrees of freedom in partitions, even if the number - * of cells is perfectly equidistributed. While this is regrettable, it is - * not a problem in practice since the number of degrees of freedom on - * partition boundaries is asymptotically vanishing as we refine the mesh as - * long as the number of partitions is kept constant. + * partition boundaries. In these cases, we assign them to the process with + * the smaller subdomain id. This may lead to different numbers of degrees + * of freedom in partitions, even if the number of cells is perfectly + * equidistributed. While this is regrettable, it is not a problem in + * practice since the number of degrees of freedom on partition boundaries + * is asymptotically vanishing as we refine the mesh as long as the number + * of partitions is kept constant. * * This function returns the association of each DoF with one subdomain. If * you are looking for the association of each @em cell with a subdomain, diff --git a/source/dofs/dof_tools.cc b/source/dofs/dof_tools.cc index d58d735a82..7ad39de9c3 100644 --- a/source/dofs/dof_tools.cc +++ b/source/dofs/dof_tools.cc @@ -1460,7 +1460,7 @@ namespace DoFTools local_dof_indices.reserve (max_dofs_per_cell(dof_handler)); // loop over all cells and record which subdomain a DoF belongs to. - // toss a coin in case it is on an interface + // give to the smaller subdomain_id in case it is on an interface typename DoFHandlerType::active_cell_iterator cell = dof_handler.begin_active(), endc = dof_handler.end(); diff --git a/tests/mpi/blockwise_parallel_01.cc b/tests/mpi/blockwise_parallel_01.cc new file mode 100644 index 0000000000..73fd9560f7 --- /dev/null +++ b/tests/mpi/blockwise_parallel_01.cc @@ -0,0 +1,78 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2015 - 2017 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 DoFRenumbering::block_wise with parallel::shared::Triangulation + +#include "../tests.h" + +// all include files you need here +#include +#include +#include +#include +#include +#include +#include +#include + +template +void test() +{ + const unsigned int this_mpi_process = + Utilities::MPI::this_mpi_process(MPI_COMM_WORLD); + + FESystem fe(FE_Q(1), 2); + + { + parallel::shared::Triangulation tria(MPI_COMM_WORLD, Triangulation::none, false, + parallel::shared::Triangulation::Settings::partition_zorder); + GridGenerator::hyper_cube(tria, -1, 1); + tria.refine_global(2); + + DoFHandler dh(tria); + dh.distribute_dofs(fe); + + DoFRenumbering::block_wise(dh); + + const std::vector locally_owned_dofs_per_subdomain = + DoFTools::locally_owned_dofs_per_subdomain(dh); + + const types::global_dof_index dofs_per_block = dh.n_dofs() / 2; + std::vector locally_owned_dofs_per_block(2); + locally_owned_dofs_per_block[0] = + locally_owned_dofs_per_subdomain[this_mpi_process].get_view(0, dofs_per_block); + locally_owned_dofs_per_block[1] = + locally_owned_dofs_per_subdomain[this_mpi_process].get_view(dofs_per_block, dh.n_dofs()); + + // This assertion fails + Assert((locally_owned_dofs_per_block[0] == locally_owned_dofs_per_block[1]), + ExcMessage("Locally owned dofs differ across blocks.")); + } + + if (Utilities::MPI::this_mpi_process (MPI_COMM_WORLD) == 0) + deallog << "OK for " << dim << "d" << std::endl; +} + +int main (int argc, char* argv[]) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv, 1); + MPILogInitAll log; + + test<2>(); + test<3>(); + + return 0; +} + diff --git a/tests/mpi/blockwise_parallel_01.mpirun=2.output b/tests/mpi/blockwise_parallel_01.mpirun=2.output new file mode 100644 index 0000000000..5637e3faa0 --- /dev/null +++ b/tests/mpi/blockwise_parallel_01.mpirun=2.output @@ -0,0 +1,3 @@ +DEAL:0::OK for 2d +DEAL:0::OK for 3d +