]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Get rid of VectorSpaceVector in the LinearAlgebra::distributed::BlockVector class.
authorWolfgang Bangerth <bangerth@colostate.edu>
Mon, 3 Jul 2023 21:11:23 +0000 (15:11 -0600)
committerWolfgang Bangerth <bangerth@colostate.edu>
Tue, 4 Jul 2023 01:48:45 +0000 (19:48 -0600)
include/deal.II/lac/la_parallel_block_vector.h
include/deal.II/lac/la_parallel_block_vector.templates.h

index f0946bfb0f64cae1a75bb251de41bded17deb8e1..693a9fe5ebce432df64aa2a80eb1481f4a8b3c46 100644 (file)
@@ -82,8 +82,7 @@ namespace LinearAlgebra
      * @ref GlossBlockLA "Block (linear algebra)"
      */
     template <typename Number>
-    class BlockVector : public BlockVectorBase<Vector<Number>>,
-                        public VectorSpaceVector<Number>
+    class BlockVector : public BlockVectorBase<Vector<Number>>
     {
     public:
       /**
@@ -200,14 +199,14 @@ namespace LinearAlgebra
        *   in a different section. The Intel compiler is prone to this
        *   behavior.
        */
-      virtual ~BlockVector() override = default;
+      ~BlockVector() = default;
 
       /**
        * Copy operator: fill all components of the vector with the given
        * scalar value.
        */
-      virtual BlockVector &
-      operator=(const value_type s) override;
+      BlockVector &
+      operator=(const value_type s);
 
       /**
        * Copy operator for arguments of the same type. Resize the present
@@ -399,8 +398,8 @@ namespace LinearAlgebra
        * the value on the owning processor and an exception is thrown if these
        * elements do not agree.
        */
-      virtual void
-      compress(VectorOperation::values operation) override;
+      void
+      compress(VectorOperation::values operation);
 
       /**
        * Fills the data field for ghost indices with the values stored in the
@@ -475,14 +474,14 @@ namespace LinearAlgebra
        * This function is mainly for internal consistency checks and should
        * seldom be used when not in debug mode since it uses quite some time.
        */
-      virtual bool
-      all_zero() const override;
+      bool
+      all_zero() const;
 
       /**
        * Compute the mean value of all the entries in the vector.
        */
-      virtual Number
-      mean_value() const override;
+      Number
+      mean_value() const;
 
       /**
        * $l_p$-norm of the vector. The pth root of the sum of the pth powers
@@ -513,41 +512,33 @@ namespace LinearAlgebra
       /** @} */
 
       /**
-       * @name 2: Implementation of VectorSpaceVector
+       * @name 2: Vector space operations
        */
       /** @{ */
 
-      /**
-       * Change the dimension to that of the vector V. The elements of V are not
-       * copied.
-       */
-      virtual void
-      reinit(const VectorSpaceVector<Number> &V,
-             const bool omit_zeroing_entries = false) override;
-
       /**
        * Multiply the entire vector by a fixed factor.
        */
-      virtual BlockVector<Number> &
-      operator*=(const Number factor) override;
+      BlockVector<Number> &
+      operator*=(const Number factor);
 
       /**
        * Divide the entire vector by a fixed factor.
        */
-      virtual BlockVector<Number> &
-      operator/=(const Number factor) override;
+      BlockVector<Number> &
+      operator/=(const Number factor);
 
       /**
        * Add the vector @p V to the present one.
        */
-      virtual BlockVector<Number> &
-      operator+=(const VectorSpaceVector<Number> &V) override;
+      BlockVector<Number> &
+      operator+=(const BlockVector<Number> &V);
 
       /**
        * Subtract the vector @p V from the present one.
        */
-      virtual BlockVector<Number> &
-      operator-=(const VectorSpaceVector<Number> &V) override;
+      BlockVector<Number> &
+      operator-=(const BlockVector<Number> &V);
 
       /**
        * Import all the elements present in the vector's IndexSet from the input
@@ -557,21 +548,21 @@ namespace LinearAlgebra
        * communication pattern is used multiple times. This can be used to
        * improve performance.
        */
-      virtual void
+      void
       import_elements(
         const LinearAlgebra::ReadWriteVector<Number> &V,
         VectorOperation::values                       operation,
         std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
-          communication_pattern = {}) override;
+          communication_pattern = {});
 
       /**
        * @deprecated Use import_elements() instead.
        */
-      DEAL_II_DEPRECATED virtual void
+      DEAL_II_DEPRECATED void
       import(const LinearAlgebra::ReadWriteVector<Number> &V,
              VectorOperation::values                       operation,
              std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
-               communication_pattern = {}) override
+               communication_pattern = {})
       {
         import_elements(V, operation, communication_pattern);
       }
