]> https://gitweb.dealii.org/ - dealii.git/commitdiff
local_size() -> locally_owned_size()
authorDavid Wells <drwells@email.unc.edu>
Sun, 24 May 2020 00:03:51 +0000 (20:03 -0400)
committerDaniel Arndt <arndtd@ornl.gov>
Thu, 11 Feb 2021 16:48:49 +0000 (11:48 -0500)
15 files changed:
include/deal.II/lac/block_vector_base.h
include/deal.II/lac/la_parallel_block_vector.templates.h
include/deal.II/lac/la_parallel_vector.h
include/deal.II/lac/petsc_vector_base.h
include/deal.II/lac/read_write_vector.h
include/deal.II/lac/read_write_vector.templates.h
include/deal.II/lac/trilinos_epetra_vector.h
include/deal.II/lac/trilinos_tpetra_vector.h
include/deal.II/lac/trilinos_tpetra_vector.templates.h
include/deal.II/lac/trilinos_vector.h
include/deal.II/lac/vector.h
source/lac/petsc_vector_base.cc
source/lac/trilinos_epetra_vector.cc
tests/mpi/local_size.cc
tests/mpi/local_size.with_trilinos_with_tpetra=on.with_petsc=on.mpirun=4.output

index fc5dca1fae70e495300ce5b7ac70b41e7587b7b5..2833ca3d075060891fad552f8dc5bd37583853eb 100644 (file)
@@ -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<VectorType>::size() const
 
 template <class VectorType>
 inline std::size_t
-BlockVectorBase<VectorType>::local_size() const
+BlockVectorBase<VectorType>::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;
 }
 
index 0274c0aa93a6290de4c688df611e6b1cfd7489d4..75bb7d2707bf857ad5606d56730416c971de4f5d 100644 (file)
@@ -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());
index 05e2ad31fa701bdbdfea2f7ffb2a12386fe575df..42b1220d0432428d26e9ef2713fd7c13d1e42705 100644 (file)
@@ -105,8 +105,9 @@ namespace LinearAlgebra
      * <li> 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()).
      * </ul>
      *
      * 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 <typename Number, typename MemorySpace>
+    inline typename Vector<Number, MemorySpace>::size_type
+    Vector<Number, MemorySpace>::locally_owned_size() const
+    {
+      return partitioner->local_size();
+    }
+
+
+
     template <typename Number, typename MemorySpace>
     inline bool
     Vector<Number, MemorySpace>::in_local_range(
index d3b4cb91e6c49784fcb68a779784795a2ee1844d..5ec8af78887f14397bbdec08b602d62752ccde0d 100644 (file)
@@ -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
index 1df77dd799144e06a2b538bb0a7201c8bcf84f79..1ce0fac946311440e95f1e90fb690e13bc4e7f8d 100644 (file)
@@ -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 <typename Number>
   inline typename ReadWriteVector<Number>::size_type
-  ReadWriteVector<Number>::local_size() const
+  ReadWriteVector<Number>::locally_owned_size() const
   {
     return stored_elements.n_elements();
   }
index fc563239d5dd7562077332154d0e8d5188bb9ff9..639becf5f4257acad799a3c9ebf1835b54a32b1e 100644 (file)
@@ -494,7 +494,7 @@ namespace LinearAlgebra
       VecGetArray(static_cast<const Vec &>(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
index 421c7865c2641cee92c0fd5d5488c0fe2f0c1559..b840c3da508428ca2b229255e62a3050cbd6d2a7 100644 (file)
@@ -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.
index 729a6d480648e04c3775d85704d592a9008c5f91..91be6bf4f497132c70969dee00ebe8c5b5f386fa 100644 (file)
@@ -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.
index 55b04638a2351abbb4cb40735d42e541e983c213..a6bcd44a1f55bd068904cff0719013d34fab730a 100644 (file)
@@ -550,7 +550,7 @@ namespace LinearAlgebra
 
     template <typename Number>
     typename Vector<Number>::size_type
-    Vector<Number>::local_size() const
+    Vector<Number>::locally_owned_size() const
     {
       return vector->getLocalLength();
     }
index 7daa481a00937ba5f4ed7d19472eabb55423eca3..d46b23f4d03b3e0cfe6458f23422734305dd2af0 100644 (file)
@@ -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::size_type, Vector::size_type>
     Vector::local_range() const
     {
index de32e87ddc0eb8c7f2ccd390d0930a1052e40dc0..d7c8eed32d32c80ab39134f046519cbc456d00c8 100644 (file)
@@ -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<Number>::size() const
 
 template <typename Number>
 inline typename Vector<Number>::size_type
-Vector<Number>::local_size() const
+Vector<Number>::locally_owned_size() const
 {
   return values.size();
 }
index 75608c295348aaa4296bf235d56f2d52dcd936e0..3b9cf737fa355f21f2d0ba3e4388a26c127a1ada 100644 (file)
@@ -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)
index e45f395f88735ff699553d1949787990acf7706f..dd506269547ae5c0d3d50237066996c2d567a6e6 100644 (file)
@@ -541,7 +541,7 @@ namespace LinearAlgebra
 
 
     Vector::size_type
-    Vector::local_size() const
+    Vector::locally_owned_size() const
     {
       return vector->MyLength();
     }
index 1fcb51bfa4b14529f99287db8df012f178dd44f2..e2c4ca0eb874dc01629050acf4d26303323a0359 100644 (file)
@@ -13,7 +13,7 @@
 //
 // ---------------------------------------------------------------------
 
-// check Vector::local_size() for all supported vector types
+// check Vector::locally_owned_size() for all supported vector types
 
 #include <deal.II/base/index_set.h>
 #include <deal.II/base/utilities.h>
@@ -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<IndexSet> local_blocks {local_indices, local_indices};
+  std::vector<IndexSet> local_blocks{local_indices, local_indices};
   // for variety do not ghost the second component
-  std::vector<IndexSet> ghost_blocks {ghost_indices, IndexSet()};
+  std::vector<IndexSet> 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<TrilinosWrappers::MPI::Vector>();
 
   // block vectors:
-  check_ghosted_parallel_block<LinearAlgebra::distributed::BlockVector<double>>();
+  check_ghosted_parallel_block<
+    LinearAlgebra::distributed::BlockVector<double>>();
   check_ghosted_parallel_block<PETScWrappers::MPI::BlockVector>();
   check_ghosted_parallel_block<TrilinosWrappers::MPI::BlockVector>();
 }
index 622b61fcb98acf02e377c043e5bc5182f8f4f976..91e505e0aacbf4c8cceb40e7de31da9d690f5f50 100644 (file)
@@ -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<double>
 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<double>
@@ -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<double>
 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<double>
 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<double>
 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
 

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.