From d023baf25c8b280d2e627eeb44f6a0c70ee5514e Mon Sep 17 00:00:00 2001 From: Marc Fehling Date: Mon, 12 Aug 2019 16:32:03 -0600 Subject: [PATCH] New test: Chains of constraints on parallel hp-adaptive applications. --- tests/mpi/hp_chains_of_constraints_01.cc | 166 ++++++++++++++++++ ...traints_01.with_p4est=true.mpirun=1.output | 9 + ...traints_01.with_p4est=true.mpirun=4.output | 39 ++++ 3 files changed, 214 insertions(+) create mode 100644 tests/mpi/hp_chains_of_constraints_01.cc create mode 100644 tests/mpi/hp_chains_of_constraints_01.with_p4est=true.mpirun=1.output create mode 100644 tests/mpi/hp_chains_of_constraints_01.with_p4est=true.mpirun=4.output diff --git a/tests/mpi/hp_chains_of_constraints_01.cc b/tests/mpi/hp_chains_of_constraints_01.cc new file mode 100644 index 0000000000..eba9f5075e --- /dev/null +++ b/tests/mpi/hp_chains_of_constraints_01.cc @@ -0,0 +1,166 @@ +// --------------------------------------------------------------------- +// +// 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. +// +// --------------------------------------------------------------------- + + + +// Check whether chains of constraints traverse over dofs on faces, +// and if this is an issue in parallel applications. +// +// Domain in 2D/3D with no normal flux constraints on the boundary: +// +---------+ +// | | +// | FE_Q(1) | +// | | +// +---------+---------+---------+ +// | | | | +// | FE_Q(1) | FE_Q(2) | FE_Q(1) | +// | | | | +// +---------+---------+---------+ + + +#include + +#include + +#include +#include + +#include + +#include +#include + +#include + +#include + +#include + +#include "../tests.h" + + +template +void +test(const unsigned int degree_center, + const unsigned int degree_other, + const bool enable_no_flux_constraints = false, + const bool print_constraints = false) +{ + Assert(dim > 1, ExcNotImplemented()); + + // ------ setup ------ + // prepare Triangulation as a hyper_cross + parallel::distributed::Triangulation tria(MPI_COMM_WORLD); + + std::vector sizes; + sizes.insert(sizes.end(), {1, 1}); // extension in in x-direction + sizes.insert(sizes.end(), {0, 1}); // extension in in y-direction + if (dim > 2) + sizes.insert(sizes.end(), {0, 0}); // extension in in z-direction + + GridGenerator::hyper_cross(tria, sizes); + + // prepare FECollection with arbitrary number of entries + hp::FECollection fe_collection; + fe_collection.push_back(FESystem(FE_Q(degree_other), dim)); + fe_collection.push_back(FESystem(FE_Q(degree_center), dim)); + + // prepare DoFHandler + hp::DoFHandler dh(tria); + + // center cell has ID 0 + const auto ¢er = dh.begin_active(); + if (center->is_locally_owned() && center->id().to_string() == "0_0:") + { + center->set_active_fe_index(1); + +#ifdef DEBUG + // verify if our scenario is initialized correctly + // by checking the number of neighbors of the center cell + unsigned int n_neighbors = 0; + for (unsigned int i = 0; i < GeometryInfo::faces_per_cell; ++i) + if (static_cast(center->neighbor_index(i)) != + numbers::invalid_unsigned_int) + ++n_neighbors; + Assert(n_neighbors == 3, ExcInternalError()); +#endif + } + + dh.distribute_dofs(fe_collection); + + // ---- constrain ---- + IndexSet locally_relevant_dofs; + DoFTools::extract_locally_relevant_dofs(dh, locally_relevant_dofs); + + AffineConstraints constraints; + constraints.clear(); + constraints.reinit(locally_relevant_dofs); + + DoFTools::make_hanging_node_constraints(dh, constraints); + + if (enable_no_flux_constraints) + VectorTools::compute_no_normal_flux_constraints( + dh, + 0, /*first component*/ + std::set{0}, + constraints); + + constraints.close(); + + // ------ verify ----- + IndexSet locally_active_dofs; + DoFTools::extract_locally_active_dofs(dh, locally_active_dofs); + + if (print_constraints) + { + deallog << "constraints:" << std::endl; + constraints.print(deallog.get_file_stream()); + } + deallog << "consistent? " + << constraints.is_consistent_in_parallel( + dh.locally_owned_dofs_per_processor(), + locally_active_dofs, + MPI_COMM_WORLD, + true) + << std::endl; + + deallog << "OK" << std::endl; +} + + +int +main(int argc, char *argv[]) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + MPILogInitAll log; + + deallog.push("2d"); + deallog.push("hanging_nodes"); + test<2>(2, 1, false); + deallog.pop(); + deallog.push("hanging_nodes+no_flux"); + test<2>(2, 1, true); + deallog.pop(); + deallog.pop(); + + deallog.push("3d"); + deallog.push("hanging_nodes"); + test<3>(2, 1, false); + deallog.pop(); + deallog.push("hanging_nodes+no_flux"); + test<3>(2, 1, true); + deallog.pop(); + deallog.pop(); +} diff --git a/tests/mpi/hp_chains_of_constraints_01.with_p4est=true.mpirun=1.output b/tests/mpi/hp_chains_of_constraints_01.with_p4est=true.mpirun=1.output new file mode 100644 index 0000000000..812c0d845f --- /dev/null +++ b/tests/mpi/hp_chains_of_constraints_01.with_p4est=true.mpirun=1.output @@ -0,0 +1,9 @@ + +DEAL:0:2d:hanging_nodes::consistent? 1 +DEAL:0:2d:hanging_nodes::OK +DEAL:0:2d:hanging_nodes+no_flux::consistent? 1 +DEAL:0:2d:hanging_nodes+no_flux::OK +DEAL:0:3d:hanging_nodes::consistent? 1 +DEAL:0:3d:hanging_nodes::OK +DEAL:0:3d:hanging_nodes+no_flux::consistent? 1 +DEAL:0:3d:hanging_nodes+no_flux::OK diff --git a/tests/mpi/hp_chains_of_constraints_01.with_p4est=true.mpirun=4.output b/tests/mpi/hp_chains_of_constraints_01.with_p4est=true.mpirun=4.output new file mode 100644 index 0000000000..b10aee54af --- /dev/null +++ b/tests/mpi/hp_chains_of_constraints_01.with_p4est=true.mpirun=4.output @@ -0,0 +1,39 @@ + +DEAL:0:2d:hanging_nodes::consistent? 1 +DEAL:0:2d:hanging_nodes::OK +DEAL:0:2d:hanging_nodes+no_flux::consistent? 1 +DEAL:0:2d:hanging_nodes+no_flux::OK +DEAL:0:3d:hanging_nodes::consistent? 1 +DEAL:0:3d:hanging_nodes::OK +DEAL:0:3d:hanging_nodes+no_flux::consistent? 1 +DEAL:0:3d:hanging_nodes+no_flux::OK + +DEAL:1:2d:hanging_nodes::consistent? 1 +DEAL:1:2d:hanging_nodes::OK +DEAL:1:2d:hanging_nodes+no_flux::consistent? 1 +DEAL:1:2d:hanging_nodes+no_flux::OK +DEAL:1:3d:hanging_nodes::consistent? 1 +DEAL:1:3d:hanging_nodes::OK +DEAL:1:3d:hanging_nodes+no_flux::consistent? 1 +DEAL:1:3d:hanging_nodes+no_flux::OK + + +DEAL:2:2d:hanging_nodes::consistent? 1 +DEAL:2:2d:hanging_nodes::OK +DEAL:2:2d:hanging_nodes+no_flux::consistent? 1 +DEAL:2:2d:hanging_nodes+no_flux::OK +DEAL:2:3d:hanging_nodes::consistent? 1 +DEAL:2:3d:hanging_nodes::OK +DEAL:2:3d:hanging_nodes+no_flux::consistent? 1 +DEAL:2:3d:hanging_nodes+no_flux::OK + + +DEAL:3:2d:hanging_nodes::consistent? 1 +DEAL:3:2d:hanging_nodes::OK +DEAL:3:2d:hanging_nodes+no_flux::consistent? 1 +DEAL:3:2d:hanging_nodes+no_flux::OK +DEAL:3:3d:hanging_nodes::consistent? 1 +DEAL:3:3d:hanging_nodes::OK +DEAL:3:3d:hanging_nodes+no_flux::consistent? 1 +DEAL:3:3d:hanging_nodes+no_flux::OK + -- 2.39.5