@@ -579,8 +570,8 @@ namespace LinearAlgebra
       /**
        * Return the scalar product of two vectors.
        */
-      virtual Number
-      operator*(const VectorSpaceVector<Number> &V) const override;
+      Number
+      operator*(const BlockVector<Number> &V) const;
 
       /**
        * Calculate the scalar product between each block of this vector and @p V
@@ -645,29 +636,29 @@ namespace LinearAlgebra
       /**
        * Add @p a to all components. Note that @p a is a scalar not a vector.
        */
-      virtual void
-      add(const Number a) override;
+      void
+      add(const Number a);
 
       /**
        * Simple addition of a multiple of a vector, i.e. <tt>*this += a*V</tt>.
        */
-      virtual void
-      add(const Number a, const VectorSpaceVector<Number> &V) override;
+      void
+      add(const Number a, const BlockVector<Number> &V);
 
       /**
        * Multiple addition of scaled vectors, i.e. <tt>*this += a*V+b*W</tt>.
        */
-      virtual void
-      add(const Number                     a,
-          const VectorSpaceVector<Number> &V,
-          const Number                     b,
-          const VectorSpaceVector<Number> &W) override;
+      void
+      add(const Number               a,
+          const BlockVector<Number> &V,
+          const Number               b,
+          const BlockVector<Number> &W);
 
       /**
        * A collective add operation: This function adds a whole set of values
        * stored in @p values to the vector components specified by @p indices.
        */
-      virtual void
+      void
       add(const std::vector<size_type> &indices,
           const std::vector<Number> &   values);
 
@@ -675,38 +666,36 @@ namespace LinearAlgebra
        * Scaling and simple addition of a multiple of a vector, i.e. <tt>*this =
        * s*(*this)+a*V</tt>.
        */
-      virtual void
-      sadd(const Number                     s,
-           const Number                     a,
-           const VectorSpaceVector<Number> &V) override;
+      void
+      sadd(const Number s, const Number a, const BlockVector<Number> &V);
 
       /**
        * Scale each element of this vector by the corresponding element in the
        * argument. This function is mostly meant to simulate multiplication (and
        * immediate re-assignment) by a diagonal scaling matrix.
        */
-      virtual void
-      scale(const VectorSpaceVector<Number> &scaling_factors) override;
+      void
+      scale(const BlockVector<Number> &scaling_factors);
 
       /**
        * Assignment <tt>*this = a*V</tt>.
        */
-      virtual void
-      equ(const Number a, const VectorSpaceVector<Number> &V) override;
+      void
+      equ(const Number a, const BlockVector<Number> &V);
 
       /**
        * Return the l<sub>1</sub> norm of the vector (i.e., the sum of the
        * absolute values of all entries among all processors).
        */
-      virtual real_type
-      l1_norm() const override;
+      real_type
+      l1_norm() const;
 
       /**
        * Return the $l_2$ norm of the vector (i.e., the square root of
        * the sum of the square of all entries among all processors).
        */
-      virtual real_type
-      l2_norm() const override;
+      real_type
+      l2_norm() const;
 
       /**
        * Return the square of the $l_2$ norm of the vector.
@@ -718,8 +707,8 @@ namespace LinearAlgebra
        * Return the maximum norm of the vector (i.e., the maximum absolute value
        * among all entries and among all processors).
        */
-      virtual real_type
-      linfty_norm() const override;
+      real_type
+      linfty_norm() const;
 
       /**
        * Perform a combined operation of a vector addition and a subsequent
@@ -741,10 +730,10 @@ namespace LinearAlgebra
        * implemented as
        * $\left<v,w\right>=\sum_i v_i \bar{w_i}$.
        */
-      virtual Number
-      add_and_dot(const Number                     a,
-                  const VectorSpaceVector<Number> &V,
-                  const VectorSpaceVector<Number> &W) override;
+      Number
+      add_and_dot(const Number               a,
+                  const BlockVector<Number> &V,
+                  const BlockVector<Number> &W);
 
       /**
        * Return the global size of the vector, equal to the sum of the number of
@@ -764,23 +753,23 @@ namespace LinearAlgebra
        *  vec.locally_owned_elements() == complete_index_set(vec.size())
        * @endcode
        */
-      virtual dealii::IndexSet
-      locally_owned_elements() const override;
+      dealii::IndexSet
+      locally_owned_elements() const;
 
       /**
        * Print the vector to the output stream @p out.
        */
