From: Wolfgang Bangerth Date: Fri, 30 Aug 2024 14:12:26 +0000 (-0600) Subject: Add the testcase from #17599. X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=856881c0442a6b21eae248fd7361f0a95b744922;p=dealii.git Add the testcase from #17599. --- diff --git a/tests/lac/constraints_make_consistent_in_parallel_06.cc b/tests/lac/constraints_make_consistent_in_parallel_06.cc new file mode 100644 index 0000000000..ab168c1daf --- /dev/null +++ b/tests/lac/constraints_make_consistent_in_parallel_06.cc @@ -0,0 +1,131 @@ +// ------------------------------------------------------------------------ +// +// SPDX-License-Identifier: LGPL-2.1-or-later +// Copyright (C) 2018 - 2023 by the deal.II authors +// +// This file is part of the deal.II library. +// +// Part of the source code is dual licensed under Apache-2.0 WITH +// LLVM-exception OR LGPL-2.1-or-later. Detailed license information +// governing the source code and code contributions can be found in +// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II. +// +// ------------------------------------------------------------------------ + + + +// test AffineConstraints::make_consistent_in_parallel for case +// where constraints need to be combined between different ranks + +#include + +#include +#include + +#include + +#include +#include + +#include + +#include + +#include "../tests.h" + +using namespace dealii; + +void +test() +{ + static const int dim = 3; + static const int spacedim = 3; + + parallel::distributed::Triangulation triangulation( + MPI_COMM_WORLD); + DoFHandler dof_handler(triangulation); + + hp::FECollection fe_collection; + + IndexSet locally_active_dofs; + IndexSet locally_relevant_dofs; + + std::vector repetitions({4, 2, 1}); + Point bottom_left(-2, -1, 0); + Point top_right(2, 1, 1); + + GridGenerator::subdivided_hyper_rectangle(triangulation, + repetitions, + bottom_left, + top_right); + + // hp-refine center part + for (const auto &cell : dof_handler.active_cell_iterators() | + IteratorFilters::LocallyOwnedCell()) + { + cell->set_active_fe_index(1); + + const auto ¢er = cell->center(); + if (std::abs(center[0]) < 1.) + { + if (center[1] > 0.) + cell->set_active_fe_index(2); + else + cell->set_refine_flag(); + } + } + + triangulation.execute_coarsening_and_refinement(); + + // set dofs + for (unsigned int degree = 1; degree <= 10; ++degree) + fe_collection.push_back(FE_Q(degree)); + dof_handler.distribute_dofs(fe_collection); + + const IndexSet &locally_owned_dofs = dof_handler.locally_owned_dofs(); + locally_active_dofs = DoFTools::extract_locally_active_dofs(dof_handler); + locally_relevant_dofs = DoFTools::extract_locally_relevant_dofs(dof_handler); + AffineConstraints constraints(locally_owned_dofs, + locally_relevant_dofs); + + const auto show_constraints_63_64 = [&]() { + deallog << "What process " + << Utilities::MPI::this_mpi_process(MPI_COMM_WORLD) + << " believes about DoF 63:" << std::endl; + for (const auto &c : constraints.get_lines()) + if (c.index == 63) + for (const auto &entry : c.entries) + deallog << " constrained against " << entry.first + << " with weight " << entry.second << std::endl; + deallog << "What process " + << Utilities::MPI::this_mpi_process(MPI_COMM_WORLD) + << " believes about DoF 64:" << std::endl; + for (const auto &c : constraints.get_lines()) + if (c.index == 64) + for (const auto &entry : c.entries) + deallog << " constrained against " << entry.first + << " with weight " << entry.second << std::endl; + }; + + + deallog << "------------- make_hanging_node_constraints():" << std::endl; + DoFTools::make_hanging_node_constraints(dof_handler, constraints); + show_constraints_63_64(); + + deallog << "------------- make_consistent_in_parallel():" << std::endl; + constraints.make_consistent_in_parallel(locally_owned_dofs, + locally_relevant_dofs, + dof_handler.get_communicator()); + show_constraints_63_64(); +} + + + +int +main(int argc, char *argv[]) +{ + Utilities::MPI::MPI_InitFinalize mpi_init(argc, argv); + MPILogInitAll all; + + test(); +} diff --git a/tests/lac/constraints_make_consistent_in_parallel_06.mpirun=2.output b/tests/lac/constraints_make_consistent_in_parallel_06.mpirun=2.output new file mode 100644 index 0000000000..ef5a5e273a --- /dev/null +++ b/tests/lac/constraints_make_consistent_in_parallel_06.mpirun=2.output @@ -0,0 +1,36 @@ + +DEAL:0::------------- make_hanging_node_constraints(): +DEAL:0::What process 0 believes about DoF 63: +DEAL:0:: constrained against 3 with weight 0.323607 +DEAL:0:: constrained against 7 with weight -0.123607 +DEAL:0:: constrained against 19 with weight 0.800000 +DEAL:0::What process 0 believes about DoF 64: +DEAL:0:: constrained against 3 with weight -0.123607 +DEAL:0:: constrained against 7 with weight 0.323607 +DEAL:0:: constrained against 19 with weight 0.800000 +DEAL:0::------------- make_consistent_in_parallel(): +DEAL:0::What process 0 believes about DoF 63: +DEAL:0:: constrained against 3 with weight 0.323607 +DEAL:0:: constrained against 7 with weight -0.123607 +DEAL:0:: constrained against 19 with weight 0.800000 +DEAL:0::What process 0 believes about DoF 64: +DEAL:0:: constrained against 3 with weight -0.123607 +DEAL:0:: constrained against 7 with weight 0.323607 +DEAL:0:: constrained against 19 with weight 0.800000 + +DEAL:1::------------- make_hanging_node_constraints(): +DEAL:1::What process 1 believes about DoF 63: +DEAL:1::What process 1 believes about DoF 64: +DEAL:1:: constrained against 3 with weight -0.447214 +DEAL:1:: constrained against 7 with weight 0.447214 +DEAL:1:: constrained against 63 with weight 1.00000 +DEAL:1::------------- make_consistent_in_parallel(): +DEAL:1::What process 1 believes about DoF 63: +DEAL:1:: constrained against 3 with weight 0.323607 +DEAL:1:: constrained against 7 with weight -0.123607 +DEAL:1:: constrained against 19 with weight 0.800000 +DEAL:1::What process 1 believes about DoF 64: +DEAL:1:: constrained against 3 with weight -0.123607 +DEAL:1:: constrained against 7 with weight 0.323607 +DEAL:1:: constrained against 19 with weight 0.800000 +