From 61b5c0ca418024d16ecf894cd6ec23cf14c0eeeb Mon Sep 17 00:00:00 2001 From: David Wells Date: Fri, 22 May 2020 21:52:17 -0400 Subject: [PATCH] consistently use local_size(). This is the name given to the function by la_parallel_vector. Defining this everywhere gives us a true check as to whether or not a vector is distributed - i.e., if size() != local_size() then the vector is distributed and otherwise it is not. --- include/deal.II/lac/block_vector_base.h | 20 +++ include/deal.II/lac/la_parallel_vector.h | 2 +- include/deal.II/lac/read_write_vector.h | 20 +++ include/deal.II/lac/trilinos_epetra_vector.h | 7 + include/deal.II/lac/trilinos_tpetra_vector.h | 7 + .../lac/trilinos_tpetra_vector.templates.h | 9 + include/deal.II/lac/vector.h | 20 +++ source/lac/trilinos_epetra_vector.cc | 8 + tests/mpi/local_size.cc | 167 ++++++++++++++++++ ...th_tpetra=on.with_petsc=on.mpirun=4.output | 147 +++++++++++++++ 10 files changed, 406 insertions(+), 1 deletion(-) create mode 100644 tests/mpi/local_size.cc create mode 100644 tests/mpi/local_size.with_trilinos_with_tpetra=on.with_petsc=on.mpirun=4.output diff --git a/include/deal.II/lac/block_vector_base.h b/include/deal.II/lac/block_vector_base.h index 870a6ae635..fc5dca1fae 100644 --- a/include/deal.II/lac/block_vector_base.h +++ b/include/deal.II/lac/block_vector_base.h @@ -544,6 +544,14 @@ public: std::size_t size() const; + /** + * Return local dimension of the vector. This is the sum of the local + * dimensions (i.e., values stored on the current processor) of all + * components. + */ + std::size_t + local_size() const; + /** * Return an index set that describes which elements of this vector are * owned by the current processor. Note that this index set does not include @@ -1430,6 +1438,18 @@ BlockVectorBase::size() const +template +inline std::size_t +BlockVectorBase::local_size() const +{ + std::size_t local_size = 0; + for (unsigned int b = 0; b < n_blocks(); ++b) + local_size += block(b).local_size(); + return local_size; +} + + + template inline IndexSet BlockVectorBase::locally_owned_elements() const diff --git a/include/deal.II/lac/la_parallel_vector.h b/include/deal.II/lac/la_parallel_vector.h index 209de68f5c..05e2ad31fa 100644 --- a/include/deal.II/lac/la_parallel_vector.h +++ b/include/deal.II/lac/la_parallel_vector.h @@ -1514,7 +1514,7 @@ namespace LinearAlgebra inline typename Vector::size_type Vector::local_size() const { - return partitioner->local_size(); + return locally_owned_size(); } diff --git a/include/deal.II/lac/read_write_vector.h b/include/deal.II/lac/read_write_vector.h index bb6f47e98d..1df77dd799 100644 --- a/include/deal.II/lac/read_write_vector.h +++ b/include/deal.II/lac/read_write_vector.h @@ -399,10 +399,21 @@ namespace LinearAlgebra * This function returns the number of elements stored. It is smaller or * equal to the dimension of the vector space that is modeled by an object * of this kind. This dimension is return by size(). + * + * @deprecated use local_size() instead. */ + DEAL_II_DEPRECATED size_type n_elements() const; + /** + * Return the local size of the vector, i.e., the number of indices + * owned locally. + */ + size_type + local_size() const; + + /** * Return the IndexSet that represents the indices of the elements stored. */ @@ -821,6 +832,15 @@ namespace LinearAlgebra + template + inline typename ReadWriteVector::size_type + ReadWriteVector::local_size() const + { + return stored_elements.n_elements(); + } + + + template inline const IndexSet & ReadWriteVector::get_stored_elements() const diff --git a/include/deal.II/lac/trilinos_epetra_vector.h b/include/deal.II/lac/trilinos_epetra_vector.h index 34ab495165..421c7865c2 100644 --- a/include/deal.II/lac/trilinos_epetra_vector.h +++ b/include/deal.II/lac/trilinos_epetra_vector.h @@ -284,6 +284,13 @@ namespace LinearAlgebra virtual size_type size() const override; + /** + * Return the local size of the vector, i.e., the number of indices + * owned locally. + */ + size_type + local_size() const; + /** * Return the MPI communicator object in use with this object. */ diff --git a/include/deal.II/lac/trilinos_tpetra_vector.h b/include/deal.II/lac/trilinos_tpetra_vector.h index 95fb600629..729a6d4806 100644 --- a/include/deal.II/lac/trilinos_tpetra_vector.h +++ b/include/deal.II/lac/trilinos_tpetra_vector.h @@ -307,6 +307,13 @@ namespace LinearAlgebra virtual size_type size() const override; + /** + * Return the local size of the vector, i.e., the number of indices + * owned locally. + */ + size_type + local_size() const; + /** * Return the MPI communicator object in use with this object. */ diff --git a/include/deal.II/lac/trilinos_tpetra_vector.templates.h b/include/deal.II/lac/trilinos_tpetra_vector.templates.h index d75d263d5a..55b04638a2 100644 --- a/include/deal.II/lac/trilinos_tpetra_vector.templates.h +++ b/include/deal.II/lac/trilinos_tpetra_vector.templates.h @@ -548,6 +548,15 @@ namespace LinearAlgebra + template + typename Vector::size_type + Vector::local_size() const + { + return vector->getLocalLength(); + } + + + template MPI_Comm Vector::get_mpi_communicator() const diff --git a/include/deal.II/lac/vector.h b/include/deal.II/lac/vector.h index 221ff506a0..de32e87ddc 100644 --- a/include/deal.II/lac/vector.h +++ b/include/deal.II/lac/vector.h @@ -953,6 +953,16 @@ public: size_type size() const; + /** + * Return local dimension of the vector. Since this vector does not support + * distributed data this is always the same value as size(). + * + * @note This function exists for compatibility with + * LinearAlgebra::ReadWriteVector. + */ + size_type + local_size() const; + /** * Return whether the vector contains only elements with value zero. This * function is mainly for internal consistency checks and should seldom be @@ -1088,6 +1098,16 @@ Vector::size() const } + +template +inline typename Vector::size_type +Vector::local_size() const +{ + return values.size(); +} + + + template inline bool Vector::in_local_range(const size_type) const diff --git a/source/lac/trilinos_epetra_vector.cc b/source/lac/trilinos_epetra_vector.cc index 54b4587642..e45f395f88 100644 --- a/source/lac/trilinos_epetra_vector.cc +++ b/source/lac/trilinos_epetra_vector.cc @@ -540,6 +540,14 @@ namespace LinearAlgebra + Vector::size_type + Vector::local_size() const + { + return vector->MyLength(); + } + + + MPI_Comm Vector::get_mpi_communicator() const { diff --git a/tests/mpi/local_size.cc b/tests/mpi/local_size.cc new file mode 100644 index 0000000000..1fcb51bfa4 --- /dev/null +++ b/tests/mpi/local_size.cc @@ -0,0 +1,167 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2020 - 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. +// +// --------------------------------------------------------------------- + +// check Vector::local_size() for all supported vector types + +#include +#include + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include "../tests.h" + + +template +void +check_serial() +{ + const auto dofs_per_proc = 4; + VEC vec(dofs_per_proc); + deallog << "type: " << Utilities::type_to_string(vec) << std::endl; + deallog << "local size: " << vec.local_size() << std::endl; + deallog << "size: " << vec.size() << std::endl; +} + + + +template +void +check_unghosted_parallel() +{ + const auto n_procs = dealii::Utilities::MPI::n_mpi_processes(MPI_COMM_WORLD); + const auto my_rank = dealii::Utilities::MPI::this_mpi_process(MPI_COMM_WORLD); + + const auto dofs_per_proc = 4; + const auto n_dofs = dofs_per_proc*n_procs; + + IndexSet local_indices(n_dofs); + const auto my_dofs_begin = dofs_per_proc*my_rank; + const auto my_dofs_end = dofs_per_proc*(my_rank + 1); + local_indices.add_range(my_dofs_begin, my_dofs_end); + local_indices.compress(); + + VEC vec(local_indices, MPI_COMM_WORLD); + deallog << "type: " << Utilities::type_to_string(vec) << std::endl; + deallog << "index set size: " << local_indices.n_elements() << std::endl; + deallog << "local size: " << vec.local_size() << std::endl; + deallog << "size: " << vec.size() << std::endl; +} + + + +template +void +check_ghosted_parallel() +{ + const auto n_procs = dealii::Utilities::MPI::n_mpi_processes(MPI_COMM_WORLD); + const auto my_rank = dealii::Utilities::MPI::this_mpi_process(MPI_COMM_WORLD); + + const auto dofs_per_proc = 4; + const auto n_dofs = dofs_per_proc*n_procs; + + IndexSet local_indices(n_dofs); + IndexSet ghost_indices(n_dofs); + const auto my_dofs_begin = dofs_per_proc*my_rank; + const auto my_dofs_end = dofs_per_proc*(my_rank + 1); + local_indices.add_range(my_dofs_begin, my_dofs_end); + local_indices.compress(); + if (my_rank == 0) + { + ghost_indices.add_index(my_dofs_end); + } + else + { + ghost_indices.add_index(my_dofs_begin - 1); + ghost_indices.add_index(my_dofs_end % n_dofs); + } + + VEC vec(local_indices, ghost_indices, MPI_COMM_WORLD); + deallog << "type: " << Utilities::type_to_string(vec) << std::endl; + deallog << "index set size: " << local_indices.n_elements() << std::endl; + deallog << "local size: " << vec.local_size() << std::endl; + deallog << "size: " << vec.size() << std::endl; +} + + + +template +void +check_ghosted_parallel_block() +{ + const auto n_procs = dealii::Utilities::MPI::n_mpi_processes(MPI_COMM_WORLD); + const auto my_rank = dealii::Utilities::MPI::this_mpi_process(MPI_COMM_WORLD); + + const auto dofs_per_proc = 4; + const auto n_dofs = dofs_per_proc*n_procs; + + IndexSet local_indices(n_dofs); + IndexSet ghost_indices(n_dofs); + const auto my_dofs_begin = dofs_per_proc*my_rank; + const auto my_dofs_end = dofs_per_proc*(my_rank + 1); + local_indices.add_range(my_dofs_begin, my_dofs_end); + local_indices.compress(); + if (my_rank == 0) + { + ghost_indices.add_index(my_dofs_end); + } + else + { + ghost_indices.add_index(my_dofs_begin - 1); + ghost_indices.add_index(my_dofs_end % n_dofs); + } + + std::vector local_blocks {local_indices, local_indices}; + // for variety do not ghost the second component + std::vector ghost_blocks {ghost_indices, IndexSet()}; + + VEC vec(local_blocks, ghost_blocks, MPI_COMM_WORLD); + deallog << "type: " << Utilities::type_to_string(vec) << std::endl; + deallog << "local size: " << vec.local_size() << std::endl; + deallog << "size: " << vec.size() << std::endl; +} + + + +int +main(int argc, char *argv[]) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + MPILogInitAll log; + + // non-block vectors: + check_serial>(); + check_serial>(); + + check_unghosted_parallel(); + check_unghosted_parallel>(); + + check_ghosted_parallel>(); + check_ghosted_parallel(); + check_ghosted_parallel(); + + // block vectors: + check_ghosted_parallel_block>(); + check_ghosted_parallel_block(); + check_ghosted_parallel_block(); +} diff --git a/tests/mpi/local_size.with_trilinos_with_tpetra=on.with_petsc=on.mpirun=4.output b/tests/mpi/local_size.with_trilinos_with_tpetra=on.with_petsc=on.mpirun=4.output new file mode 100644 index 0000000000..622b61fcb9 --- /dev/null +++ b/tests/mpi/local_size.with_trilinos_with_tpetra=on.with_petsc=on.mpirun=4.output @@ -0,0 +1,147 @@ + +DEAL:0::type: dealii::Vector +DEAL:0::local size: 4 +DEAL:0::size: 4 +DEAL:0::type: dealii::LinearAlgebra::Vector +DEAL:0::local size: 4 +DEAL:0::size: 4 +DEAL:0::type: dealii::LinearAlgebra::EpetraWrappers::Vector +DEAL:0::index set size: 4 +DEAL:0::local size: 4 +DEAL:0::size: 16 +DEAL:0::type: dealii::LinearAlgebra::TpetraWrappers::Vector +DEAL:0::index set size: 4 +DEAL:0::local size: 4 +DEAL:0::size: 16 +DEAL:0::type: dealii::LinearAlgebra::distributed::Vector +DEAL:0::index set size: 4 +DEAL:0::local size: 4 +DEAL:0::size: 16 +DEAL:0::type: dealii::PETScWrappers::MPI::Vector +DEAL:0::index set size: 4 +DEAL:0::local size: 4 +DEAL:0::size: 16 +DEAL:0::type: dealii::TrilinosWrappers::MPI::Vector +DEAL:0::index set size: 4 +DEAL:0::local size: 5 +DEAL:0::size: 16 +DEAL:0::type: dealii::LinearAlgebra::distributed::BlockVector +DEAL:0::local size: 8 +DEAL:0::size: 32 +DEAL:0::type: dealii::PETScWrappers::MPI::BlockVector +DEAL:0::local size: 8 +DEAL:0::size: 32 +DEAL:0::type: dealii::TrilinosWrappers::MPI::BlockVector +DEAL:0::local size: 9 +DEAL:0::size: 32 + +DEAL:1::type: dealii::Vector +DEAL:1::local size: 4 +DEAL:1::size: 4 +DEAL:1::type: dealii::LinearAlgebra::Vector +DEAL:1::local size: 4 +DEAL:1::size: 4 +DEAL:1::type: dealii::LinearAlgebra::EpetraWrappers::Vector +DEAL:1::index set size: 4 +DEAL:1::local size: 4 +DEAL:1::size: 16 +DEAL:1::type: dealii::LinearAlgebra::TpetraWrappers::Vector +DEAL:1::index set size: 4 +DEAL:1::local size: 4 +DEAL:1::size: 16 +DEAL:1::type: dealii::LinearAlgebra::distributed::Vector +DEAL:1::index set size: 4 +DEAL:1::local size: 4 +DEAL:1::size: 16 +DEAL:1::type: dealii::PETScWrappers::MPI::Vector +DEAL:1::index set size: 4 +DEAL:1::local size: 4 +DEAL:1::size: 16 +DEAL:1::type: dealii::TrilinosWrappers::MPI::Vector +DEAL:1::index set size: 4 +DEAL:1::local size: 6 +DEAL:1::size: 16 +DEAL:1::type: dealii::LinearAlgebra::distributed::BlockVector +DEAL:1::local size: 8 +DEAL:1::size: 32 +DEAL:1::type: dealii::PETScWrappers::MPI::BlockVector +DEAL:1::local size: 8 +DEAL:1::size: 32 +DEAL:1::type: dealii::TrilinosWrappers::MPI::BlockVector +DEAL:1::local size: 10 +DEAL:1::size: 32 + + +DEAL:2::type: dealii::Vector +DEAL:2::local size: 4 +DEAL:2::size: 4 +DEAL:2::type: dealii::LinearAlgebra::Vector +DEAL:2::local size: 4 +DEAL:2::size: 4 +DEAL:2::type: dealii::LinearAlgebra::EpetraWrappers::Vector +DEAL:2::index set size: 4 +DEAL:2::local size: 4 +DEAL:2::size: 16 +DEAL:2::type: dealii::LinearAlgebra::TpetraWrappers::Vector +DEAL:2::index set size: 4 +DEAL:2::local size: 4 +DEAL:2::size: 16 +DEAL:2::type: dealii::LinearAlgebra::distributed::Vector +DEAL:2::index set size: 4 +DEAL:2::local size: 4 +DEAL:2::size: 16 +DEAL:2::type: dealii::PETScWrappers::MPI::Vector +DEAL:2::index set size: 4 +DEAL:2::local size: 4 +DEAL:2::size: 16 +DEAL:2::type: dealii::TrilinosWrappers::MPI::Vector +DEAL:2::index set size: 4 +DEAL:2::local size: 6 +DEAL:2::size: 16 +DEAL:2::type: dealii::LinearAlgebra::distributed::BlockVector +DEAL:2::local size: 8 +DEAL:2::size: 32 +DEAL:2::type: dealii::PETScWrappers::MPI::BlockVector +DEAL:2::local size: 8 +DEAL:2::size: 32 +DEAL:2::type: dealii::TrilinosWrappers::MPI::BlockVector +DEAL:2::local size: 10 +DEAL:2::size: 32 + + +DEAL:3::type: dealii::Vector +DEAL:3::local size: 4 +DEAL:3::size: 4 +DEAL:3::type: dealii::LinearAlgebra::Vector +DEAL:3::local size: 4 +DEAL:3::size: 4 +DEAL:3::type: dealii::LinearAlgebra::EpetraWrappers::Vector +DEAL:3::index set size: 4 +DEAL:3::local size: 4 +DEAL:3::size: 16 +DEAL:3::type: dealii::LinearAlgebra::TpetraWrappers::Vector +DEAL:3::index set size: 4 +DEAL:3::local size: 4 +DEAL:3::size: 16 +DEAL:3::type: dealii::LinearAlgebra::distributed::Vector +DEAL:3::index set size: 4 +DEAL:3::local size: 4 +DEAL:3::size: 16 +DEAL:3::type: dealii::PETScWrappers::MPI::Vector +DEAL:3::index set size: 4 +DEAL:3::local size: 4 +DEAL:3::size: 16 +DEAL:3::type: dealii::TrilinosWrappers::MPI::Vector +DEAL:3::index set size: 4 +DEAL:3::local size: 6 +DEAL:3::size: 16 +DEAL:3::type: dealii::LinearAlgebra::distributed::BlockVector +DEAL:3::local size: 8 +DEAL:3::size: 32 +DEAL:3::type: dealii::PETScWrappers::MPI::BlockVector +DEAL:3::local size: 8 +DEAL:3::size: 32 +DEAL:3::type: dealii::TrilinosWrappers::MPI::BlockVector +DEAL:3::local size: 10 +DEAL:3::size: 32 + -- 2.39.5