]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Hanging nodes: Loosen assertion for artificial cells. 8365/head
authormarcfehling <marc.fehling@gmx.net>
Wed, 10 Jul 2019 16:12:29 +0000 (18:12 +0200)
committerMarc Fehling <marc.fehling@gmx.net>
Wed, 10 Jul 2019 16:29:13 +0000 (18:29 +0200)
source/dofs/dof_tools_constraints.cc
tests/mpi/hp_hanging_node_constraints.cc [new file with mode: 0644]
tests/mpi/hp_hanging_node_constraints.with_p4est=true.mpirun=2.output [new file with mode: 0644]

index 32d01d2ee47b33e8cc38da89aff5a6562717e774..ab6b067f96ebaa7605f8cc733cd6893d2f5eaee4 100644 (file)
@@ -840,8 +840,10 @@ namespace DoFTools
                        ExcInternalError());
                 for (unsigned int c = 0; c < cell->face(face)->n_children();
                      ++c)
-                  AssertDimension(
-                    cell->face(face)->child(c)->n_active_fe_indices(), 1);
+                  if (!cell->neighbor_child_on_subface(face, c)
+                         ->is_artificial())
+                    AssertDimension(
+                      cell->face(face)->child(c)->n_active_fe_indices(), 1);
 
                 // right now, all that is implemented is the case that both
                 // sides use the same fe, and not only that but also that all
@@ -1115,8 +1117,11 @@ namespace DoFTools
                        ExcInternalError());
                 for (unsigned int c = 0; c < cell->face(face)->n_children();
                      ++c)
-                  Assert(cell->face(face)->child(c)->n_active_fe_indices() == 1,
-                         ExcInternalError());
+                  if (!cell->neighbor_child_on_subface(face, c)
+                         ->is_artificial())
+                    Assert(cell->face(face)->child(c)->n_active_fe_indices() ==
+                             1,
+                           ExcInternalError());
 
                 // first find out whether we can constrain each of the subfaces
                 // to the mother face. in the lingo of the hp paper, this would
diff --git a/tests/mpi/hp_hanging_node_constraints.cc b/tests/mpi/hp_hanging_node_constraints.cc
new file mode 100644 (file)
index 0000000..d8e0512
--- /dev/null
@@ -0,0 +1,102 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2019 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.
+//
+// ---------------------------------------------------------------------
+
+
+
+// this test failed at some stage:
+// DoFTools::make_hanging_node_constraints() has not worked with an
+// hp::DoFHandler on ghost cells that neighbor artificial cells.
+//
+// Geometry in 2D divided on two subdomains:
+// Subdomain0:          Subdomain1:
+// +-------+-------+    +-------+-------+
+// |       |       |    |       |       |
+// |   g   |   g   |    |   o   |   o   |
+// |       |       |    |       |       |
+// +---+---+---+!!!+    +---+---+---+---+
+// | o | o | g | a |    | g | g | o | o |
+// +---+---+---+---+    +---+---+---+---+
+// | o | o | g | a |    | a | g | o | o |
+// +---+---+---+---+    +---+---+---+---+
+
+
+#include <deal.II/base/index_set.h>
+
+#include <deal.II/distributed/tria.h>
+
+#include <deal.II/dofs/dof_tools.h>
+
+#include <deal.II/fe/fe_q.h>
+
+#include <deal.II/grid/grid_generator.h>
+
+#include <deal.II/hp/dof_handler.h>
+#include <deal.II/hp/fe_collection.h>
+
+#include <deal.II/lac/affine_constraints.h>
+
+#include "../tests.h"
+
+
+template <int dim>
+void
+test()
+{
+  // setup triangulation
+  parallel::distributed::Triangulation<dim> tria(MPI_COMM_WORLD);
+  GridGenerator::hyper_cube(tria);
+  tria.refine_global(1);
+
+  // refine half of all cells
+  unsigned int i = 0;
+  for (const auto &cell : tria.active_cell_iterators())
+    if (i++ < .5 * tria.n_active_cells())
+      cell->set_refine_flag();
+  tria.execute_coarsening_and_refinement();
+
+  // setup finite elemets
+  hp::FECollection<dim> fes;
+  fes.push_back(FE_Q<dim>(1));
+
+  hp::DoFHandler<dim> dh(tria);
+  dh.distribute_dofs(fes);
+
+  // make constraints
+  IndexSet locally_relevant_dofs;
+  DoFTools::extract_locally_relevant_dofs(dh, locally_relevant_dofs);
+
+  AffineConstraints<double> constraints;
+  constraints.reinit(locally_relevant_dofs);
+  DoFTools::make_hanging_node_constraints(dh, constraints);
+
+  deallog << "OK" << std::endl;
+}
+
+
+int
+main(int argc, char *argv[])
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
+
+  MPILogInitAll all;
+
+  deallog.push("2d");
+  test<2>();
+  deallog.pop();
+  MPI_Barrier(MPI_COMM_WORLD);
+  deallog.push("3d");
+  test<3>();
+  deallog.pop();
+}
diff --git a/tests/mpi/hp_hanging_node_constraints.with_p4est=true.mpirun=2.output b/tests/mpi/hp_hanging_node_constraints.with_p4est=true.mpirun=2.output
new file mode 100644 (file)
index 0000000..4aa9062
--- /dev/null
@@ -0,0 +1,7 @@
+
+DEAL:0:2d::OK
+DEAL:0:3d::OK
+
+DEAL:1:2d::OK
+DEAL:1:3d::OK
+

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.