]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Get rid of VectorSpaceVector in the LinearAlgebra::TpetraWrappers::Vector class. 15631/head
authorWolfgang Bangerth <bangerth@colostate.edu>
Mon, 3 Jul 2023 20:41:56 +0000 (14:41 -0600)
committerWolfgang Bangerth <bangerth@colostate.edu>
Mon, 3 Jul 2023 20:48:24 +0000 (14:48 -0600)
include/deal.II/lac/trilinos_tpetra_vector.h
include/deal.II/lac/trilinos_tpetra_vector.templates.h

index 3ad3cd01f2e85cdec57915f39ca63f580ca2db04..6cd715bd1a7b879e2ec5fdd0cce7b0486373f509 100644 (file)
@@ -95,8 +95,7 @@ namespace LinearAlgebra
   {
     /**
      * This class implements a wrapper to the Trilinos distributed vector
-     * class Tpetra::Vector. This class is derived from the
-     * LinearAlgebra::VectorSpaceVector class and requires Trilinos to be
+     * class Tpetra::Vector. This class requires Trilinos to be
      * compiled with MPI support.
      *
      * Tpetra uses Kokkos for thread-parallelism and chooses the execution and
@@ -114,14 +113,12 @@ namespace LinearAlgebra
      * @ingroup Vectors
      */
     template <typename Number>
-    class Vector : public VectorSpaceVector<Number>,
-                   public ReadVector<Number>,
-                   public Subscriptor
+    class Vector : public ReadVector<Number>, public Subscriptor
     {
     public:
       using value_type = Number;
-
-      using size_type = typename VectorSpaceVector<Number>::size_type;
+      using real_type  = typename numbers::NumberTraits<Number>::real_type;
+      using size_type  = types::global_dof_index;
 
       /**
        * Constructor. Create a vector of dimension zero.
@@ -158,9 +155,8 @@ namespace LinearAlgebra
        * Change the dimension to that of the vector @p V. The elements of @p V are not
        * copied.
        */
-      virtual void
-      reinit(const VectorSpaceVector<Number> &V,
-             const bool omit_zeroing_entries = false) override;
+      void
+      reinit(const Vector<Number> &V, const bool omit_zeroing_entries = false);
 
       /**
        * Extract a range of elements all at once.
@@ -182,8 +178,8 @@ namespace LinearAlgebra
        * Sets all elements of the vector to the scalar @p s. This operation is
        * only allowed if @p s is equal to zero.
        */
-      virtual Vector &
-      operator=(const Number s) override;
+      Vector &
+      operator=(const Number s);
 
       /**
        * Imports all the elements present in the vector's IndexSet from the
@@ -193,22 +189,22 @@ namespace LinearAlgebra
        * communication pattern is used multiple times. This can be used to
        * improve performance.
        */
-      virtual void
+      void
       import_elements(
         const 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
+      void
       import(const ReadWriteVector<Number> &V,
              VectorOperation::values        operation,
              std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
-               communication_pattern = {}) override
+               communication_pattern = {})
       {
         import_elements(V, operation, communication_pattern);
       }
@@ -216,65 +212,63 @@ namespace LinearAlgebra
       /**
        * Multiply the entire vector by a fixed factor.
        */
-      virtual Vector &
-      operator*=(const Number factor) override;
+      Vector &
+      operator*=(const Number factor);
 
       /**
        * Divide the entire vector by a fixed factor.
        */
-      virtual Vector &
-      operator/=(const Number factor) override;
+      Vector &
+      operator/=(const Number factor);
 
       /**
        * Add the vector @p V to the present one.
        */
-      virtual Vector &
-      operator+=(const VectorSpaceVector<Number> &V) override;
+      Vector &
+      operator+=(const Vector<Number> &V);
 
       /**
        * Subtract the vector @p V from the present one.
        */
-      virtual Vector &
-      operator-=(const VectorSpaceVector<Number> &V) override;
+      Vector &
+      operator-=(const Vector<Number> &V);
 
       /**
        * Return the scalar product of two vectors. The vectors need to have the
        * same layout.
        */
-      virtual Number
-      operator*(const VectorSpaceVector<Number> &V) const override;
+      Number
+      operator*(const Vector<Number> &V) const;
 
       /**
        * Add @p a to all components. Note that @p 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>. The vectors need to have the same layout.
        */
-      virtual void
-      add(const Number a, const VectorSpaceVector<Number> &V) override;
+      void
+      add(const Number a, const Vector<Number> &V);
 
       /**
        * Multiple addition of multiple of a vector, i.e. <tt>*this> +=
        * a*V+b*W</tt>. The vectors need to have the same layout.
        */
-      virtual void
-      add(const Number                     a,
-          const VectorSpaceVector<Number> &V,
-          const Number                     b,
-          const VectorSpaceVector<Number> &W) override;
+      void
+      add(const Number          a,
+          const Vector<Number> &V,
+          const Number          b,
+          const Vector<Number> &W);
 
       /**
        * 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 Vector<Number> &V);
 
       /**
        * Scale each element of this vector by the corresponding element in the
@@ -282,47 +276,47 @@ namespace LinearAlgebra
        * (and immediate re-assignment) by a diagonal scaling matrix. The
        * vectors need to have the same layout.
        */
-      virtual void
-      scale(const VectorSpaceVector<Number> &scaling_factors) override;
+      void
+      scale(const Vector<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 Vector<Number> &V);
 
       /**
        * Return whether the vector contains only elements with value zero.
        */
-      virtual bool
-      all_zero() const override;
+      bool
+      all_zero() const;
 
       /**
        * Return the mean value of the element of this vector.
        */
-      virtual Number
-      mean_value() const override;
+      Number
+      mean_value() const;
 
       /**
        * 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 typename LinearAlgebra::VectorSpaceVector<Number>::real_type
-      l1_norm() const override;
+      real_type
+      l1_norm() const;
 
       /**
        * Return the l<sub>2</sub> norm of the vector (i.e., the square root of
        * the sum of the square of all entries among all processors).
        */
-      virtual typename LinearAlgebra::VectorSpaceVector<Number>::real_type
-      l2_norm() const override;
+      real_type
+      l2_norm() const;
 
       /**
        * Return the maximum norm of the vector (i.e., the maximum absolute value
        * among all entries and among all processors).
        */
-      virtual typename LinearAlgebra::VectorSpaceVector<Number>::real_type
-      linfty_norm() const override;
+      real_type
+      linfty_norm() const;
 
       /**
        * Performs a combined operation of a vector addition and a subsequent
@@ -346,10 +340,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 Vector<Number> &V,
+                  const Vector<Number> &W);
       /**
        * This function always returns false and is present only for backward
        * compatibility.
@@ -388,9 +382,21 @@ 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;
 
+      /**
+       * Compress the underlying representation of the Trilinos object, i.e.
+       * flush the buffers of the vector object if it has any. This function is
+       * necessary after writing into a vector element-by-element and before
+       * anything else can be done on it.
+       *
+       * See
+       * @ref GlossCompress "Compressing distributed objects"
+       * for more information.
+       */
+      void
+      compress(const VectorOperation::values operation);
       /**
        * Return a const reference to the underlying Trilinos
        * Tpetra::Vector class.
@@ -408,17 +414,17 @@ namespace LinearAlgebra
       /**
        * Prints 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;
 
       /**
        * The vectors have different partitioning, i.e. their IndexSet objects
index e1f1dc06b5ef011f8b86e33510ca77c0b048111b..fc8fb5c05747386296c925bf16a8ccacc4a34266 100644 (file)
@@ -101,18 +101,11 @@ namespace LinearAlgebra
 
     template <typename Number>
     void
-    Vector<Number>::reinit(const VectorSpaceVector<Number> &V,
-                           const bool omit_zeroing_entries)
+    Vector<Number>::reinit(const Vector<Number> &V,
+                           const bool            omit_zeroing_entries)
     {
-      // Check that casting will work.
-      Assert(dynamic_cast<const Vector<Number> *>(&V) != nullptr,
-             ExcVectorTypeNotCompatible());
-
-      // Downcast V. If fails, throws an exception.
-      const Vector<Number> &down_V = dynamic_cast<const Vector<Number> &>(V);
-
-      reinit(down_V.locally_owned_elements(),
-             down_V.get_mpi_communicator(),
+      reinit(V.locally_owned_elements(),
+             V.get_mpi_communicator(),
              omit_zeroing_entries);
     }
 
@@ -296,34 +289,26 @@ namespace LinearAlgebra
 
     template <typename Number>
     Vector<Number> &
-    Vector<Number>::operator+=(const VectorSpaceVector<Number> &V)
+    Vector<Number>::operator+=(const Vector<Number> &V)
     {
-      // Check that casting will work.
-      Assert(dynamic_cast<const Vector<Number> *>(&V) != nullptr,
-             ExcVectorTypeNotCompatible());
-
-      // Downcast V. If fails, throws an exception.
-      const Vector<Number> &down_V = dynamic_cast<const Vector<Number> &>(V);
       // If the maps are the same we can update right away.
-      if (vector->getMap()->isSameAs(*(down_V.trilinos_vector().getMap())))
+      if (vector->getMap()->isSameAs(*(V.trilinos_vector().getMap())))
         {
-          vector->update(1., down_V.trilinos_vector(), 1.);
+          vector->update(1., V.trilinos_vector(), 1.);
         }
       else
         {
-          Assert(this->size() == down_V.size(),
-                 ExcDimensionMismatch(this->size(), down_V.size()));
+          Assert(this->size() == V.size(),
+                 ExcDimensionMismatch(this->size(), V.size()));
 
           // TODO: Tpetra doesn't have a combine mode that also updates local
           // elements, maybe there is a better workaround.
           Tpetra::Vector<Number, int, types::signed_global_dof_index> dummy(
             vector->getMap(), false);
           Tpetra::Import<int, types::signed_global_dof_index> data_exchange(
-            down_V.trilinos_vector().getMap(), dummy.getMap());
+            V.trilinos_vector().getMap(), dummy.getMap());
 
-          dummy.doImport(down_V.trilinos_vector(),
-                         data_exchange,
-                         Tpetra::INSERT);
+          dummy.doImport(V.trilinos_vector(), data_exchange, Tpetra::INSERT);
 
           vector->update(1.0, dummy, 1.0);
         }
@@ -335,7 +320,7 @@ namespace LinearAlgebra
 
     template <typename Number>
     Vector<Number> &
-    Vector<Number>::operator-=(const VectorSpaceVector<Number> &V)
+    Vector<Number>::operator-=(const Vector<Number> &V)
     {
       this->add(-1., V);
 
@@ -346,20 +331,14 @@ namespace LinearAlgebra
 
     template <typename Number>
     Number
-    Vector<Number>::operator*(const VectorSpaceVector<Number> &V) const
+    Vector<Number>::operator*(const Vector<Number> &V) const
     {
-      // Check that casting will work.
-      Assert(dynamic_cast<const Vector<Number> *>(&V) != nullptr,
-             ExcVectorTypeNotCompatible());
-
-      // Downcast V. If fails, throws an exception.
-      const Vector<Number> &down_V = dynamic_cast<const Vector<Number> &>(V);
-      Assert(this->size() == down_V.size(),
-             ExcDimensionMismatch(this->size(), down_V.size()));
-      Assert(vector->getMap()->isSameAs(*down_V.trilinos_vector().getMap()),
+      Assert(this->size() == V.size(),
+             ExcDimensionMismatch(this->size(), V.size()));
+      Assert(vector->getMap()->isSameAs(*V.trilinos_vector().getMap()),
              ExcDifferentParallelPartitioning());
 
-      return vector->dot(down_V.trilinos_vector());
+      return vector->dot(V.trilinos_vector());
     }
 
 
@@ -397,68 +376,45 @@ namespace LinearAlgebra
 
     template <typename Number>
     void
-    Vector<Number>::add(const Number a, const VectorSpaceVector<Number> &V)
+    Vector<Number>::add(const Number a, const Vector<Number> &V)
     {
-      // Check that casting will work.
-      Assert(dynamic_cast<const Vector<Number> *>(&V) != nullptr,
-             ExcVectorTypeNotCompatible());
-
-      // Downcast V. If fails, throws an exception.
-      const Vector<Number> &down_V = dynamic_cast<const Vector<Number> &>(V);
       AssertIsFinite(a);
-      Assert(vector->getMap()->isSameAs(*(down_V.trilinos_vector().getMap())),
+      Assert(vector->getMap()->isSameAs(*(V.trilinos_vector().getMap())),
              ExcDifferentParallelPartitioning());
 
-      vector->update(a, down_V.trilinos_vector(), 1.);
+      vector->update(a, V.trilinos_vector(), 1.);
     }
 
 
 
     template <typename Number>
     void
-    Vector<Number>::add(const Number                     a,
-                        const VectorSpaceVector<Number> &V,
-                        const Number                     b,
-                        const VectorSpaceVector<Number> &W)
-    {
-      // Check that casting will work.
-      Assert(dynamic_cast<const Vector<Number> *>(&V) != nullptr,
-             ExcVectorTypeNotCompatible());
-      // Check that casting will work.
-      Assert(dynamic_cast<const Vector<Number> *>(&W) != nullptr,
-             ExcVectorTypeNotCompatible());
-
-      // Downcast V. If fails, throws an exception.
-      const Vector<Number> &down_V = dynamic_cast<const Vector<Number> &>(V);
-      // Downcast W. If fails, throws an exception.
-      const Vector<Number> &down_W = dynamic_cast<const Vector<Number> &>(W);
-      Assert(vector->getMap()->isSameAs(*(down_V.trilinos_vector().getMap())),
+    Vector<Number>::add(const Number          a,
+                        const Vector<Number> &V,
+                        const Number          b,
+                        const Vector<Number> &W)
+    {
+      Assert(vector->getMap()->isSameAs(*(V.trilinos_vector().getMap())),
              ExcDifferentParallelPartitioning());
-      Assert(vector->getMap()->isSameAs(*(down_W.trilinos_vector().getMap())),
+      Assert(vector->getMap()->isSameAs(*(W.trilinos_vector().getMap())),
              ExcDifferentParallelPartitioning());
       AssertIsFinite(a);
       AssertIsFinite(b);
 
-      vector->update(
-        a, down_V.trilinos_vector(), b, down_W.trilinos_vector(), 1.);
+      vector->update(a, V.trilinos_vector(), b, W.trilinos_vector(), 1.);
     }
 
 
 
     template <typename Number>
     void
-    Vector<Number>::sadd(const Number                     s,
-                         const Number                     a,
-                         const VectorSpaceVector<Number> &V)
+    Vector<Number>::sadd(const Number          s,
+                         const Number          a,
+                         const Vector<Number> &V)
     {
-      // Check that casting will work.
-      Assert(dynamic_cast<const Vector<Number> *>(&V) != nullptr,
-             ExcVectorTypeNotCompatible());
-
       *this *= s;
-      // Downcast V. It fails, throws an exception.
-      const Vector<Number> &down_V = dynamic_cast<const Vector<Number> &>(V);
-      Vector<Number>        tmp(down_V);
+
+      Vector<Number> tmp(V);
       tmp *= a;
       *this += tmp;
     }
@@ -467,45 +423,28 @@ namespace LinearAlgebra
 
     template <typename Number>
     void
-    Vector<Number>::scale(const VectorSpaceVector<Number> &scaling_factors)
+    Vector<Number>::scale(const Vector<Number> &scaling_factors)
     {
-      // Check that casting will work.
-      Assert(dynamic_cast<const Vector<Number> *>(&scaling_factors) != nullptr,
-             ExcVectorTypeNotCompatible());
-
-      // Downcast scaling_factors. If fails, throws an exception.
-      const Vector<Number> &down_scaling_factors =
-        dynamic_cast<const Vector<Number> &>(scaling_factors);
       Assert(vector->getMap()->isSameAs(
-               *(down_scaling_factors.trilinos_vector().getMap())),
+               *(scaling_factors.trilinos_vector().getMap())),
              ExcDifferentParallelPartitioning());
 
-      vector->elementWiseMultiply(1.,
-                                  *down_scaling_factors.vector,
-                                  *vector,
-                                  0.);
+      vector->elementWiseMultiply(1., *scaling_factors.vector, *vector, 0.);
     }
 
 
 
     template <typename Number>
     void
-    Vector<Number>::equ(const Number a, const VectorSpaceVector<Number> &V)
+    Vector<Number>::equ(const Number a, const Vector<Number> &V)
     {
-      // Check that casting will work.
-      Assert(dynamic_cast<const Vector<Number> *>(&V) != nullptr,
-             ExcVectorTypeNotCompatible());
-
-      // Downcast V. If fails, throws an exception.
-      const Vector<Number> &down_V = dynamic_cast<const Vector<Number> &>(V);
       // If we don't have the same map, copy.
-      if (vector->getMap()->isSameAs(*down_V.trilinos_vector().getMap()) ==
-          false)
+      if (vector->getMap()->isSameAs(*V.trilinos_vector().getMap()) == false)
         this->sadd(0., a, V);
       else
         {
           // Otherwise, just update
-          vector->update(a, down_V.trilinos_vector(), 0.);
+          vector->update(a, V.trilinos_vector(), 0.);
         }
     }
 
@@ -554,7 +493,7 @@ namespace LinearAlgebra
 
 
     template <typename Number>
-    typename LinearAlgebra::VectorSpaceVector<Number>::real_type
+    typename Vector<Number>::real_type
     Vector<Number>::l1_norm() const
     {
       return vector->norm1();
@@ -563,7 +502,7 @@ namespace LinearAlgebra
 
 
     template <typename Number>
-    typename LinearAlgebra::VectorSpaceVector<Number>::real_type
+    typename Vector<Number>::real_type
     Vector<Number>::l2_norm() const
     {
       return vector->norm2();
@@ -572,7 +511,7 @@ namespace LinearAlgebra
 
 
     template <typename Number>
-    typename LinearAlgebra::VectorSpaceVector<Number>::real_type
+    typename Vector<Number>::real_type
     Vector<Number>::linfty_norm() const
     {
       return vector->normInf();
@@ -582,9 +521,9 @@ namespace LinearAlgebra
 
     template <typename Number>
     Number
-    Vector<Number>::add_and_dot(const Number                     a,
-                                const VectorSpaceVector<Number> &V,
-                                const VectorSpaceVector<Number> &W)
+    Vector<Number>::add_and_dot(const Number          a,
+                                const Vector<Number> &V,
+                                const Vector<Number> &W)
     {
       this->add(a, V);
 
@@ -653,6 +592,12 @@ namespace LinearAlgebra
 
 
 
+    template <typename Number>
+    void
+    Vector<Number>::compress(const VectorOperation::values /*operation*/)
+    {}
+
+
     template <typename Number>
     const Tpetra::Vector<Number, int, types::signed_global_dof_index> &
     Vector<Number>::trilinos_vector() const

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.