From: Marc Fehling <mafehling.git@gmail.com>
Date: Sat, 26 Nov 2022 02:23:51 +0000 (-0700)
Subject: Simplify reinit() of PETSc/Trilinos BlockVector.
X-Git-Tag: v9.5.0-rc1~796^2~1
X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=931b47237e306a143524f473239d82b2fc978f62;p=dealii.git

Simplify reinit() of PETSc/Trilinos BlockVector.
---

diff --git a/include/deal.II/lac/petsc_block_vector.h b/include/deal.II/lac/petsc_block_vector.h
index 97c80245a6..ebe458fbe1 100644
--- a/include/deal.II/lac/petsc_block_vector.h
+++ b/include/deal.II/lac/petsc_block_vector.h
@@ -312,11 +312,13 @@ namespace PETScWrappers
     inline BlockVector::BlockVector(const BlockVector &v)
       : BlockVectorBase<Vector>()
     {
-      this->components.resize(v.n_blocks());
       this->block_indices = v.block_indices;
 
+      this->components.resize(this->n_blocks());
       for (unsigned int i = 0; i < this->n_blocks(); ++i)
         this->components[i] = v.components[i];
+
+      this->collect_sizes();
     }
 
     inline BlockVector::BlockVector(
@@ -346,16 +348,17 @@ namespace PETScWrappers
     {
       // we only allow assignment to vectors with the same number of blocks
       // or to an empty BlockVector
-      Assert(n_blocks() == 0 || n_blocks() == v.n_blocks(),
-             ExcDimensionMismatch(n_blocks(), v.n_blocks()));
+      Assert(this->n_blocks() == 0 || this->n_blocks() == v.n_blocks(),
+             ExcDimensionMismatch(this->n_blocks(), v.n_blocks()));
 
       if (this->n_blocks() != v.n_blocks())
-        reinit(v.n_blocks());
+        this->block_indices = v.block_indices;
 
-      for (size_type i = 0; i < this->n_blocks(); ++i)
-        this->components[i] = v.block(i);
+      this->components.resize(this->n_blocks());
+      for (unsigned int i = 0; i < this->n_blocks(); ++i)
+        this->components[i] = v.components[i];
 
-      collect_sizes();
+      this->collect_sizes();
 
       return *this;
     }
@@ -384,42 +387,45 @@ namespace PETScWrappers
                         const bool                    omit_zeroing_entries)
     {
       this->block_indices.reinit(block_sizes);
-      if (this->components.size() != this->n_blocks())
-        this->components.resize(this->n_blocks());
 
+      this->components.resize(this->n_blocks());
       for (unsigned int i = 0; i < this->n_blocks(); ++i)
         this->components[i].reinit(communicator,
                                    block_sizes[i],
                                    locally_owned_sizes[i],
                                    omit_zeroing_entries);
+
+      this->collect_sizes();
     }
 
 
     inline void
     BlockVector::reinit(const BlockVector &v, const bool omit_zeroing_entries)
     {
-      this->block_indices = v.get_block_indices();
-      if (this->components.size() != this->n_blocks())
-        this->components.resize(this->n_blocks());
+      if (this->n_blocks() != v.n_blocks())
+        this->block_indices = v.get_block_indices();
 
+      this->components.resize(this->n_blocks());
       for (unsigned int i = 0; i < this->n_blocks(); ++i)
-        block(i).reinit(v.block(i), omit_zeroing_entries);
+        this->components[i].reinit(v.components[i], omit_zeroing_entries);
+
+      this->collect_sizes();
     }
 
     inline void
     BlockVector::reinit(const std::vector<IndexSet> &parallel_partitioning,
                         const MPI_Comm &             communicator)
     {
-      std::vector<size_type> sizes(parallel_partitioning.size());
-      for (unsigned int i = 0; i < parallel_partitioning.size(); ++i)
-        sizes[i] = parallel_partitioning[i].size();
-
-      this->block_indices.reinit(sizes);
-      if (this->components.size() != this->n_blocks())
-        this->components.resize(this->n_blocks());
+      // update the number of blocks
+      this->block_indices.reinit(parallel_partitioning.size(), 0);
 
+      // initialize each block
+      this->components.resize(this->n_blocks());
       for (unsigned int i = 0; i < this->n_blocks(); ++i)
-        block(i).reinit(parallel_partitioning[i], communicator);
+        this->components[i].reinit(parallel_partitioning[i], communicator);
+
+      // update block_indices content
+      this->collect_sizes();
     }
 
     inline void
