From 5252bdb341e04ace0ad9f0bd31e2894b4c98a912 Mon Sep 17 00:00:00 2001 From: Marc Fehling Date: Fri, 2 Jun 2023 21:35:00 -0600 Subject: [PATCH] Choose to initialize BlockVector with ghost elements using reinit(partitioners). --- doc/news/changes/minor/20230602Fehling | 6 ++ .../deal.II/lac/la_parallel_block_vector.h | 13 +++ .../lac/la_parallel_block_vector.templates.h | 13 +++ include/deal.II/lac/petsc_block_vector.h | 33 ++++++ .../lac/trilinos_parallel_block_vector.h | 16 +++ source/lac/trilinos_block_vector.cc | 23 ++++ ....cc => parallel_block_vector_reinit_01.cc} | 4 +- ...l_block_vector_reinit_01.mpirun=10.output} | 2 +- .../petsc/parallel_block_vector_reinit_01.cc | 101 ++++++++++++++++++ ...el_block_vector_reinit_01.mpirun=10.output | 5 + .../parallel_block_vector_reinit_01.cc | 100 +++++++++++++++++ ...el_block_vector_reinit_01.mpirun=10.output | 5 + 12 files changed, 318 insertions(+), 3 deletions(-) create mode 100644 doc/news/changes/minor/20230602Fehling rename tests/mpi/{parallel_block_vector_14.cc => parallel_block_vector_reinit_01.cc} (96%) rename tests/mpi/{parallel_block_vector_14.mpirun=10.output => parallel_block_vector_reinit_01.mpirun=10.output} (71%) create mode 100644 tests/petsc/parallel_block_vector_reinit_01.cc create mode 100644 tests/petsc/parallel_block_vector_reinit_01.mpirun=10.output create mode 100644 tests/trilinos/parallel_block_vector_reinit_01.cc create mode 100644 tests/trilinos/parallel_block_vector_reinit_01.mpirun=10.output diff --git a/doc/news/changes/minor/20230602Fehling b/doc/news/changes/minor/20230602Fehling new file mode 100644 index 0000000000..1b1fe6f2f5 --- /dev/null +++ b/doc/news/changes/minor/20230602Fehling @@ -0,0 +1,6 @@ +New: All parallel BlockVector classes now have a reinit function that takes +collections of Utilities::MPI::Partitioner objects as an argument. This +affects TrilinosWrappers::MPI::BlockVector, PETScWrappers::MPI::BlockVector, +and LinearAlgebra::distributed::BlockVector. The interface is compatible. +
+(Marc Fehling, 2023/06/02) diff --git a/include/deal.II/lac/la_parallel_block_vector.h b/include/deal.II/lac/la_parallel_block_vector.h index 321ea64071..08b44c824f 100644 --- a/include/deal.II/lac/la_parallel_block_vector.h +++ b/include/deal.II/lac/la_parallel_block_vector.h @@ -363,6 +363,19 @@ namespace LinearAlgebra & partitioners, const MPI_Comm &comm_sm = MPI_COMM_SELF); + /** + * This function exists purely for reasons of compatibility with the + * PETScWrappers::MPI::Vector and TrilinosWrappers::MPI::Vector classes. + * + * It calls the function above, and ignores the parameter @p make_ghosted. + */ + void + reinit( + const std::vector> + & partitioners, + const bool make_ghosted, + const MPI_Comm &comm_sm = MPI_COMM_SELF); + /** * This function copies the data that has accumulated in the data buffer * for ghost indices to the owning processor. For the meaning of the diff --git a/include/deal.II/lac/la_parallel_block_vector.templates.h b/include/deal.II/lac/la_parallel_block_vector.templates.h index da7ef55ac0..990c2950d3 100644 --- a/include/deal.II/lac/la_parallel_block_vector.templates.h +++ b/include/deal.II/lac/la_parallel_block_vector.templates.h @@ -220,6 +220,19 @@ namespace LinearAlgebra + template + void + BlockVector::reinit( + const std::vector> + &partitioners, + const bool /*make_ghosted*/, + const MPI_Comm &comm_sm) + { + reinit(partitioners, comm_sm); + } + + + template BlockVector & BlockVector::operator=(const value_type s) diff --git a/include/deal.II/lac/petsc_block_vector.h b/include/deal.II/lac/petsc_block_vector.h index 94e4a1aa76..2c7694da10 100644 --- a/include/deal.II/lac/petsc_block_vector.h +++ b/include/deal.II/lac/petsc_block_vector.h @@ -245,6 +245,19 @@ namespace PETScWrappers const std::vector &ghost_entries, const MPI_Comm communicator); + /** + * Initialize each block given to each parallel partitioning described in + * @p partitioners. + * + * You can decide whether your vector will contain ghost elements with + * @p make_ghosted. + */ + void + reinit( + const std::vector> + & partitioners, + const bool make_ghosted = true); + /** * This function collects the sizes of the sub-objects and stores them * in internal arrays, in order to be able to relay global indices into @@ -569,6 +582,26 @@ namespace PETScWrappers + inline void + BlockVector::reinit( + const std::vector> + & partitioners, + const bool make_ghosted) + { + // update the number of blocks + this->block_indices.reinit(partitioners.size(), 0); + + // initialize each block + this->components.resize(this->n_blocks()); + for (unsigned int i = 0; i < this->n_blocks(); ++i) + this->components[i].reinit(partitioners[i], make_ghosted); + + // update block_indices content + this->collect_sizes(); + } + + + inline MPI_Comm BlockVector::get_mpi_communicator() const { diff --git a/include/deal.II/lac/trilinos_parallel_block_vector.h b/include/deal.II/lac/trilinos_parallel_block_vector.h index 160792b7e6..33c3ac309f 100644 --- a/include/deal.II/lac/trilinos_parallel_block_vector.h +++ b/include/deal.II/lac/trilinos_parallel_block_vector.h @@ -215,6 +215,22 @@ namespace TrilinosWrappers const MPI_Comm communicator = MPI_COMM_WORLD, const bool vector_writable = false); + /** + * Initialize each block given to each parallel partitioning described in + * @p partitioners. + * + * You can decide whether your vector will contain ghost elements with + * @p make_ghosted. + * + * The parameter @p vector_writable only has effect on ghosted vectors + * and is ignored for non-ghosted vectors. + */ + void + reinit( + const std::vector> + & partitioners, + const bool make_ghosted = true, + const bool vector_writable = false); /** * Change the dimension to that of the vector V. The same diff --git a/source/lac/trilinos_block_vector.cc b/source/lac/trilinos_block_vector.cc index d72464c69c..4ce2b52c91 100644 --- a/source/lac/trilinos_block_vector.cc +++ b/source/lac/trilinos_block_vector.cc @@ -113,6 +113,29 @@ namespace TrilinosWrappers + void + BlockVector::reinit( + const std::vector> + & partitioners, + const bool make_ghosted, + const bool vector_writable) + { + // update the number of blocks + this->block_indices.reinit(partitioners.size(), 0); + + // initialize each block + this->components.resize(this->n_blocks()); + for (unsigned int i = 0; i < this->n_blocks(); ++i) + this->components[i].reinit(partitioners[i], + make_ghosted, + vector_writable); + + // update block_indices content + this->collect_sizes(); + } + + + void BlockVector::reinit(const BlockVector &v, const bool omit_zeroing_entries) { diff --git a/tests/mpi/parallel_block_vector_14.cc b/tests/mpi/parallel_block_vector_reinit_01.cc similarity index 96% rename from tests/mpi/parallel_block_vector_14.cc rename to tests/mpi/parallel_block_vector_reinit_01.cc index 071c0d44d9..1fd91de201 100644 --- a/tests/mpi/parallel_block_vector_14.cc +++ b/tests/mpi/parallel_block_vector_reinit_01.cc @@ -1,6 +1,6 @@ // --------------------------------------------------------------------- // -// Copyright (C) 2022 by the deal.II authors +// Copyright (C) 2022 - 2023 by the deal.II authors // // This file is part of the deal.II library. // @@ -79,7 +79,7 @@ test() { LinearAlgebra::distributed::BlockVector block_vector(partitioners); - deallog << "partitioner: OK" << std::endl; + deallog << "partitioners: OK" << std::endl; } } diff --git a/tests/mpi/parallel_block_vector_14.mpirun=10.output b/tests/mpi/parallel_block_vector_reinit_01.mpirun=10.output similarity index 71% rename from tests/mpi/parallel_block_vector_14.mpirun=10.output rename to tests/mpi/parallel_block_vector_reinit_01.mpirun=10.output index 9d10939766..2edc592cfe 100644 --- a/tests/mpi/parallel_block_vector_14.mpirun=10.output +++ b/tests/mpi/parallel_block_vector_reinit_01.mpirun=10.output @@ -1,4 +1,4 @@ DEAL::w/o ghost indices: OK DEAL::w/ ghost indices: OK -DEAL::partitioner: OK +DEAL::partitioners: OK diff --git a/tests/petsc/parallel_block_vector_reinit_01.cc b/tests/petsc/parallel_block_vector_reinit_01.cc new file mode 100644 index 0000000000..de1885f4fa --- /dev/null +++ b/tests/petsc/parallel_block_vector_reinit_01.cc @@ -0,0 +1,101 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2022 - 2023 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. +// +// --------------------------------------------------------------------- + + +// test various constructors and reinit functions + +#include +#include + +#include + +#include + +#include "../tests.h" + + +void +test() +{ + MPI_Comm comm = MPI_COMM_WORLD; + + const unsigned int myid = Utilities::MPI::this_mpi_process(comm); + const unsigned int n_procs = Utilities::MPI::n_mpi_processes(comm); + + constexpr unsigned int n_blocks = 3; + constexpr unsigned int n_indices_per_proc_and_block = 2; + + const unsigned int n_indices_per_block = + n_indices_per_proc_and_block * n_procs; + + // set up partitioning + std::vector owned_indexsets; + std::vector relevant_indexsets; + std::vector> partitioners; + for (unsigned int b = 0; b < n_blocks; ++b) + { + // For PETSc, the IndexSet for each block must be complete over all + // processes, i.e., their union has to start from zero and must cover all + // indices. + + IndexSet owned(n_indices_per_block); + owned.add_range(n_indices_per_proc_and_block * myid, + n_indices_per_proc_and_block * (myid + 1)); + + IndexSet relevant(n_indices_per_block); + relevant.add_range(1, 2); + + partitioners.push_back( + std::make_shared(owned, + relevant, + comm)); + owned_indexsets.push_back(std::move(owned)); + relevant_indexsets.push_back(std::move(relevant)); + } + + // create block vectors using different constructors + PETScWrappers::MPI::BlockVector block_vector; + + block_vector.reinit(owned_indexsets, comm); + AssertThrow(block_vector.has_ghost_elements() == false, ExcInternalError()); + deallog << "w/o ghost indices: OK" << std::endl; + + block_vector.reinit(owned_indexsets, relevant_indexsets, comm); + AssertThrow(block_vector.has_ghost_elements() == true, ExcInternalError()); + deallog << "w/ ghost indices: OK" << std::endl; + + block_vector.reinit(partitioners, /*make_ghosted=*/false); + AssertThrow(block_vector.has_ghost_elements() == false, ExcInternalError()); + deallog << "partitioners w/o ghost indices: OK" << std::endl; + + block_vector.reinit(partitioners, /*make_ghosted=*/true); + AssertThrow(block_vector.has_ghost_elements() == true, ExcInternalError()); + deallog << "partitioners w/ ghost indices: OK" << std::endl; +} + + + +int +main(int argc, char **argv) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + + const unsigned int myid = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD); + + if (myid == 0) + initlog(); + + test(); +} diff --git a/tests/petsc/parallel_block_vector_reinit_01.mpirun=10.output b/tests/petsc/parallel_block_vector_reinit_01.mpirun=10.output new file mode 100644 index 0000000000..0daae16348 --- /dev/null +++ b/tests/petsc/parallel_block_vector_reinit_01.mpirun=10.output @@ -0,0 +1,5 @@ + +DEAL::w/o ghost indices: OK +DEAL::w/ ghost indices: OK +DEAL::partitioners w/o ghost indices: OK +DEAL::partitioners w/ ghost indices: OK diff --git a/tests/trilinos/parallel_block_vector_reinit_01.cc b/tests/trilinos/parallel_block_vector_reinit_01.cc new file mode 100644 index 0000000000..9ccece07f1 --- /dev/null +++ b/tests/trilinos/parallel_block_vector_reinit_01.cc @@ -0,0 +1,100 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2022 - 2023 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. +// +// --------------------------------------------------------------------- + + +// test various constructors and reinit functions + +#include +#include + +#include + +#include + +#include "../tests.h" + + +void +test() +{ + MPI_Comm comm = MPI_COMM_WORLD; + + const unsigned int myid = Utilities::MPI::this_mpi_process(comm); + const unsigned int n_procs = Utilities::MPI::n_mpi_processes(comm); + + constexpr unsigned int n_blocks = 3; + constexpr unsigned int n_indices_per_proc_and_block = 2; + + const unsigned int n_indices_per_block = + n_indices_per_proc_and_block * n_procs; + const unsigned int n_indices = n_blocks * n_indices_per_block; + + // set up partitioning + std::vector owned_indexsets; + std::vector relevant_indexsets; + std::vector> partitioners; + for (unsigned int b = 0; b < n_blocks; ++b) + { + const unsigned int begin = b * n_indices_per_block; + + IndexSet owned(n_indices); + owned.add_range(begin + n_indices_per_proc_and_block * myid, + begin + n_indices_per_proc_and_block * (myid + 1)); + + IndexSet relevant(n_indices); + relevant.add_range(begin + 1, begin + 2); + + partitioners.push_back( + std::make_shared(owned, + relevant, + comm)); + owned_indexsets.push_back(std::move(owned)); + relevant_indexsets.push_back(std::move(relevant)); + } + + // create block vectors using different constructors + TrilinosWrappers::MPI::BlockVector block_vector; + + block_vector.reinit(owned_indexsets, comm); + AssertThrow(block_vector.has_ghost_elements() == false, ExcInternalError()); + deallog << "w/o ghost indices: OK" << std::endl; + + block_vector.reinit(owned_indexsets, relevant_indexsets, comm); + AssertThrow(block_vector.has_ghost_elements() == true, ExcInternalError()); + deallog << "w/ ghost indices: OK" << std::endl; + + block_vector.reinit(partitioners, /*make_ghosted=*/false); + AssertThrow(block_vector.has_ghost_elements() == false, ExcInternalError()); + deallog << "partitioners w/o ghost indices: OK" << std::endl; + + block_vector.reinit(partitioners, /*make_ghosted=*/true); + AssertThrow(block_vector.has_ghost_elements() == true, ExcInternalError()); + deallog << "partitioners w/ ghost indices: OK" << std::endl; +} + + + +int +main(int argc, char **argv) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + + const unsigned int myid = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD); + + if (myid == 0) + initlog(); + + test(); +} diff --git a/tests/trilinos/parallel_block_vector_reinit_01.mpirun=10.output b/tests/trilinos/parallel_block_vector_reinit_01.mpirun=10.output new file mode 100644 index 0000000000..0daae16348 --- /dev/null +++ b/tests/trilinos/parallel_block_vector_reinit_01.mpirun=10.output @@ -0,0 +1,5 @@ + +DEAL::w/o ghost indices: OK +DEAL::w/ ghost indices: OK +DEAL::partitioners w/o ghost indices: OK +DEAL::partitioners w/ ghost indices: OK -- 2.39.5