From: Wolfgang Bangerth Date: Thu, 5 Oct 2023 20:59:23 +0000 (-0600) Subject: Add a test. X-Git-Tag: relicensing~425^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=e69c0282c03bbfdda6d88be27c7cbde5fe9d32c4;p=dealii.git Add a test. --- diff --git a/tests/mpi/affine_constraints_get_view_01.cc b/tests/mpi/affine_constraints_get_view_01.cc new file mode 100644 index 0000000000..def3035761 --- /dev/null +++ b/tests/mpi/affine_constraints_get_view_01.cc @@ -0,0 +1,210 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2009 - 2018 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 AffineConstraints.get_view() +// +// The way this test works is by setting up a Stokes system with a +// Q2xQ1 element on an adaptively refined mesh for which we compute +// the corresponding constraints (which will contain constraints for +// both the velocity and pressure degrees of freedom). Then we create +// a random vector to which we apply distribute(). This is all +// existing functionality. +// +// Next, we create views into the constraints object for both the +// velocity variable only and the pressure variable, and apply these +// views to the individual blocks of the (original) random +// vector. This must result in the same output vector as above. + +#include +#include + +#include +#include +#include +#include + +#include +#include + +#include +#include +#include + +#include +#include + +#include + +#include "../tests.h" + + + +template +void +test() +{ + deallog << "dim=" << dim << std::endl; + + const unsigned int myid = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD); + const unsigned int n_processes = + Utilities::MPI::n_mpi_processes(MPI_COMM_WORLD); + + // Set up of triangulation and dof_handler + + parallel::distributed::Triangulation tr(MPI_COMM_WORLD); + + const double R0 = 6371000. - 2890000.; + const double R1 = 6371000. - 35000.; + GridGenerator::hyper_shell(tr, Point(), R0, R1, 12, true); + + tr.refine_global(1); + for (unsigned int step = 0; step < 5; ++step) + { + typename Triangulation::active_cell_iterator cell = + tr.begin_active(), + endc = tr.end(); + + for (; cell != endc; ++cell) + if (Testing::rand() % 42 == 1) + cell->set_refine_flag(); + + tr.execute_coarsening_and_refinement(); + } + deallog << "Cells on process #" << myid << ": " << tr.n_active_cells() + << std::endl; + + FESystem fe(FE_Q(2), dim, FE_Q(1), 1); + DoFHandler dof_handler(tr); + + dof_handler.distribute_dofs(fe); + + std::vector stokes_sub_blocks(dim + 1, 0); + stokes_sub_blocks[dim] = 1; + DoFRenumbering::component_wise(dof_handler, stokes_sub_blocks); + + deallog << "Total number of DoFs: " << dof_handler.n_dofs() << std::endl; + + const std::vector stokes_dofs_per_block = + DoFTools::count_dofs_per_fe_block(dof_handler, stokes_sub_blocks); + + const types::global_dof_index n_u = stokes_dofs_per_block[0], + n_p = stokes_dofs_per_block[1]; + + const IndexSet &stokes_locally_owned_index_set = + dof_handler.locally_owned_dofs(); + const IndexSet stokes_locally_relevant_set = + DoFTools::extract_locally_relevant_dofs(dof_handler); + + std::vector stokes_partitioning; + stokes_partitioning.push_back( + stokes_locally_owned_index_set.get_view(0, n_u)); + stokes_partitioning.push_back( + stokes_locally_owned_index_set.get_view(n_u, n_u + n_p)); + + std::vector stokes_relevant_partitioning; + stokes_relevant_partitioning.push_back( + stokes_locally_relevant_set.get_view(0, n_u)); + stokes_relevant_partitioning.push_back( + stokes_locally_relevant_set.get_view(n_u, n_u + n_p)); + + // Set up a constraints object from hanging node constraints and + // boundary conditions. + AffineConstraints stokes_constraints(stokes_locally_owned_index_set, + stokes_locally_relevant_set); + + DoFTools::make_hanging_node_constraints(dof_handler, stokes_constraints); + + const FEValuesExtractors::Vector velocity_components(0); + VectorTools::interpolate_boundary_values( + dof_handler, + 0, + Functions::ZeroFunction(dim + 1), + stokes_constraints, + fe.component_mask(velocity_components)); + + std::set no_normal_flux_boundaries; + no_normal_flux_boundaries.insert(1); + VectorTools::compute_no_normal_flux_constraints(dof_handler, + 0, + no_normal_flux_boundaries, + stokes_constraints); + stokes_constraints.close(); + deallog << "Number of constrains: " << stokes_constraints.n_constraints() + << std::endl; + + // Now set up a block vector for this situation. Fill it with random numbers. + TrilinosWrappers::MPI::BlockVector random_vector; + random_vector.reinit(stokes_partitioning, MPI_COMM_WORLD); + for (const unsigned int i : stokes_locally_owned_index_set) + random_vector(i) = random_value(-100, 100); + random_vector.compress(VectorOperation::insert); + + + // Check 1: Call distribute() on a vector that contains + // everything, with an AffineConstraint object that contains all + // constraints. + TrilinosWrappers::MPI::BlockVector check_vector_1; + { + check_vector_1.reinit(stokes_partitioning, MPI_COMM_WORLD); + check_vector_1 = random_vector; + + stokes_constraints.distribute(check_vector_1); + } + + // Check 2: Call distribute() on the two blocks individually, with + // AffineConstraint objects that only contain views. + TrilinosWrappers::MPI::BlockVector check_vector_2; + { + check_vector_2.reinit(stokes_partitioning, MPI_COMM_WORLD); + check_vector_2 = random_vector; + + + IndexSet velocity_dofs(dof_handler.n_dofs()); + velocity_dofs.add_range(0, n_u); + velocity_dofs.compress(); + + const AffineConstraints velocity_constraints = + stokes_constraints.get_view(velocity_dofs); + + + IndexSet pressure_dofs(dof_handler.n_dofs()); + pressure_dofs.add_range(n_u, n_u + n_p); + pressure_dofs.compress(); + + const AffineConstraints pressure_constraints = + stokes_constraints.get_view(pressure_dofs); + + velocity_constraints.distribute(check_vector_2.block(0)); + pressure_constraints.distribute(check_vector_2.block(1)); + } + + // Now the assertion part: These vectors are the same: + Assert(check_vector_1 == check_vector_2, ExcInternalError()); + + deallog << "OK" << std::endl; +} + + +int +main(int argc, char *argv[]) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + mpi_initlog(); + + test<2>(); + test<3>(); +} diff --git a/tests/mpi/affine_constraints_get_view_01.with_trilinos=true.mpirun=1.output b/tests/mpi/affine_constraints_get_view_01.with_trilinos=true.mpirun=1.output new file mode 100644 index 0000000000..8369f6974b --- /dev/null +++ b/tests/mpi/affine_constraints_get_view_01.with_trilinos=true.mpirun=1.output @@ -0,0 +1,11 @@ + +DEAL::dim=2 +DEAL::Cells on process #0: 72 +DEAL::Total number of DoFs: 878 +DEAL::Number of constrains: 306 +DEAL::OK +DEAL::dim=3 +DEAL::Cells on process #0: 509 +DEAL::Total number of DoFs: 18124 +DEAL::Number of constrains: 7698 +DEAL::OK diff --git a/tests/mpi/affine_constraints_get_view_01.with_trilinos=true.mpirun=2.output b/tests/mpi/affine_constraints_get_view_01.with_trilinos=true.mpirun=2.output new file mode 100644 index 0000000000..9cd86c8a47 --- /dev/null +++ b/tests/mpi/affine_constraints_get_view_01.with_trilinos=true.mpirun=2.output @@ -0,0 +1,11 @@ + +DEAL::dim=2 +DEAL::Cells on process #0: 51 +DEAL::Total number of DoFs: 853 +DEAL::Number of constrains: 170 +DEAL::OK +DEAL::dim=3 +DEAL::Cells on process #0: 481 +DEAL::Total number of DoFs: 20074 +DEAL::Number of constrains: 6280 +DEAL::OK diff --git a/tests/mpi/affine_constraints_get_view_01.with_trilinos=true.mpirun=4.output b/tests/mpi/affine_constraints_get_view_01.with_trilinos=true.mpirun=4.output new file mode 100644 index 0000000000..6ea1c7537b --- /dev/null +++ b/tests/mpi/affine_constraints_get_view_01.with_trilinos=true.mpirun=4.output @@ -0,0 +1,11 @@ + +DEAL::dim=2 +DEAL::Cells on process #0: 45 +DEAL::Total number of DoFs: 903 +DEAL::Number of constrains: 123 +DEAL::OK +DEAL::dim=3 +DEAL::Cells on process #0: 215 +DEAL::Total number of DoFs: 9429 +DEAL::Number of constrains: 2582 +DEAL::OK