@@ -427,18 +433,20 @@ namespace PETScWrappers
                         const std::vector<IndexSet> &ghost_entries,
                         const MPI_Comm &             communicator)
     {
-      std::vector<types::global_dof_index> sizes(parallel_partitioning.size());
-      for (unsigned int i = 0; i < parallel_partitioning.size(); ++i)
-        sizes[i] = parallel_partitioning[i].size();
+      AssertDimension(parallel_partitioning.size(), ghost_entries.size());
 
-      this->block_indices.reinit(sizes);
-      if (this->components.size() != this->n_blocks())
-        this->components.resize(this->n_blocks());
+      // update the number of blocks
+      this->block_indices.reinit(parallel_partitioning.size(), 0);
 
+      // initialize each block
+      this->components.resize(this->n_blocks());
       for (unsigned int i = 0; i < this->n_blocks(); ++i)
-        block(i).reinit(parallel_partitioning[i],
-                        ghost_entries[i],
-                        communicator);
+        this->components[i].reinit(parallel_partitioning[i],
+                                   ghost_entries[i],
+                                   communicator);
+
+      // update block_indices content
+      this->collect_sizes();
     }
 
 
diff --git a/include/deal.II/lac/trilinos_parallel_block_vector.h b/include/deal.II/lac/trilinos_parallel_block_vector.h
index 3d32bbbd9d..fe42d5a113 100644
--- a/include/deal.II/lac/trilinos_parallel_block_vector.h
+++ b/include/deal.II/lac/trilinos_parallel_block_vector.h
@@ -348,11 +348,13 @@ namespace TrilinosWrappers
     inline BlockVector::BlockVector(const BlockVector &v)
       : dealii::BlockVectorBase<MPI::Vector>()
     {
-      this->components.resize(v.n_blocks());
       this->block_indices = v.block_indices;
 
-      for (size_type i = 0; i < this->n_blocks(); ++i)
+      this->components.resize(this->n_blocks());
+      for (unsigned int i = 0; i < this->n_blocks(); ++i)
         this->components[i] = v.components[i];
+
+      this->collect_sizes();
     }
 
 
@@ -370,18 +372,19 @@ namespace TrilinosWrappers
     BlockVector &
     BlockVector::operator=(const ::dealii::BlockVector<Number> &v)
     {
-      if (n_blocks() != v.n_blocks())
-        {
-          std::vector<size_type> block_sizes(v.n_blocks(), 0);
-          block_indices.reinit(block_sizes);
-          if (components.size() != n_blocks())
-            components.resize(n_blocks());
-        }
-
-      for (size_type i = 0; i < this->n_blocks(); ++i)
+      // we only allow assignment to vectors with the same number of blocks
+      // or to an empty BlockVector
+      Assert(this->n_blocks() == 0 || this->n_blocks() == v.n_blocks(),
+             ExcDimensionMismatch(this->n_blocks(), v.n_blocks()));
+
+      if (this->n_blocks() != v.n_blocks())
+        this->block_indices = v.get_block_indices();
+
+      this->components.resize(this->n_blocks());
+      for (unsigned int i = 0; i < this->n_blocks(); ++i)
         this->components[i] = v.block(i);
 
-      collect_sizes();
+      this->collect_sizes();
 
       return *this;
     }
diff --git a/source/lac/trilinos_block_vector.cc b/source/lac/trilinos_block_vector.cc
index cd19d600bb..486fa0c418 100644
--- a/source/lac/trilinos_block_vector.cc
+++ b/source/lac/trilinos_block_vector.cc
@@ -41,16 +41,17 @@ namespace TrilinosWrappers
     {
       // we only allow assignment to vectors with the same number of blocks
       // or to an empty BlockVector
-      Assert(n_blocks() == 0 || n_blocks() == v.n_blocks(),
-             ExcDimensionMismatch(n_blocks(), v.n_blocks()));
+      Assert(this->n_blocks() == 0 || this->n_blocks() == v.n_blocks(),
+             ExcDimensionMismatch(this->n_blocks(), v.n_blocks()));
 
       if (this->n_blocks() != v.n_blocks())
-        reinit(v.n_blocks());
+        this->block_indices = v.block_indices;
 
-      for (size_type i = 0; i < this->n_blocks(); ++i)
-        this->components[i] = v.block(i);
+      this->components.resize(this->n_blocks());
+      for (unsigned int i = 0; i < this->n_blocks(); ++i)
+        this->components[i] = v.components[i];
 
-      collect_sizes();
+      this->collect_sizes();
 
       return *this;
     }
