]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add additional test with hanging nodes. Fix bug in identification of periodic vertices.
authorMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Wed, 27 Jan 2016 16:39:05 +0000 (17:39 +0100)
committerMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Thu, 28 Jan 2016 08:08:21 +0000 (09:08 +0100)
source/distributed/tria.cc
tests/mpi/distribute_mg_dofs_periodicity_01.cc
tests/mpi/distribute_mg_dofs_periodicity_02.cc [new file with mode: 0644]
tests/mpi/distribute_mg_dofs_periodicity_02.mpirun=7.output [new file with mode: 0644]

index b0f1149c6a96fde86375614dddfb88c0f2ff47e0..fe6b63bd7ac947a0d3e7785ff37dfc73a724228a 100644 (file)
@@ -3187,7 +3187,7 @@ namespace parallel
             const unsigned int vi0 = topological_vertex_numbering[periodic.cell[0]->face(periodic.face_idx[0])->vertex_index(v)];
             const unsigned int vi1 = topological_vertex_numbering[periodic.cell[1]->face(periodic.face_idx[1])->vertex_index(v)];
             const unsigned int min_index = std::min(vi0, vi1);
-            topological_vertex_numbering[vi0] = topological_vertex_numbering[vi1] = min_index;
+            topological_vertex_numbering[periodic.cell[0]->face(periodic.face_idx[0])->vertex_index(v)] = topological_vertex_numbering[periodic.cell[1]->face(periodic.face_idx[1])->vertex_index(v)] = min_index;
           }
       }
 
@@ -3203,9 +3203,10 @@ namespace parallel
         for (unsigned int i=0; i<topological_vertex_numbering.size(); ++i)
           topological_vertex_numbering[i] = i;
         for (unsigned int i=0; i<periodic_face_pairs_level_0.size(); ++i)
+          {
           identify_periodic_vertices_recursively(periodic_face_pairs_level_0[i],
                                                  topological_vertex_numbering);
-
+          }
         // this code is replicated from grid/tria.cc but using an indirection
         // for periodic boundary conditions
         bool continue_iterating = true;
index c6aacdc0b1b0c20c05259fd8b47457b2eb546aac..be6e5d328fab67455b1c0bf3468cd079136dbfb0 100644 (file)
 #include <deal.II/fe/fe_q.h>
 #include <deal.II/distributed/tria.h>
 
-#include <deal.II/numerics/data_out.h>
-
 // This is C++:
 #include <fstream>
 #include <sstream>
 
 
-template <int dim>
-void
-set_owner_recursively(typename Triangulation<dim>::cell_iterator cell,
-                      const int owner,
-                      Vector<double> &vec)
-{
-  if (cell->active())
-    {
-      vec(cell->active_cell_index()) = owner;
-    }
-  else
-    for (unsigned int c=0; c<cell->n_children(); ++c)
-      set_owner_recursively<dim>(cell->child(c), owner, vec);
-}
 
 template <int dim>
 void test()
