From 42ec1be70205b7b7da1e83ab326a82ea6a9b5036 Mon Sep 17 00:00:00 2001 From: Marc Fehling Date: Thu, 20 Feb 2020 13:50:05 +0100 Subject: [PATCH] New test: p::d::tria + hp::dofh + hanging_node_constraints. --- ...s.cc => hp_hanging_node_constraints_01.cc} | 0 ...raints_01.with_p4est=true.mpirun=2.output} | 0 tests/mpi/hp_hanging_node_constraints_02.cc | 123 ++++++++++++++++++ ...traints_02.with_p4est=true.mpirun=2.output | 7 + 4 files changed, 130 insertions(+) rename tests/mpi/{hp_hanging_node_constraints.cc => hp_hanging_node_constraints_01.cc} (100%) rename tests/mpi/{hp_hanging_node_constraints.with_p4est=true.mpirun=2.output => hp_hanging_node_constraints_01.with_p4est=true.mpirun=2.output} (100%) create mode 100644 tests/mpi/hp_hanging_node_constraints_02.cc create mode 100644 tests/mpi/hp_hanging_node_constraints_02.with_p4est=true.mpirun=2.output diff --git a/tests/mpi/hp_hanging_node_constraints.cc b/tests/mpi/hp_hanging_node_constraints_01.cc similarity index 100% rename from tests/mpi/hp_hanging_node_constraints.cc rename to tests/mpi/hp_hanging_node_constraints_01.cc diff --git a/tests/mpi/hp_hanging_node_constraints.with_p4est=true.mpirun=2.output b/tests/mpi/hp_hanging_node_constraints_01.with_p4est=true.mpirun=2.output similarity index 100% rename from tests/mpi/hp_hanging_node_constraints.with_p4est=true.mpirun=2.output rename to tests/mpi/hp_hanging_node_constraints_01.with_p4est=true.mpirun=2.output diff --git a/tests/mpi/hp_hanging_node_constraints_02.cc b/tests/mpi/hp_hanging_node_constraints_02.cc new file mode 100644 index 0000000000..e14b5cde85 --- /dev/null +++ b/tests/mpi/hp_hanging_node_constraints_02.cc @@ -0,0 +1,123 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2020 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. +// +// see: https://github.com/dealii/dealii/issues/9517 + + +#include + +#include + +#include + +#include +#include +#include + +#include + +#include +#include + +#include + +#include "../tests.h" + + +template +void +test() +{ + // setup triangulation + parallel::distributed::Triangulation triangulation( + MPI_COMM_WORLD, + typename Triangulation::MeshSmoothing( + Triangulation::smoothing_on_refinement | + Triangulation::smoothing_on_coarsening)); + GridGenerator::hyper_rectangle( + triangulation, + (dim == 3 ? Point(0.0, 0.0, 0.0) : Point(0.0, 0.0)), + (dim == 3 ? Point(1.0, 1.0, 1.0) : Point(1.0, 1.0)), + true); + triangulation.refine_global(3); + + // Refine all the cells with y < 0.5 twice + for (unsigned int n_ref = 0; n_ref < 2; ++n_ref) + { + for (const auto &cell : triangulation.active_cell_iterators()) + if (cell->is_locally_owned() && cell->center()(1) < 0.5) + cell->set_refine_flag(); + + triangulation.prepare_coarsening_and_refinement(); + triangulation.execute_coarsening_and_refinement(); + } + + // setup finite elements + FESystem void_fe(FE_Nothing(), dim); + FESystem solid_fe(FE_Q(1), dim); + + hp::FECollection fe_collection; + fe_collection.push_back(void_fe); + fe_collection.push_back(solid_fe); + + hp::DoFHandler dof_handler(triangulation); + dof_handler.set_fe(fe_collection); + + // Assign void_fe to all the cells with x < 0.5 + for (const auto &cell : dof_handler.active_cell_iterators()) + if (cell->is_locally_owned()) + { + if (cell->center()(0) < 0.5) + cell->set_active_fe_index(0); + else + cell->set_active_fe_index(1); + } + + dof_handler.distribute_dofs(fe_collection); + + // make constraints + IndexSet locally_relevant_dofs; + DoFTools::extract_locally_relevant_dofs(dof_handler, locally_relevant_dofs); + + AffineConstraints hanging_node_constraints; + hanging_node_constraints.reinit(locally_relevant_dofs); + DoFTools::make_hanging_node_constraints(dof_handler, + hanging_node_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_02.with_p4est=true.mpirun=2.output b/tests/mpi/hp_hanging_node_constraints_02.with_p4est=true.mpirun=2.output new file mode 100644 index 0000000000..4aa906268b --- /dev/null +++ b/tests/mpi/hp_hanging_node_constraints_02.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