From a9cab0d741bbc41e7e8144e4d7b67b7ae7eaf343 Mon Sep 17 00:00:00 2001 From: marcfehling Date: Wed, 10 Jul 2019 18:12:29 +0200 Subject: [PATCH] Hanging nodes: Loosen assertion for artificial cells. --- source/dofs/dof_tools_constraints.cc | 13 ++- tests/mpi/hp_hanging_node_constraints.cc | 102 ++++++++++++++++++ ...onstraints.with_p4est=true.mpirun=2.output | 7 ++ 3 files changed, 118 insertions(+), 4 deletions(-) create mode 100644 tests/mpi/hp_hanging_node_constraints.cc create mode 100644 tests/mpi/hp_hanging_node_constraints.with_p4est=true.mpirun=2.output diff --git a/source/dofs/dof_tools_constraints.cc b/source/dofs/dof_tools_constraints.cc index 32d01d2ee4..ab6b067f96 100644 --- a/source/dofs/dof_tools_constraints.cc +++ b/source/dofs/dof_tools_constraints.cc @@ -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 index 0000000000..d8e051231e --- /dev/null +++ b/tests/mpi/hp_hanging_node_constraints.cc @@ -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 + +#include + +#include + +#include + +#include + +#include +#include + +#include + +#include "../tests.h" + + +template +void +test() +{ + // setup triangulation + parallel::distributed::Triangulation 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 fes; + fes.push_back(FE_Q(1)); + + hp::DoFHandler dh(tria); + dh.distribute_dofs(fes); + + // make constraints + IndexSet locally_relevant_dofs; + DoFTools::extract_locally_relevant_dofs(dh, locally_relevant_dofs); + + AffineConstraints 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 index 0000000000..4aa906268b --- /dev/null +++ b/tests/mpi/hp_hanging_node_constraints.with_p4est=true.mpirun=2.output @@ -0,0 +1,7 @@ + +DEAL:0:2d::OK +DEAL:0:3d::OK + +DEAL:1:2d::OK +DEAL:1:3d::OK + -- 2.39.5