From: tcclevenger Date: Tue, 19 Feb 2019 23:48:34 +0000 (-0700) Subject: add random renumbering X-Git-Tag: v9.1.0-rc1~328^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=refs%2Fpull%2F7739%2Fhead;p=dealii.git add random renumbering --- diff --git a/doc/news/changes/minor/20190219ConradClevenger b/doc/news/changes/minor/20190219ConradClevenger new file mode 100644 index 0000000000..aa560d10d9 --- /dev/null +++ b/doc/news/changes/minor/20190219ConradClevenger @@ -0,0 +1,4 @@ +New: Function DoFRenumbering::random(dof_handler, level) which allows for a random renumbering +of degrees of freedom on a level in a mutilevel hierarchy. +
+(Conrad Clevenger, 2019/02/19) diff --git a/include/deal.II/dofs/dof_renumbering.h b/include/deal.II/dofs/dof_renumbering.h index d94420c66e..445fedd078 100644 --- a/include/deal.II/dofs/dof_renumbering.h +++ b/include/deal.II/dofs/dof_renumbering.h @@ -1133,6 +1133,16 @@ namespace DoFRenumbering void random(DoFHandlerType &dof_handler); + /** + * Renumber the degrees of freedom in a random way. It does the same thing as + * the above function, only that it does this for one single level of a + * multilevel discretization. The non-multigrid part of the DoFHandler + * is not touched. + */ + template + void + random(DoFHandlerType &dof_handler, const unsigned int level); + /** * Compute the renumbering vector needed by the random() function. See * there for more information on the computed random renumbering. @@ -1145,6 +1155,17 @@ namespace DoFRenumbering compute_random(std::vector &new_dof_indices, const DoFHandlerType & dof_handler); + /** + * Compute the renumbering vector needed by the random() function. Same + * as the above function but for a single level of a multilevel + * discretization. + */ + template + void + compute_random(std::vector &new_dof_indices, + const DoFHandlerType & dof_handler, + const unsigned int level); + /** * @} */ diff --git a/source/dofs/dof_renumbering.cc b/source/dofs/dof_renumbering.cc index bd96818aca..adc401e0f7 100644 --- a/source/dofs/dof_renumbering.cc +++ b/source/dofs/dof_renumbering.cc @@ -2129,6 +2129,24 @@ namespace DoFRenumbering + template + void + random(DoFHandlerType &dof_handler, const unsigned int level) + { + Assert(dof_handler.n_dofs(level) != numbers::invalid_dof_index, + ExcDoFHandlerNotInitialized()); + + std::vector renumbering( + dof_handler.locally_owned_mg_dofs(level).n_elements(), + numbers::invalid_dof_index); + + compute_random(renumbering, dof_handler, level); + + dof_handler.renumber_dofs(level, renumbering); + } + + + template void compute_random(std::vector &new_indices, @@ -2159,6 +2177,37 @@ namespace DoFRenumbering + template + void + compute_random(std::vector &new_indices, + const DoFHandlerType & dof_handler, + const unsigned int level) + { + const types::global_dof_index n_dofs = dof_handler.n_dofs(level); + Assert(new_indices.size() == n_dofs, + ExcDimensionMismatch(new_indices.size(), n_dofs)); + + for (unsigned int i = 0; i < n_dofs; ++i) + new_indices[i] = i; + + // shuffle the elements; the following is essentially std::shuffle (which + // is new in C++11) but with a boost URNG + ::boost::mt19937 random_number_generator; + for (unsigned int i = 1; i < n_dofs; ++i) + { + // get a random number between 0 and i (inclusive) + const unsigned int j = + ::boost::random::uniform_int_distribution<>(0, i)( + random_number_generator); + + // if possible, swap the elements + if (i != j) + std::swap(new_indices[i], new_indices[j]); + } + } + + + template void subdomain_wise(DoFHandlerType &dof_handler) diff --git a/source/dofs/dof_renumbering.inst.in b/source/dofs/dof_renumbering.inst.in index f7667f4b83..7ff6480386 100644 --- a/source/dofs/dof_renumbering.inst.in +++ b/source/dofs/dof_renumbering.inst.in @@ -293,6 +293,16 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS) std::vector &, const DoFHandler &); + template void + random>(DoFHandler &, + const unsigned int); + + template void + compute_random>( + std::vector &, + const DoFHandler &, + const unsigned int); + template void sort_selected_dofs_back>( DoFHandler &, diff --git a/tests/dofs/dof_renumbering_09.cc b/tests/dofs/dof_renumbering_09.cc new file mode 100644 index 0000000000..7d3e0dcc3b --- /dev/null +++ b/tests/dofs/dof_renumbering_09.cc @@ -0,0 +1,108 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2000 - 2018 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. +// +// --------------------------------------------------------------------- + + + +// tests DoFRenumbering::random function for levels + + + +#include + +#include +#include +#include + +#include +#include +#include +#include + +#include +#include +#include + +#include + +#include "../tests.h" + + + +template +void +print_dofs(const DoFHandler &dof, unsigned int level) +{ + std::vector v(dof.get_fe().dofs_per_cell); + for (typename DoFHandler::cell_iterator cell = dof.begin(level); + cell != dof.end(level); + ++cell) + { + deallog << "Cell " << cell << " -- "; + cell->get_mg_dof_indices(v); + for (unsigned int i = 0; i < v.size(); ++i) + deallog << v[i] << ' '; + deallog << std::endl; + } +} + + +template +void +test() +{ + Triangulation tr(Triangulation::limit_level_difference_at_vertices); + GridGenerator::hyper_cube(tr); + tr.refine_global((dim == 1 ? 2 : 1)); + + DoFHandler dof_handler(tr); + + FE_Q fe(1); + dof_handler.distribute_dofs(fe); + dof_handler.distribute_mg_dofs(fe); + + + // Print dofs before reordering + deallog << "Before reorder: " << std::endl; + for (unsigned int level = 0; level < tr.n_levels(); ++level) + { + deallog << "Level " << level << std::endl; + print_dofs(dof_handler, level); + } + + for (unsigned int level = 0; level < tr.n_levels(); ++level) + DoFRenumbering::random(dof_handler, level); + + // Print dofs after reordering + deallog << "After reorder: " << std::endl; + for (unsigned int level = 0; level < tr.n_levels(); ++level) + { + deallog << "Level " << level << std::endl; + print_dofs(dof_handler, level); + } +} + + +int +main() +{ + initlog(); + + deallog << "1D" << std::endl; + test<1>(); + deallog << std::endl << "2D" << std::endl; + test<2>(); + deallog << std::endl << "3D" << std::endl; + test<3>(); +} diff --git a/tests/dofs/dof_renumbering_09.output b/tests/dofs/dof_renumbering_09.output new file mode 100644 index 0000000000..a82919a5f9 --- /dev/null +++ b/tests/dofs/dof_renumbering_09.output @@ -0,0 +1,68 @@ + +DEAL::1D +DEAL::Before reorder: +DEAL::Level 0 +DEAL::Cell 0.0 -- 0 1 +DEAL::Level 1 +DEAL::Cell 1.0 -- 0 1 +DEAL::Cell 1.1 -- 1 2 +DEAL::Level 2 +DEAL::Cell 2.0 -- 0 1 +DEAL::Cell 2.1 -- 1 2 +DEAL::Cell 2.2 -- 2 3 +DEAL::Cell 2.3 -- 3 4 +DEAL::After reorder: +DEAL::Level 0 +DEAL::Cell 0.0 -- 0 1 +DEAL::Level 1 +DEAL::Cell 1.0 -- 2 1 +DEAL::Cell 1.1 -- 1 0 +DEAL::Level 2 +DEAL::Cell 2.0 -- 2 1 +DEAL::Cell 2.1 -- 1 0 +DEAL::Cell 2.2 -- 0 3 +DEAL::Cell 2.3 -- 3 4 +DEAL:: +DEAL::2D +DEAL::Before reorder: +DEAL::Level 0 +DEAL::Cell 0.0 -- 0 1 2 3 +DEAL::Level 1 +DEAL::Cell 1.0 -- 0 1 2 3 +DEAL::Cell 1.1 -- 1 4 3 5 +DEAL::Cell 1.2 -- 2 3 6 7 +DEAL::Cell 1.3 -- 3 5 7 8 +DEAL::After reorder: +DEAL::Level 0 +DEAL::Cell 0.0 -- 2 1 0 3 +DEAL::Level 1 +DEAL::Cell 1.0 -- 5 8 0 3 +DEAL::Cell 1.1 -- 8 4 3 2 +DEAL::Cell 1.2 -- 0 3 6 7 +DEAL::Cell 1.3 -- 3 2 7 1 +DEAL:: +DEAL::3D +DEAL::Before reorder: +DEAL::Level 0 +DEAL::Cell 0.0 -- 0 1 2 3 4 5 6 7 +DEAL::Level 1 +DEAL::Cell 1.0 -- 0 1 2 3 4 5 6 7 +DEAL::Cell 1.1 -- 1 8 3 9 5 10 7 11 +DEAL::Cell 1.2 -- 2 3 12 13 6 7 14 15 +DEAL::Cell 1.3 -- 3 9 13 16 7 11 15 17 +DEAL::Cell 1.4 -- 4 5 6 7 18 19 20 21 +DEAL::Cell 1.5 -- 5 10 7 11 19 22 21 23 +DEAL::Cell 1.6 -- 6 7 14 15 20 21 24 25 +DEAL::Cell 1.7 -- 7 11 15 17 21 23 25 26 +DEAL::After reorder: +DEAL::Level 0 +DEAL::Cell 0.0 -- 5 1 0 3 4 2 6 7 +DEAL::Level 1 +DEAL::Cell 1.0 -- 5 11 26 21 4 2 9 12 +DEAL::Cell 1.1 -- 11 15 21 6 2 3 12 8 +DEAL::Cell 1.2 -- 26 21 7 10 9 12 0 1 +DEAL::Cell 1.3 -- 21 6 10 22 12 8 1 17 +DEAL::Cell 1.4 -- 4 2 9 12 18 19 20 13 +DEAL::Cell 1.5 -- 2 3 12 8 19 16 13 23 +DEAL::Cell 1.6 -- 9 12 0 1 20 13 25 24 +DEAL::Cell 1.7 -- 12 8 1 17 13 23 24 14