]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Get rid of VectorSpaceVector in the Cuda Vector class. 15622/head
authorWolfgang Bangerth <bangerth@colostate.edu>
Mon, 3 Jul 2023 16:56:02 +0000 (10:56 -0600)
committerWolfgang Bangerth <bangerth@colostate.edu>
Wed, 5 Jul 2023 16:26:20 +0000 (10:26 -0600)
include/deal.II/lac/cuda_vector.h
source/lac/cuda_vector.cc

index 1eae5e3bd670bab3876ea4acedd483da504b6fb6..8032d9007d8331dfa9538409a9a2a74283cfaaa2 100644 (file)
@@ -43,8 +43,7 @@ namespace LinearAlgebra
   namespace CUDAWrappers
   {
     /**
-     * This class implements a vector using CUDA for use on Nvidia GPUs. This
-     * class is derived from the LinearAlgebra::VectorSpaceVector class.
+     * This class implements a vector using CUDA for use on Nvidia GPUs.
      *
      * @note Only float and double are supported.
      *
@@ -52,12 +51,12 @@ namespace LinearAlgebra
      * @ingroup Vectors
      */
     template <typename Number>
-    class Vector : public VectorSpaceVector<Number>
+    class Vector
     {
     public:
-      using value_type = typename VectorSpaceVector<Number>::value_type;
-      using size_type  = typename VectorSpaceVector<Number>::size_type;
-      using real_type  = typename VectorSpaceVector<Number>::real_type;
+      using value_type = Number;
+      using size_type  = types::global_dof_index;
+      using real_type  = typename numbers::NumberTraits<Number>::real_type;
 
       /**
        * Constructor. Create a vector of dimension zero.
@@ -109,11 +108,8 @@ namespace LinearAlgebra
        * standard containers. Also, there is a global function
        * <tt>swap(u,v)</tt> that simply calls <tt>u.swap(v)</tt>, again in
        * analogy to standard functions.
-       *
-       * This function is virtual in order to allow for derived classes to
-       * handle memory separately.
        */
-      virtual void
+      void
       swap(Vector<Number> &v);
 
       /**
@@ -128,9 +124,8 @@ namespace LinearAlgebra
        * 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;
+      void
+      reinit(const Vector<Number> &V, const bool omit_zeroing_entries = false);
 
       /**
        * Import all the element from the input vector @p V.
@@ -140,21 +135,21 @@ namespace LinearAlgebra
        * for distributed vectors. This is the function that should be used to
        * copy a vector to the GPU.
        */
-      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
+      DEAL_II_DEPRECATED 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);
       }
@@ -163,108 +158,106 @@ 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<Number> &
-      operator=(const Number s) override;
+      Vector<Number> &
+      operator=(const Number s);
 
       /**
        * Multiply the entive vector by a fixed factor.
        */
-      virtual Vector<Number> &
-      operator*=(const Number factor) override;
+      Vector<Number> &
+      operator*=(const Number factor);
 
       /**
        * Divide the entire vector by a fixed factor.
        */
-      virtual Vector<Number> &
-      operator/=(const Number factor) override;
+      Vector<Number> &
+      operator/=(const Number factor);
 
       /**
        * Add the vector @p V to the present one.
        */
-      virtual Vector<Number> &
-      operator+=(const VectorSpaceVector<Number> &V) override;
+      Vector<Number> &
+      operator+=(const Vector<Number> &V);
 
       /**
        * Subtract the vector @p V from the present one.
        */
-      virtual Vector<Number> &
-      operator-=(const VectorSpaceVector<Number> &V) override;
+      Vector<Number> &
+      operator-=(const Vector<Number> &V);
 
       /**
        * Return the scalar product of two vectors.
        */
-      virtual Number
-      operator*(const VectorSpaceVector<Number> &V) const override;
+      Number
+      operator*(const Vector<Number> &V) const;
 
       /**
        * Add @p 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 Vector<Number> &V);
 
       /**
        * Multiple additions 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 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
        * 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 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 all the entries of this vector.
        */
-      virtual value_type
-      mean_value() const override;
+      value_type
+      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 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 real_type
-      l2_norm() const override;
+      real_type
+      l2_norm() const;
 
       /**
        * Return the square of the $l_2$-norm.
@@ -276,8 +269,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
@@ -298,10 +291,10 @@ namespace LinearAlgebra
        * For complex-valued vectors, the scalar product in the second step is
        * 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);
 
       /**
        * Return the pointer to the underlying array. Ownership still resides
@@ -313,30 +306,30 @@ namespace LinearAlgebra
       /**
        * Return the size of the vector.
        */
-      virtual size_type
-      size() const override;
+      size_type
+      size() const;
 
       /**
        * Return an index set that describe which elements of this vector are
        * owned by the current processor, i.e. [0, size).
        */