-      virtual void
+      void
       print(std::ostream &     out,
             const unsigned int precision  = 3,
             const bool         scientific = true,
-            const bool         across     = true) const override;
+            const bool         across     = true) const;
 
       /**
        * Return the memory consumption of this class in bytes.
        */
-      virtual std::size_t
-      memory_consumption() const override;
+      std::size_t
+      memory_consumption() const;
       /** @} */
 
       /**
index 923c297a886656b6a0270e1cc8ccd56892988e0a..d4155126b5490f7151f5958cad01e1fba116931d 100644 (file)
@@ -503,19 +503,6 @@ namespace LinearAlgebra
 
 
 
-    template <typename Number>
-    void
-    BlockVector<Number>::reinit(const VectorSpaceVector<Number> &V,
-                                const bool omit_zeroing_entries)
-    {
-      Assert(dynamic_cast<const BlockVector<Number> *>(&V) != nullptr,
-             ExcVectorTypeNotCompatible());
-      const BlockVector<Number> &down_V =
-        dynamic_cast<const BlockVector<Number> &>(V);
-      reinit(down_V, omit_zeroing_entries);
-    }
-
-
     template <typename Number>
     BlockVector<Number> &
     BlockVector<Number>::operator*=(const Number factor)
@@ -539,13 +526,8 @@ namespace LinearAlgebra
 
     template <typename Number>
     void
-    BlockVector<Number>::scale(const VectorSpaceVector<Number> &vv)
+    BlockVector<Number>::scale(const BlockVector<Number> &v)
     {
-      // Downcast. Throws an exception if invalid.
-      Assert(dynamic_cast<const BlockVector<Number> *>(&vv) != nullptr,
-             ExcVectorTypeNotCompatible());
-      const BlockVector<Number> &v =
-        dynamic_cast<const BlockVector<Number> &>(vv);
       AssertDimension(this->n_blocks(), v.n_blocks());
       for (unsigned int block = 0; block < this->n_blocks(); ++block)
         this->block(block).scale(v.block(block));
@@ -555,14 +537,8 @@ namespace LinearAlgebra
 
     template <typename Number>
     void
-    BlockVector<Number>::equ(const Number                     a,
-                             const VectorSpaceVector<Number> &vv)
+    BlockVector<Number>::equ(const Number a, const BlockVector<Number> &v)
     {
-      // Downcast. Throws an exception if invalid.
-      Assert(dynamic_cast<const BlockVector<Number> *>(&vv) != nullptr,
-             ExcVectorTypeNotCompatible());
-      const BlockVector<Number> &v =
-        dynamic_cast<const BlockVector<Number> &>(vv);
       AssertDimension(this->n_blocks(), v.n_blocks());
       for (unsigned int block = 0; block < this->n_blocks(); ++block)
         this->block(block).equ(a, v.block(block));
@@ -572,13 +548,8 @@ namespace LinearAlgebra
 
     template <typename Number>
     BlockVector<Number> &
-    BlockVector<Number>::operator+=(const VectorSpaceVector<Number> &vv)
+    BlockVector<Number>::operator+=(const BlockVector<Number> &v)
     {
-      // Downcast. Throws an exception if invalid.
-      Assert(dynamic_cast<const BlockVector<Number> *>(&vv) != nullptr,
-             ExcVectorTypeNotCompatible());
-      const BlockVector<Number> &v =
-        dynamic_cast<const BlockVector<Number> &>(vv);
       AssertDimension(this->n_blocks(), v.n_blocks());
       for (unsigned int block = 0; block < this->n_blocks(); ++block)
         this->block(block) += v.block(block);
@@ -590,13 +561,8 @@ namespace LinearAlgebra
 
     template <typename Number>
     BlockVector<Number> &
-    BlockVector<Number>::operator-=(const VectorSpaceVector<Number> &vv)
+    BlockVector<Number>::operator-=(const BlockVector<Number> &v)
     {
-      // Downcast. Throws an exception if invalid.
-      Assert(dynamic_cast<const BlockVector<Number> *>(&vv) != nullptr,
-             ExcVectorTypeNotCompatible());
-      const BlockVector<Number> &v =
-        dynamic_cast<const BlockVector<Number> &>(vv);
       AssertDimension(this->n_blocks(), v.n_blocks());
       for (unsigned int block = 0; block < this->n_blocks(); ++block)
         this->block(block) -= v.block(block);
@@ -618,14 +584,8 @@ namespace LinearAlgebra
 
     template <typename Number>
     void
-    BlockVector<Number>::add(const Number                     a,
-                             const VectorSpaceVector<Number> &vv)
+    BlockVector<Number>::add(const Number a, const BlockVector<Number> &v)
     {
-      // Downcast. Throws an exception if invalid.
-      Assert(dynamic_cast<const BlockVector<Number> *>(&vv) != nullptr,
-             ExcVectorTypeNotCompatible());
-      const BlockVector<Number> &v =
-        dynamic_cast<const BlockVector<Number> &>(vv);
       AssertDimension(this->n_blocks(), v.n_blocks());
       for (unsigned int block = 0; block < this->n_blocks(); ++block)
         this->block(block).add(a, v.block(block));
@@ -635,22 +595,13 @@ namespace LinearAlgebra
 
     template <typename Number>
     void
-    BlockVector<Number>::add(const Number                     a,
-                             const VectorSpaceVector<Number> &vv,
-                             const Number                     b,
-                             const VectorSpaceVector<Number> &ww)
+    BlockVector<Number>::add(const Number               a,
+                             const BlockVector<Number> &v,
+                             const Number               b,
+                             const BlockVector<Number> &w)
     {
-      // Downcast. Throws an exception if invalid.
-      Assert(dynamic_cast<const BlockVector<Number> *>(&vv) != nullptr,
-             ExcVectorTypeNotCompatible());
-      const BlockVector<Number> &v =
-        dynamic_cast<const BlockVector<Number> &>(vv);
-      AssertDimension(this->n_blocks(), v.n_blocks());
-      Assert(dynamic_cast<const BlockVector<Number> *>(&ww) != nullptr,
-             ExcVectorTypeNotCompatible());
-      const BlockVector<Number> &w =
-        dynamic_cast<const BlockVector<Number> &>(ww);
       AssertDimension(this->n_blocks(), v.n_blocks());
+      AssertDimension(this->n_blocks(), w.n_blocks());
 
       for (unsigned int block = 0; block < this->n_blocks(); ++block)
         this->block(block).add(a, v.block(block), b, w.block(block));
@@ -660,15 +611,10 @@ namespace LinearAlgebra
 
     template <typename Number>
     void
-    BlockVector<Number>::sadd(const Number                     x,
-                              const Number                     a,
-                              const VectorSpaceVector<Number> &vv)
+    BlockVector<Number>::sadd(const Number               x,
+                              const Number               a,
+                              const BlockVector<Number> &v)
     {
-      // Downcast. Throws an exception if invalid.
-      Assert(dynamic_cast<const BlockVector<Number> *>(&vv) != nullptr,
-             ExcVectorTypeNotCompatible());
-      const BlockVector<Number> &v =
-        dynamic_cast<const BlockVector<Number> &>(vv);
       AssertDimension(this->n_blocks(), v.n_blocks());
       for (unsigned int block = 0; block < this->n_blocks(); ++block)
         this->block(block).sadd(x, a, v.block(block));
@@ -734,15 +680,9 @@ namespace LinearAlgebra
 
     template <typename Number>
     Number
-    BlockVector<Number>::operator*(const VectorSpaceVector<Number> &vv) const
+    BlockVector<Number>::operator*(const BlockVector<Number> &v) const
     {
       Assert(this->n_blocks() > 0, ExcEmptyObject());
-
-      // Downcast. Throws an exception if invalid.
-      Assert(dynamic_cast<const BlockVector<Number> *>(&vv) != nullptr,
-             ExcVectorTypeNotCompatible());
-      const BlockVector<Number> &v =
-        dynamic_cast<const BlockVector<Number> &>(vv);
       AssertDimension(this->n_blocks(), v.n_blocks());
 
       Number local_result = Number();
@@ -871,22 +811,11 @@ namespace LinearAlgebra
 
     template <typename Number>
     inline Number
-    BlockVector<Number>::add_and_dot(const Number                     a,
-                                     const VectorSpaceVector<Number> &vv,
-                                     const VectorSpaceVector<Number> &ww)
+    BlockVector<Number>::add_and_dot(const Number               a,
+                                     const BlockVector<Number> &v,
+                                     const BlockVector<Number> &w)
     {
-      // Downcast. Throws an exception if invalid.
-      Assert(dynamic_cast<const BlockVector<Number> *>(&vv) != nullptr,
-             ExcVectorTypeNotCompatible());
-      const BlockVector<Number> &v =
-        dynamic_cast<const BlockVector<Number> &>(vv);
       AssertDimension(this->n_blocks(), v.n_blocks());
-
-      // Downcast. Throws an exception if invalid.
-      Assert(dynamic_cast<const BlockVector<Number> *>(&ww) != nullptr,
-             ExcVectorTypeNotCompatible());
-      const BlockVector<Number> &w =
-        dynamic_cast<const BlockVector<Number> &>(ww);
       AssertDimension(this->n_blocks(), w.n_blocks());
 
       Number local_result = Number();

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.