From: David Wells Date: Sun, 24 May 2020 00:03:51 +0000 (-0400) Subject: local_size() -> locally_owned_size() X-Git-Tag: v9.3.0-rc1~449^2~13 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=6e76252b941b8ce65db9d6ce6e5578b238766e8c;p=dealii.git local_size() -> locally_owned_size() --- diff --git a/include/deal.II/lac/block_vector_base.h b/include/deal.II/lac/block_vector_base.h index fc5dca1fae..2833ca3d07 100644 --- a/include/deal.II/lac/block_vector_base.h +++ b/include/deal.II/lac/block_vector_base.h @@ -550,7 +550,7 @@ public: * components. */ std::size_t - local_size() const; + locally_owned_size() const; /** * Return an index set that describes which elements of this vector are @@ -1440,11 +1440,11 @@ BlockVectorBase::size() const template inline std::size_t -BlockVectorBase::local_size() const +BlockVectorBase::locally_owned_size() const { std::size_t local_size = 0; for (unsigned int b = 0; b < n_blocks(); ++b) - local_size += block(b).local_size(); + local_size += block(b).locally_owned_size(); return local_size; } 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 0274c0aa93..75bb7d2707 100644 --- a/include/deal.II/lac/la_parallel_block_vector.templates.h +++ b/include/deal.II/lac/la_parallel_block_vector.templates.h @@ -269,7 +269,7 @@ namespace LinearAlgebra &start_ptr); AssertThrow(ierr == 0, ExcPETScError(ierr)); - const size_type vec_size = this->block(i).local_size(); + const size_type vec_size = this->block(i).locally_owned_size(); petsc_helpers::copy_petsc_vector(start_ptr, start_ptr + vec_size, this->block(i).begin()); diff --git a/include/deal.II/lac/la_parallel_vector.h b/include/deal.II/lac/la_parallel_vector.h index 05e2ad31fa..42b1220d04 100644 --- a/include/deal.II/lac/la_parallel_vector.h +++ b/include/deal.II/lac/la_parallel_vector.h @@ -105,8 +105,9 @@ namespace LinearAlgebra *
  • Besides the usual global access operator() it is also possible to * access vector entries in the local index space with the function @p * local_element(). Locally owned indices are placed first, [0, - * local_size()), and then all ghost indices follow after them - * contiguously, [local_size(), local_size()+n_ghost_entries()). + * locally_owned_size()), and then all ghost indices follow after them + * contiguously, [locally_owned_size(), + * locally_owned_size()+n_ghost_entries()). * * * Functions related to parallel functionality: @@ -923,10 +924,20 @@ namespace LinearAlgebra /** * Return the local size of the vector, i.e., the number of indices * owned locally. + * + * @deprecated use locally_owned_size() instead. */ + DEAL_II_DEPRECATED size_type local_size() const; + /** + * Return the local size of the vector, i.e., the number of indices + * owned locally. + */ + size_type + locally_owned_size() const; + /** * Return true if the given global index is in the local range of this * processor. @@ -1519,6 +1530,15 @@ namespace LinearAlgebra + template + inline typename Vector::size_type + Vector::locally_owned_size() const + { + return partitioner->local_size(); + } + + + template inline bool Vector::in_local_range( diff --git a/include/deal.II/lac/petsc_vector_base.h b/include/deal.II/lac/petsc_vector_base.h index d3b4cb91e6..5ec8af7888 100644 --- a/include/deal.II/lac/petsc_vector_base.h +++ b/include/deal.II/lac/petsc_vector_base.h @@ -357,10 +357,24 @@ namespace PETScWrappers * * To figure out which elements exactly are stored locally, use * local_range(). + * + * @deprecated use locally_owned_size() instead. */ + DEAL_II_DEPRECATED size_type local_size() const; + /** + * Return the local dimension of the vector, i.e. the number of elements + * stored on the present MPI process. For sequential vectors, this number + * is the same as size(), but for parallel vectors it may be smaller. + * + * To figure out which elements exactly are stored locally, use + * local_range(). + */ + size_type + locally_owned_size() const; + /** * Return a pair of indices indicating which elements of this vector are * stored locally. The first number is the index of the first element diff --git a/include/deal.II/lac/read_write_vector.h b/include/deal.II/lac/read_write_vector.h index 1df77dd799..1ce0fac946 100644 --- a/include/deal.II/lac/read_write_vector.h +++ b/include/deal.II/lac/read_write_vector.h @@ -400,7 +400,7 @@ namespace LinearAlgebra * 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. + * @deprecated use locally_owned_size() instead. */ DEAL_II_DEPRECATED size_type @@ -411,8 +411,7 @@ namespace LinearAlgebra * owned locally. */ size_type - local_size() const; - + locally_owned_size() const; /** * Return the IndexSet that represents the indices of the elements stored. @@ -834,7 +833,7 @@ namespace LinearAlgebra template inline typename ReadWriteVector::size_type - ReadWriteVector::local_size() const + ReadWriteVector::locally_owned_size() const { return stored_elements.n_elements(); } diff --git a/include/deal.II/lac/read_write_vector.templates.h b/include/deal.II/lac/read_write_vector.templates.h index fc563239d5..639becf5f4 100644 --- a/include/deal.II/lac/read_write_vector.templates.h +++ b/include/deal.II/lac/read_write_vector.templates.h @@ -494,7 +494,7 @@ namespace LinearAlgebra VecGetArray(static_cast(petsc_vec), &start_ptr); AssertThrow(ierr == 0, ExcPETScError(ierr)); - const size_type vec_size = petsc_vec.local_size(); + const size_type vec_size = petsc_vec.locally_owned_size(); internal::copy_petsc_vector(start_ptr, start_ptr + vec_size, begin()); // restore the representation of the vector diff --git a/include/deal.II/lac/trilinos_epetra_vector.h b/include/deal.II/lac/trilinos_epetra_vector.h index 421c7865c2..b840c3da50 100644 --- a/include/deal.II/lac/trilinos_epetra_vector.h +++ b/include/deal.II/lac/trilinos_epetra_vector.h @@ -289,7 +289,7 @@ namespace LinearAlgebra * owned locally. */ size_type - local_size() const; + locally_owned_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 729a6d4806..91be6bf4f4 100644 --- a/include/deal.II/lac/trilinos_tpetra_vector.h +++ b/include/deal.II/lac/trilinos_tpetra_vector.h @@ -312,7 +312,7 @@ namespace LinearAlgebra * owned locally. */ size_type - local_size() const; + locally_owned_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 55b04638a2..a6bcd44a1f 100644 --- a/include/deal.II/lac/trilinos_tpetra_vector.templates.h +++ b/include/deal.II/lac/trilinos_tpetra_vector.templates.h @@ -550,7 +550,7 @@ namespace LinearAlgebra template typename Vector::size_type - Vector::local_size() const + Vector::locally_owned_size() const { return vector->getLocalLength(); } diff --git a/include/deal.II/lac/trilinos_vector.h b/include/deal.II/lac/trilinos_vector.h index 7daa481a00..d46b23f4d0 100644 --- a/include/deal.II/lac/trilinos_vector.h +++ b/include/deal.II/lac/trilinos_vector.h @@ -738,10 +738,20 @@ namespace TrilinosWrappers * * If the vector contains ghost elements, they are included in this * number. + * + * @deprecated This function is deprecated. */ + DEAL_II_DEPRECATED size_type local_size() const; + /** + * Return the local size of the vector, i.e., the number of indices + * owned locally. + */ + size_type + locally_owned_size() const; + /** * Return a pair of indices indicating which elements of this vector are * stored locally. The first number is the index of the first element @@ -1739,6 +1749,14 @@ namespace TrilinosWrappers + inline Vector::size_type + Vector::locally_owned_size() const + { + return owned_elements.n_elements(); + } + + + inline std::pair Vector::local_range() const { diff --git a/include/deal.II/lac/vector.h b/include/deal.II/lac/vector.h index de32e87ddc..d7c8eed32d 100644 --- a/include/deal.II/lac/vector.h +++ b/include/deal.II/lac/vector.h @@ -961,7 +961,7 @@ public: * LinearAlgebra::ReadWriteVector. */ size_type - local_size() const; + locally_owned_size() const; /** * Return whether the vector contains only elements with value zero. This @@ -1101,7 +1101,7 @@ Vector::size() const template inline typename Vector::size_type -Vector::local_size() const +Vector::locally_owned_size() const { return values.size(); } diff --git a/source/lac/petsc_vector_base.cc b/source/lac/petsc_vector_base.cc index 75608c2953..3b9cf737fa 100644 --- a/source/lac/petsc_vector_base.cc +++ b/source/lac/petsc_vector_base.cc @@ -257,6 +257,18 @@ namespace PETScWrappers + VectorBase::size_type + VectorBase::locally_owned_size() const + { + PetscInt sz; + const PetscErrorCode ierr = VecGetLocalSize(vector, &sz); + AssertThrow(ierr == 0, ExcPETScError(ierr)); + + return sz; + } + + + VectorBase::size_type VectorBase::local_size() const { @@ -869,10 +881,10 @@ namespace PETScWrappers out.setf(std::ios::fixed, std::ios::floatfield); if (across) - for (size_type i = 0; i < local_size(); ++i) + for (size_type i = 0; i < locally_owned_size(); ++i) out << val[i] << ' '; else - for (size_type i = 0; i < local_size(); ++i) + for (size_type i = 0; i < locally_owned_size(); ++i) out << val[i] << std::endl; out << std::endl; @@ -915,7 +927,7 @@ namespace PETScWrappers // TH: I am relatively sure that PETSc is // storing the local data in a contiguous // block without indices: - mem += local_size() * sizeof(PetscScalar); + mem += locally_owned_size() * sizeof(PetscScalar); // assume that PETSc is storing one index // and one double per ghost element if (ghosted) diff --git a/source/lac/trilinos_epetra_vector.cc b/source/lac/trilinos_epetra_vector.cc index e45f395f88..dd50626954 100644 --- a/source/lac/trilinos_epetra_vector.cc +++ b/source/lac/trilinos_epetra_vector.cc @@ -541,7 +541,7 @@ namespace LinearAlgebra Vector::size_type - Vector::local_size() const + Vector::locally_owned_size() const { return vector->MyLength(); } diff --git a/tests/mpi/local_size.cc b/tests/mpi/local_size.cc index 1fcb51bfa4..e2c4ca0eb8 100644 --- a/tests/mpi/local_size.cc +++ b/tests/mpi/local_size.cc @@ -13,7 +13,7 @@ // // --------------------------------------------------------------------- -// check Vector::local_size() for all supported vector types +// check Vector::locally_owned_size() for all supported vector types #include #include @@ -37,9 +37,9 @@ void check_serial() { const auto dofs_per_proc = 4; - VEC vec(dofs_per_proc); + VEC vec(dofs_per_proc); deallog << "type: " << Utilities::type_to_string(vec) << std::endl; - deallog << "local size: " << vec.local_size() << std::endl; + deallog << "local size: " << vec.locally_owned_size() << std::endl; deallog << "size: " << vec.size() << std::endl; } @@ -53,18 +53,18 @@ check_unghosted_parallel() 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; + 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); + 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 << "local size: " << vec.locally_owned_size() << std::endl; deallog << "size: " << vec.size() << std::endl; } @@ -78,12 +78,12 @@ check_ghosted_parallel() 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; + 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); + 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) @@ -99,7 +99,7 @@ check_ghosted_parallel() 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 << "local size: " << vec.locally_owned_size() << std::endl; deallog << "size: " << vec.size() << std::endl; } @@ -113,12 +113,12 @@ check_ghosted_parallel_block() 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; + 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); + 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) @@ -131,13 +131,13 @@ check_ghosted_parallel_block() ghost_indices.add_index(my_dofs_end % n_dofs); } - std::vector local_blocks {local_indices, local_indices}; + std::vector local_blocks{local_indices, local_indices}; // for variety do not ghost the second component - std::vector ghost_blocks {ghost_indices, IndexSet()}; + 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 << "local size: " << vec.locally_owned_size() << std::endl; deallog << "size: " << vec.size() << std::endl; } @@ -161,7 +161,8 @@ main(int argc, char *argv[]) check_ghosted_parallel(); // block vectors: - check_ghosted_parallel_block>(); + check_ghosted_parallel_block< + LinearAlgebra::distributed::BlockVector>(); 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 index 622b61fcb9..91e505e0aa 100644 --- 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 @@ -23,7 +23,7 @@ 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::local size: 4 DEAL:0::size: 16 DEAL:0::type: dealii::LinearAlgebra::distributed::BlockVector DEAL:0::local size: 8 @@ -32,7 +32,7 @@ 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::local size: 8 DEAL:0::size: 32 DEAL:1::type: dealii::Vector @@ -59,7 +59,7 @@ 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::local size: 4 DEAL:1::size: 16 DEAL:1::type: dealii::LinearAlgebra::distributed::BlockVector DEAL:1::local size: 8 @@ -68,7 +68,7 @@ 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::local size: 8 DEAL:1::size: 32 @@ -96,7 +96,7 @@ 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::local size: 4 DEAL:2::size: 16 DEAL:2::type: dealii::LinearAlgebra::distributed::BlockVector DEAL:2::local size: 8 @@ -105,7 +105,7 @@ 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::local size: 8 DEAL:2::size: 32 @@ -133,7 +133,7 @@ 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::local size: 4 DEAL:3::size: 16 DEAL:3::type: dealii::LinearAlgebra::distributed::BlockVector DEAL:3::local size: 8 @@ -142,6 +142,6 @@ 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::local size: 8 DEAL:3::size: 32