Make GpuVector support default constructing vectors and resizing

opm/simulators/linalg/gpuistl/GpuVector.cpp

@@ -67,13 +67,15 @@ template <class T>
 GpuVector<T>&
 GpuVector<T>::operator=(const GpuVector<T>& other)
 {
-    assertHasElements();
-    assertSameSize(other);
-
-    OPM_GPU_SAFE_CALL(cudaMemcpy(m_dataOnDevice,
-                                  other.m_dataOnDevice,
-                                  detail::to_size_t(m_numberOfElements) * sizeof(T),
-                                  cudaMemcpyDeviceToDevice));
+    // Only copy data if both vectors have elements and same size
+    if (m_numberOfElements > 0 && other.m_numberOfElements > 0) {
+        assertSameSize(other);
+        OPM_GPU_SAFE_CALL(cudaMemcpy(m_dataOnDevice,
+                                      other.m_dataOnDevice,
+                                      detail::to_size_t(m_numberOfElements) * sizeof(T),
+                                      cudaMemcpyDeviceToDevice));
+    }
+    // If both are zero-sized, assignment is trivial (do nothing)
     return *this;
 }
 
@@ -81,12 +83,15 @@ template <class T>
 GpuVector<T>::GpuVector(const GpuVector<T>& other)
     : GpuVector(other.m_numberOfElements)
 {
-    assertHasElements();
-    assertSameSize(other);
-    OPM_GPU_SAFE_CALL(cudaMemcpy(m_dataOnDevice,
-                                  other.m_dataOnDevice,
-                                  detail::to_size_t(m_numberOfElements) * sizeof(T),
-                                  cudaMemcpyDeviceToDevice));
+    // Only copy data if both vectors have elements and same size
+    if (m_numberOfElements > 0) {
+        assertSameSize(other);
+        OPM_GPU_SAFE_CALL(cudaMemcpy(m_dataOnDevice,
+                                      other.m_dataOnDevice,
+                                      detail::to_size_t(m_numberOfElements) * sizeof(T),
+                                      cudaMemcpyDeviceToDevice));
+    }
+    // If other is zero-sized, assignment is trivial (do nothing)
 }
 
 template <class T>
@@ -114,6 +119,47 @@ GpuVector<T>::dim() const
     return detail::to_size_t(m_numberOfElements);
 }
 
+template <typename T>
+void
+GpuVector<T>::resize(size_t new_size)
+{
+    const int new_elements = detail::to_int(new_size);
+
+    if (new_elements == m_numberOfElements) {
+        return;
+    }
+
+    if (new_elements == 0) {
+        // Free existing memory and set to empty state
+        if (m_dataOnDevice != nullptr) {
+            OPM_GPU_WARN_IF_ERROR(cudaFree(m_dataOnDevice));
+            m_dataOnDevice = nullptr;
+        }
+        m_numberOfElements = 0;
+        return;
+    }
+
+    // Allocate new memory
+    T* new_data = nullptr;
+    OPM_GPU_SAFE_CALL(cudaMalloc(&new_data, sizeof(T) * new_size));
+
+    if (m_dataOnDevice != nullptr && m_numberOfElements > 0) {
+        // Copy existing data (up to the minimum of old and new size)
+        const size_t copy_elements = std::min(detail::to_size_t(m_numberOfElements), new_size);
+        if (copy_elements > 0) {
+            OPM_GPU_SAFE_CALL(cudaMemcpy(new_data, m_dataOnDevice,
+                                          sizeof(T) * copy_elements,
+                                          cudaMemcpyDeviceToDevice));
+        }
+
+        // Free old memory
+        OPM_GPU_WARN_IF_ERROR(cudaFree(m_dataOnDevice));
+    }
+
+    m_dataOnDevice = new_data;
+    m_numberOfElements = new_elements;
+}
+
 template <typename T>
 std::vector<T>
 GpuVector<T>::asStdVector() const

opm/simulators/linalg/gpuistl/GpuVector.hpp

@@ -130,6 +130,11 @@ public:
      */
     GpuVector& operator=(T scalar);
 
+    /**
+     * @brief GpuVector default constructor creates a zero-sized vector (no GPU memory allocated)
+     */
+    GpuVector() : m_dataOnDevice(nullptr), m_numberOfElements(0), m_cuBlasHandle(detail::CuBlasHandle::getInstance()) {}
+
     /**
      * @brief GpuVector allocates new GPU memory of size numberOfElements * sizeof(T)
      *
@@ -449,6 +454,15 @@ public:
      */
     size_type dim() const;
 
+    /**
+     * @brief resize changes the size of the vector, preserving existing data if new size is larger
+     * @param new_size the new number of elements
+     * @note If new_size is larger, existing data is preserved and new elements are uninitialized
+     * @note If new_size is smaller, data is truncated
+     * @note If new_size equals current size, no operation is performed
+     * @note For now new_size needs to be within the limits of int due to restrictions in cublas
+     */
+    void resize(size_t new_size);
 
     /**
      * @brief creates an std::vector of the same size and copies the GPU data to this std::vector
@@ -457,8 +471,8 @@ public:
     std::vector<T> asStdVector() const;
 
     /**
-     * @brief creates an std::vector of the same size and copies the GPU data to this std::vector
-     * @return an std::vector containing the elements copied from the GPU.
+     * @brief creates a Dune::BlockVector of the same size and copies the GPU data to this Dune::BlockVector
+     * @return a Dune::BlockVector containing the elements copied from the GPU.
      */
     template <int blockSize>
     Dune::BlockVector<Dune::FieldVector<T, blockSize>> asDuneBlockVector() const
@@ -506,7 +520,7 @@ private:
 
     // Note that we store this as int to make sure we are always cublas compatible.
     // This gives the added benefit that a size_t to int conversion error occurs during construction.
-    const int m_numberOfElements;
+    int m_numberOfElements;
     detail::CuBlasHandle& m_cuBlasHandle;
 
     void assertSameSize(const GpuVector<T>& other) const;