diff --git a/tests/mpi/distribute_mg_dofs_periodicity_02.cc b/tests/mpi/distribute_mg_dofs_periodicity_02.cc
new file mode 100644 (file)
index 0000000..043d732
--- /dev/null
@@ -0,0 +1,102 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2008 - 2015 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.
+//
+// ---------------------------------------------------------------------
+
+// this tests distribute_mg_dofs on a mesh with periodic boundary
+// conditions. we used to forget to enforce the 2:1 level balance over
+// periodic boundaries in artificial cells, which lead to the case where the
+// coarser MG level did not provide the correct ghost layer. This test
+// explicitly looks into the ghost layer on an adaptive case for processor 0
+
+#include "../tests.h"
+
+#include <deal.II/base/quadrature_lib.h>
+#include <deal.II/base/logstream.h>
+#include <deal.II/base/utilities.h>
+
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/tria_iterator.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_tools.h>
+#include <deal.II/grid/grid_out.h>
+
+#include <deal.II/distributed/tria.h>
+
+// This is C++:
+#include <fstream>
+#include <sstream>
+
+
+
+template <int dim>
+void test()
+{
+  const unsigned int myid = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD);
+
+  parallel::distributed::Triangulation<dim>
+  tria(MPI_COMM_WORLD, dealii::Triangulation<dim>::limit_level_difference_at_vertices,
+       parallel::distributed::Triangulation<dim>::construct_multigrid_hierarchy);
+  GridGenerator::subdivided_hyper_cube (tria, 3);
+  // set periodic boundary conditions in x and z directions
+  for (typename Triangulation<dim>::cell_iterator cell=tria.begin();
+       cell != tria.end(); ++cell)
+    for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
+      if (cell->at_boundary(f))
+        cell->face(f)->set_all_boundary_ids(f+10);
+
+  std::vector<GridTools::PeriodicFacePair<typename Triangulation<dim>::cell_iterator> > periodic_faces;
+  GridTools::collect_periodic_faces(tria, 0+10, 1+10, 0, periodic_faces);
+  GridTools::collect_periodic_faces(tria, 2+10, 3+10, 1, periodic_faces);
+  tria.add_periodicity(periodic_faces);
+  //tria.refine_global(3);
+
+  // adaptively refine into the corner on proc 0
+  for (unsigned int i=0; i<2; ++i)
+    {
+      for (typename Triangulation<dim>::active_cell_iterator cell = tria.begin_active();
+           cell != tria.end(); ++cell)
+        if (cell->at_boundary(0) && cell->at_boundary(2) && cell->is_locally_owned())
+          cell->set_refine_flag();
+      tria.execute_coarsening_and_refinement();
+    }
+
+  for (typename Triangulation<dim>::cell_iterator cell=tria.begin();
+       cell != tria.end(); ++cell)
+    deallog << cell->id().to_string() << " " << cell->level_subdomain_id() << std::endl;
+  if (1)
+    {
+      std::ofstream grid_output (("out"+Utilities::to_string(myid)+".svg").c_str());
+      GridOut grid_out;
+      GridOutFlags::Svg flags;
+      flags.label_level_subdomain_id = true;
+      flags.coloring = GridOutFlags::Svg::level_subdomain_id;
+      flags.convert_level_number_to_height = true;
+      grid_out.set_flags(flags);
+
+      grid_out.write_svg (tria, grid_output);
+    }
+}
+
+
+
+int main (int argc, char *argv[])
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv, 1);
+  MPILogInitAll log;
+
+  test<2>();
+
+  return 0;
+}
diff --git a/tests/mpi/distribute_mg_dofs_periodicity_02.mpirun=7.output b/tests/mpi/distribute_mg_dofs_periodicity_02.mpirun=7.output
new file mode 100644 (file)
index 0000000..b74e1de
--- /dev/null
@@ -0,0 +1,208 @@
+
+DEAL:0::0_0: 0
+DEAL:0::1_0: 2
+DEAL:0::2_0: 3
+DEAL:0::3_0: 2
+DEAL:0::4_0: 2
+DEAL:0::5_0: 4
+DEAL:0::6_0: 4
+DEAL:0::7_0: 5
+DEAL:0::8_0: 6
+DEAL:0::0_1:0 0
+DEAL:0::0_1:1 1
+DEAL:0::0_1:2 1
+DEAL:0::0_1:3 2
+DEAL:0::2_1:0 4294967294
+DEAL:0::2_1:1 3
+DEAL:0::2_1:2 4294967294
+DEAL:0::2_1:3 3
+DEAL:0::6_1:0 4294967294
+DEAL:0::6_1:1 4294967294
+DEAL:0::6_1:2 4
+DEAL:0::6_1:3 4
+DEAL:0::8_1:0 4294967294
+DEAL:0::8_1:1 4294967294
+DEAL:0::8_1:2 4294967294
+DEAL:0::8_1:3 6
+DEAL:0::0_2:00 0
+DEAL:0::0_2:01 0
+DEAL:0::0_2:02 0
+DEAL:0::0_2:03 0
+
+DEAL:1::0_0: 0
+DEAL:1::1_0: 2
+DEAL:1::2_0: 4294967294
+DEAL:1::3_0: 2
+DEAL:1::4_0: 4294967294
+DEAL:1::5_0: 4294967294
+DEAL:1::6_0: 4294967294
+DEAL:1::7_0: 4294967294
+DEAL:1::8_0: 4294967294
+DEAL:1::0_1:0 0
+DEAL:1::0_1:1 1
+DEAL:1::0_1:2 1
+DEAL:1::0_1:3 2
+DEAL:1::2_1:0 4294967294
+DEAL:1::2_1:1 3
+DEAL:1::2_1:2 4294967294
+DEAL:1::2_1:3 3
+DEAL:1::6_1:0 4294967294
+DEAL:1::6_1:1 4294967294
+DEAL:1::6_1:2 4
+DEAL:1::6_1:3 4
+DEAL:1::0_2:00 4294967294
+DEAL:1::0_2:01 0
+DEAL:1::0_2:02 0
+DEAL:1::0_2:03 0
+
+
+DEAL:2::0_0: 0
+DEAL:2::1_0: 2
+DEAL:2::2_0: 3
+DEAL:2::3_0: 2
+DEAL:2::4_0: 2
+DEAL:2::5_0: 4
+DEAL:2::6_0: 4
+DEAL:2::7_0: 5
+DEAL:2::8_0: 6
+DEAL:2::0_1:0 0
+DEAL:2::0_1:1 1
+DEAL:2::0_1:2 1
+DEAL:2::0_1:3 2
+DEAL:2::2_1:0 3
+DEAL:2::2_1:1 4294967294
+DEAL:2::2_1:2 3
+DEAL:2::2_1:3 4294967294
+DEAL:2::6_1:0 4
+DEAL:2::6_1:1 4
+DEAL:2::6_1:2 4294967294
+DEAL:2::6_1:3 4294967294
+DEAL:2::8_1:0 4294967294
+DEAL:2::8_1:1 4294967294
+DEAL:2::8_1:2 4294967294
+DEAL:2::8_1:3 4294967294
+DEAL:2::0_2:00 4294967294
+DEAL:2::0_2:01 4294967294
+DEAL:2::0_2:02 4294967294
+DEAL:2::0_2:03 4294967294
+
+
+DEAL:3::0_0: 0
+DEAL:3::1_0: 2
+DEAL:3::2_0: 3
+DEAL:3::3_0: 2
+DEAL:3::4_0: 2
+DEAL:3::5_0: 4
+DEAL:3::6_0: 4
+DEAL:3::7_0: 5
+DEAL:3::8_0: 6
+DEAL:3::0_1:0 0
+DEAL:3::0_1:1 4294967294
+DEAL:3::0_1:2 1
+DEAL:3::0_1:3 4294967294
+DEAL:3::2_1:0 3
+DEAL:3::2_1:1 3
+DEAL:3::2_1:2 3
+DEAL:3::2_1:3 3
+DEAL:3::6_1:0 4294967294
+DEAL:3::6_1:1 4294967294
+DEAL:3::6_1:2 4
+DEAL:3::6_1:3 4294967294
+DEAL:3::8_1:0 4294967294
+DEAL:3::8_1:1 4294967294
+DEAL:3::8_1:2 6
+DEAL:3::8_1:3 6
+DEAL:3::0_2:00 4294967294
+DEAL:3::0_2:01 4294967294
+DEAL:3::0_2:02 4294967294
+DEAL:3::0_2:03 4294967294
+
+
+DEAL:4::0_0: 0
+DEAL:4::1_0: 2
+DEAL:4::2_0: 3
+DEAL:4::3_0: 2
+DEAL:4::4_0: 2
+DEAL:4::5_0: 4
+DEAL:4::6_0: 4
+DEAL:4::7_0: 5
+DEAL:4::8_0: 6
+DEAL:4::0_1:0 0
+DEAL:4::0_1:1 1
+DEAL:4::0_1:2 4294967294
+DEAL:4::0_1:3 4294967294
+DEAL:4::2_1:0 4294967294
+DEAL:4::2_1:1 3
+DEAL:4::2_1:2 3
+DEAL:4::2_1:3 3
+DEAL:4::6_1:0 4
+DEAL:4::6_1:1 4
+DEAL:4::6_1:2 4
+DEAL:4::6_1:3 4
+DEAL:4::8_1:0 6
+DEAL:4::8_1:1 6
+DEAL:4::8_1:2 4294967294
+DEAL:4::8_1:3 6
+DEAL:4::0_2:00 4294967294
+DEAL:4::0_2:01 4294967294
+DEAL:4::0_2:02 4294967294
+DEAL:4::0_2:03 4294967294
+
+
+DEAL:5::0_0: 0
+DEAL:5::1_0: 2
+DEAL:5::2_0: 3
+DEAL:5::3_0: 2
+DEAL:5::4_0: 2
+DEAL:5::5_0: 4
+DEAL:5::6_0: 4
+DEAL:5::7_0: 5
+DEAL:5::8_0: 6
+DEAL:5::0_1:0 4294967294
+DEAL:5::0_1:1 4294967294
+DEAL:5::0_1:2 4294967294
+DEAL:5::0_1:3 4294967294
+DEAL:5::2_1:0 4294967294
+DEAL:5::2_1:1 4294967294
+DEAL:5::2_1:2 4294967294
+DEAL:5::2_1:3 4294967294
+DEAL:5::6_1:0 4294967294
+DEAL:5::6_1:1 4
+DEAL:5::6_1:2 4294967294
+DEAL:5::6_1:3 4
+DEAL:5::8_1:0 6
+DEAL:5::8_1:1 4294967294
+DEAL:5::8_1:2 6
+DEAL:5::8_1:3 4294967294
+
+
+DEAL:6::0_0: 0
+DEAL:6::1_0: 2
+DEAL:6::2_0: 3
+DEAL:6::3_0: 2
+DEAL:6::4_0: 2
+DEAL:6::5_0: 4
+DEAL:6::6_0: 4
+DEAL:6::7_0: 5
+DEAL:6::8_0: 6
+DEAL:6::0_1:0 0
+DEAL:6::0_1:1 4294967294
+DEAL:6::0_1:2 4294967294
+DEAL:6::0_1:3 4294967294
+DEAL:6::2_1:0 3
+DEAL:6::2_1:1 3
+DEAL:6::2_1:2 4294967294
+DEAL:6::2_1:3 4294967294
+DEAL:6::6_1:0 4
+DEAL:6::6_1:1 4294967294
+DEAL:6::6_1:2 4
+DEAL:6::6_1:3 4294967294
+DEAL:6::8_1:0 6
+DEAL:6::8_1:1 6
+DEAL:6::8_1:2 6
+DEAL:6::8_1:3 6
+DEAL:6::0_2:00 4294967294
+DEAL:6::0_2:01 4294967294
+DEAL:6::0_2:02 4294967294
+DEAL:6::0_2:03 4294967294
+

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.