]> https://gitweb.dealii.org/ - dealii.git/commitdiff
add random renumbering 7739/head
authortcclevenger <tcleven@clemson.edu>
Tue, 19 Feb 2019 23:48:34 +0000 (16:48 -0700)
committertcclevenger <tcleven@clemson.edu>
Wed, 20 Feb 2019 00:47:50 +0000 (17:47 -0700)
doc/news/changes/minor/20190219ConradClevenger [new file with mode: 0644]
include/deal.II/dofs/dof_renumbering.h
source/dofs/dof_renumbering.cc
source/dofs/dof_renumbering.inst.in
tests/dofs/dof_renumbering_09.cc [new file with mode: 0644]
tests/dofs/dof_renumbering_09.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20190219ConradClevenger b/doc/news/changes/minor/20190219ConradClevenger
new file mode 100644 (file)
index 0000000..aa560d1
--- /dev/null
@@ -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. 
+<br>
+(Conrad Clevenger, 2019/02/19)
index d94420c66e7858921c577959d2c98c4641f413bb..445fedd07883d1a0d89441e9a60196daa0d131a1 100644 (file)
@@ -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 <typename DoFHandlerType>
+  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<types::global_dof_index> &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 <typename DoFHandlerType>
+  void
+  compute_random(std::vector<types::global_dof_index> &new_dof_indices,
+                 const DoFHandlerType &                dof_handler,
+                 const unsigned int                    level);
+
   /**
    * @}
    */
index bd96818aca415b5c31bab31f328fe22608c9145c..adc401e0f70d37b2ff5f4c434e4546ef018fa613 100644 (file)
@@ -2129,6 +2129,24 @@ namespace DoFRenumbering
 
 
 
+  template <typename DoFHandlerType>
+  void
+  random(DoFHandlerType &dof_handler, const unsigned int level)
+  {
+    Assert(dof_handler.n_dofs(level) != numbers::invalid_dof_index,
+           ExcDoFHandlerNotInitialized());
+
+    std::vector<types::global_dof_index> 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 <typename DoFHandlerType>
   void
   compute_random(std::vector<types::global_dof_index> &new_indices,
@@ -2159,6 +2177,37 @@ namespace DoFRenumbering
 
 
 
+  template <typename DoFHandlerType>
+  void
+  compute_random(std::vector<types::global_dof_index> &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 <typename DoFHandlerType>
   void
   subdomain_wise(DoFHandlerType &dof_handler)
index f7667f4b8316643cb37496cc0fd4bb5d5b80099b..7ff6480386efc267f7fade25e4b63f59c4c6bea5 100644 (file)
@@ -293,6 +293,16 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS)
         std::vector<types::global_dof_index> &,
         const DoFHandler<deal_II_dimension> &);
 
+      template void
+      random<DoFHandler<deal_II_dimension>>(DoFHandler<deal_II_dimension> &,
+                                            const unsigned int);
+
+      template void
+      compute_random<DoFHandler<deal_II_dimension>>(
+        std::vector<types::global_dof_index> &,
+        const DoFHandler<deal_II_dimension> &,
+        const unsigned int);
+
       template void
       sort_selected_dofs_back<DoFHandler<deal_II_dimension>>(
         DoFHandler<deal_II_dimension> &,
diff --git a/tests/dofs/dof_renumbering_09.cc b/tests/dofs/dof_renumbering_09.cc
new file mode 100644 (file)
index 0000000..7d3e0dc
--- /dev/null
@@ -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 <deal.II/base/function_lib.h>
+
+#include <deal.II/dofs/dof_accessor.h>
+#include <deal.II/dofs/dof_handler.h>
+#include <deal.II/dofs/dof_renumbering.h>
+
+#include <deal.II/fe/fe_dgp.h>
+#include <deal.II/fe/fe_dgq.h>
+#include <deal.II/fe/fe_q.h>
+#include <deal.II/fe/fe_system.h>
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/tria_iterator.h>
+
+#include <deal.II/lac/vector.h>
+
+#include "../tests.h"
+
+
+
+template <int dim>
+void
+print_dofs(const DoFHandler<dim> &dof, unsigned int level)
+{
+  std::vector<types::global_dof_index> v(dof.get_fe().dofs_per_cell);
+  for (typename DoFHandler<dim>::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 <int dim>
+void
+test()
+{
+  Triangulation<dim> tr(Triangulation<dim>::limit_level_difference_at_vertices);
+  GridGenerator::hyper_cube(tr);
+  tr.refine_global((dim == 1 ? 2 : 1));
+
+  DoFHandler<dim> dof_handler(tr);
+
+  FE_Q<dim> 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 (file)
index 0000000..a82919a
--- /dev/null
@@ -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

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.