-      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  = 2,
             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;
 
       /**
        * Attempt to perform an operation between two incompatible vector types.
index 43f62b1de51acf2c2baab100a2ada0d2ae50f6c6..9e3e54ef37688972045cf9775367813163c940ee 100644 (file)
@@ -116,8 +116,8 @@ 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)
     {
       reinit(V.size(), omit_zeroing_entries);
     }
@@ -219,22 +219,16 @@ 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 it fails, it throw an exception.
-      const Vector<Number> &down_V = dynamic_cast<const Vector<Number> &>(V);
-      Assert(down_V.size() == this->size(),
+      Assert(V.size() == this->size(),
              ExcMessage(
                "Cannot add two vectors with different numbers of elements"));
 
       const int n_blocks = 1 + (n_elements - 1) / (chunk_size * block_size);
 
       kernel::vector_bin_op<Number, kernel::Binop_Addition>
-        <<<n_blocks, block_size>>>(val.get(), down_V.val.get(), n_elements);
+        <<<n_blocks, block_size>>>(val.get(), V.val.get(), n_elements);
       AssertCudaKernel();
 
       return *this;
@@ -244,22 +238,16 @@ 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);
-      Assert(down_V.size() == this->size(),
+      Assert(V.size() == this->size(),
              ExcMessage(
                "Cannot add two vectors with different numbers of elements."));
 
       const int n_blocks = 1 + (n_elements - 1) / (chunk_size * block_size);
 
       kernel::vector_bin_op<Number, kernel::Binop_Subtraction>
-        <<<n_blocks, block_size>>>(val.get(), down_V.val.get(), n_elements);
+        <<<n_blocks, block_size>>>(val.get(), V.val.get(), n_elements);
       AssertCudaKernel();
 
       return *this;
@@ -269,15 +257,9 @@ 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(down_V.size() == this->size(),
+      Assert(V.size() == this->size(),
              ExcMessage(
                "Cannot add two vectors with different numbers of elements"));
 
@@ -291,7 +273,7 @@ namespace LinearAlgebra
       kernel::double_vector_reduction<Number, kernel::DotProduct<Number>>
         <<<dim3(n_blocks, 1), dim3(block_size)>>>(result_device,
                                                   val.get(),
-                                                  down_V.val.get(),
+                                                  V.val.get(),
                                                   static_cast<unsigned int>(
                                                     n_elements));
 
@@ -325,23 +307,17 @@ 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)
     {
       AssertIsFinite(a);
 
-      // Check that casting will work.
-      Assert(dynamic_cast<const Vector<Number> *>(&V) != nullptr,
-             ExcVectorTypeNotCompatible());
-
-      // Downcast V. If fails, throw an exception.
-      const Vector<Number> &down_V = dynamic_cast<const Vector<Number> &>(V);
-      Assert(down_V.size() == this->size(),
+      Assert(V.size() == this->size(),
              ExcMessage(
                "Cannot add two vectors with different numbers of elements."));
 
       const int n_blocks = 1 + (n_elements - 1) / (chunk_size * block_size);
       kernel::add_aV<Number><<<dim3(n_blocks, 1), dim3(block_size)>>>(
-        val.get(), a, down_V.val.get(), n_elements);
+        val.get(), a, V.val.get(), n_elements);
       AssertCudaKernel();
     }
 
@@ -349,37 +325,25 @@ namespace LinearAlgebra
 
     template <typename Number>
     void
-    Vector<Number>::add(const Number                     a,
-                        const VectorSpaceVector<Number> &V,
-                        const Number                     b,
-                        const VectorSpaceVector<Number> &W)
+    Vector<Number>::add(const Number          a,
+                        const Vector<Number> &V,
+                        const Number          b,
+                        const Vector<Number> &W)
     {
       AssertIsFinite(a);
       AssertIsFinite(b);
 
-      // Check that casting will work.
-      Assert(dynamic_cast<const Vector<Number> *>(&V) != nullptr,
-             ExcVectorTypeNotCompatible());
-
-      // Downcast V. If fails, throw an exception.
-      const Vector<Number> &down_V = dynamic_cast<const Vector<Number> &>(V);
-      Assert(down_V.size() == this->size(),
+      Assert(V.size() == this->size(),
              ExcMessage(
                "Cannot add two vectors with different numbers of elements."));
 
-      // Check that casting will work.
-      Assert(dynamic_cast<const Vector<Number> *>(&W) != nullptr,
-             ExcVectorTypeNotCompatible());
-
-      // Downcast V. If fails, throw an exception.
-      const Vector<Number> &down_W = dynamic_cast<const Vector<Number> &>(W);
-      Assert(down_W.size() == this->size(),
+      Assert(W.size() == this->size(),
              ExcMessage(
                "Cannot add two vectors with different numbers of elements."));
 
       const int n_blocks = 1 + (n_elements - 1) / (chunk_size * block_size);
       kernel::add_aVbW<Number><<<dim3(n_blocks, 1), dim3(block_size)>>>(
-        val.get(), a, down_V.val.get(), b, down_W.val.get(), n_elements);
+        val.get(), a, V.val.get(), b, W.val.get(), n_elements);
       AssertCudaKernel();
     }
 
@@ -387,26 +351,20 @@ namespace LinearAlgebra
 
     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)
     {
       AssertIsFinite(s);
       AssertIsFinite(a);
 
-      // Check that casting will work.
-      Assert(dynamic_cast<const Vector<Number> *>(&V) != nullptr,
-             ExcVectorTypeNotCompatible());
-
-      // Downcast V. If fails, throw an exception.
-      const Vector<Number> &down_V = dynamic_cast<const Vector<Number> &>(V);
-      Assert(down_V.size() == this->size(),
+      Assert(V.size() == this->size(),
              ExcMessage(
                "Cannot add two vectors with different numbers of elements."));
 
       const int n_blocks = 1 + (n_elements - 1) / (chunk_size * block_size);
       kernel::sadd<Number><<<dim3(n_blocks, 1), dim3(block_size)>>>(
-        s, val.get(), a, down_V.val.get(), n_elements);
+        s, val.get(), a, V.val.get(), n_elements);
       AssertCudaKernel();
     }
 
@@ -414,22 +372,17 @@ 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 V. If fails, throw an exception.
-      const Vector<Number> &down_scaling_factors =
-        dynamic_cast<const Vector<Number> &>(scaling_factors);
-      Assert(down_scaling_factors.size() == this->size(),
+      Assert(scaling_factors.size() == this->size(),
              ExcMessage(
                "Cannot scale two vectors with different numbers of elements."));
 
       const int n_blocks = 1 + (n_elements - 1) / (chunk_size * block_size);
-      kernel::scale<Number><<<dim3(n_blocks, 1), dim3(block_size)>>>(
-        val.get(), down_scaling_factors.val.get(), n_elements);
+      kernel::scale<Number>
+        <<<dim3(n_blocks, 1), dim3(block_size)>>>(val.get(),
+                                                  scaling_factors.val.get(),
+                                                  n_elements);
       AssertCudaKernel();
     }
 
@@ -437,24 +390,20 @@ namespace LinearAlgebra
 
     template <typename Number>
     void
-    Vector<Number>::equ(const Number a, const VectorSpaceVector<Number> &V)
+    Vector<Number>::equ(const Number a, const Vector<Number> &V)
     {
       AssertIsFinite(a);
 
-      // Check that casting will work.
-      Assert(dynamic_cast<const Vector<Number> *>(&V) != nullptr,
-             ExcVectorTypeNotCompatible());
-
-      // Downcast V. If fails, throw an exception.
-      const Vector<Number> &down_V = dynamic_cast<const Vector<Number> &>(V);
       Assert(
-        down_V.size() == this->size(),
+        V.size() == this->size(),
         ExcMessage(
           "Cannot assign two vectors with different numbers of elements."));
 
       const int n_blocks = 1 + (n_elements - 1) / (chunk_size * block_size);
-      kernel::equ<Number><<<dim3(n_blocks, 1), dim3(block_size)>>>(
-        val.get(), a, down_V.val.get(), n_elements);
+      kernel::equ<Number><<<dim3(n_blocks, 1), dim3(block_size)>>>(val.get(),
+                                                                   a,
+                                                                   V.val.get(),
+                                                                   n_elements);
       AssertCudaKernel();
     }
 
@@ -580,24 +529,15 @@ 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)
     {
       AssertIsFinite(a);
 
-      // Check that casting will work
-      Assert(dynamic_cast<const Vector<Number> *>(&V) != nullptr,
-             ExcVectorTypeNotCompatible());
-      Assert(dynamic_cast<const Vector<Number> *>(&W) != nullptr,
-             ExcVectorTypeNotCompatible());
-
-      // Downcast V and W. If it fails, throw an exception.
-      const Vector<Number> &down_V = dynamic_cast<const Vector<Number> &>(V);
-      Assert(down_V.size() == this->size(),
+      Assert(V.size() == this->size(),
              ExcMessage("Vector V has the wrong size."));
-      const Vector<Number> &down_W = dynamic_cast<const Vector<Number> &>(W);
-      Assert(down_W.size() == this->size(),
+      Assert(W.size() == this->size(),
              ExcMessage("Vector W has the wrong size."));
 
       Number *    result_device;
@@ -607,13 +547,8 @@ namespace LinearAlgebra
       AssertCuda(error_code);
 
       const int n_blocks = 1 + (n_elements - 1) / (chunk_size * block_size);
-      kernel::add_and_dot<Number>
-        <<<dim3(n_blocks, 1), dim3(block_size)>>>(result_device,
-                                                  val.get(),
-                                                  down_V.val.get(),
-                                                  down_W.val.get(),
-                                                  a,
-                                                  n_elements);
+      kernel::add_and_dot<Number><<<dim3(n_blocks, 1), dim3(block_size)>>>(
+        result_device, val.get(), V.val.get(), W.val.get(), a, n_elements);
 
       Number result;
       error_code = cudaMemcpy(&result,

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.