@@ -71,24 +72,18 @@ namespace TrilinosWrappers
                         const MPI_Comm &             communicator,
                         const bool                   omit_zeroing_entries)
     {
-      const size_type        no_blocks = parallel_partitioning.size();
-      std::vector<size_type> block_sizes(no_blocks);
-
-      for (size_type i = 0; i < no_blocks; ++i)
-        {
-          block_sizes[i] = parallel_partitioning[i].size();
-        }
-
-      this->block_indices.reinit(block_sizes);
-      if (components.size() != n_blocks())
-        components.resize(n_blocks());
-
-      for (size_type i = 0; i < n_blocks(); ++i)
-        components[i].reinit(parallel_partitioning[i],
-                             communicator,
-                             omit_zeroing_entries);
-
-      collect_sizes();
+      // update the number of blocks
+      this->block_indices.reinit(parallel_partitioning.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(parallel_partitioning[i],
+                                   communicator,
+                                   omit_zeroing_entries);
+
+      // update block_indices content
+      this->collect_sizes();
     }
 
 
@@ -99,25 +94,21 @@ namespace TrilinosWrappers
                         const MPI_Comm &             communicator,
                         const bool                   vector_writable)
     {
-      const size_type        no_blocks = parallel_partitioning.size();
-      std::vector<size_type> block_sizes(no_blocks);
-
-      for (size_type i = 0; i < no_blocks; ++i)
-        {
-          block_sizes[i] = parallel_partitioning[i].size();
-        }
+      AssertDimension(parallel_partitioning.size(), ghost_values.size());
 
-      this->block_indices.reinit(block_sizes);
-      if (components.size() != n_blocks())
-        components.resize(n_blocks());
+      // update the number of blocks
+      this->block_indices.reinit(parallel_partitioning.size(), 0);
 
-      for (size_type i = 0; i < n_blocks(); ++i)
-        components[i].reinit(parallel_partitioning[i],
-                             ghost_values[i],
-                             communicator,
-                             vector_writable);
+      // initialize each block
+      this->components.resize(this->n_blocks());
+      for (unsigned int i = 0; i < this->n_blocks(); ++i)
+        this->components[i].reinit(parallel_partitioning[i],
+                                   ghost_values[i],
+                                   communicator,
+                                   vector_writable);
 
-      collect_sizes();
+      // update block_indices content
+      this->collect_sizes();
     }
 
 
@@ -125,14 +116,16 @@ namespace TrilinosWrappers
     void
     BlockVector::reinit(const BlockVector &v, const bool omit_zeroing_entries)
     {
-      block_indices = v.get_block_indices();
-      if (components.size() != n_blocks())
-        components.resize(n_blocks());
+      if (this->n_blocks() != v.n_blocks())
+        this->block_indices = v.get_block_indices();
 
-      for (size_type i = 0; i < n_blocks(); ++i)
-        components[i].reinit(v.block(i), omit_zeroing_entries, false);
+      this->components.resize(this->n_blocks());
+      for (unsigned int i = 0; i < this->n_blocks(); ++i)
+        this->components[i].reinit(v.components[i],
+                                   omit_zeroing_entries,
+                                   false);
 
-      collect_sizes();
+      this->collect_sizes();
     }
 
 
@@ -140,15 +133,11 @@ namespace TrilinosWrappers
     void
     BlockVector::reinit(const size_type num_blocks)
     {
-      std::vector<size_type> block_sizes(num_blocks, 0);
-      this->block_indices.reinit(block_sizes);
-      if (this->components.size() != this->n_blocks())
-        this->components.resize(this->n_blocks());
-
-      for (size_type i = 0; i < this->n_blocks(); ++i)
-        components[i].clear();
+      this->block_indices.reinit(num_blocks, 0);
 
-      collect_sizes();
+      this->components.resize(this->n_blocks());
+      for (auto &c : this->components)
+        c.clear();
     }
 
 
@@ -158,21 +147,18 @@ namespace TrilinosWrappers
       const TrilinosWrappers::BlockSparseMatrix &m,
       const BlockVector &                        v)
     {
-      Assert(m.n_block_rows() == v.n_blocks(),
-             ExcDimensionMismatch(m.n_block_rows(), v.n_blocks()));
-      Assert(m.n_block_cols() == v.n_blocks(),
-             ExcDimensionMismatch(m.n_block_cols(), v.n_blocks()));
+      AssertDimension(m.n_block_rows(), v.n_blocks());
+      AssertDimension(m.n_block_cols(), v.n_blocks());
 
-      if (v.n_blocks() != n_blocks())
-        {
-          block_indices = v.get_block_indices();
-          components.resize(v.n_blocks());
-        }
+      if (this->n_blocks() != v.n_blocks())
+        this->block_indices = v.get_block_indices();
 
-      for (size_type i = 0; i < this->n_blocks(); ++i)
-        components[i].import_nonlocal_data_for_fe(m.block(i, i), v.block(i));
+      this->components.resize(this->n_blocks());
+      for (unsigned int i = 0; i < this->n_blocks(); ++i)
+        this->components[i].import_nonlocal_data_for_fe(m.block(i, i),
+                                                        v.block(i));
 
-      collect_sizes();
+      this->collect_sizes();
     }
 
 
@@ -183,7 +169,7 @@ namespace TrilinosWrappers
                        const bool         scientific,
                        const bool         across) const
     {
-      for (size_type i = 0; i < this->n_blocks(); ++i)
+      for (unsigned int i = 0; i < this->n_blocks(); ++i)
         {
           if (across)
             out << 'C' << i << ':';