Fix alignment errors from cuSPARSE

cmake/ExternalHYPRE.cmake

@@ -141,7 +141,7 @@ if(PALACE_WITH_CUDA)
     "--with-cuda-home=${CUDA_DIR}"
     "--enable-curand"
     "--enable-cusparse"
-    "--enable-device-memory-pool"
+    # "--enable-device-memory-pool"
   )
   if(NOT "${CMAKE_CUDA_ARCHITECTURES}" STREQUAL "")
     list(GET CMAKE_CUDA_ARCHITECTURES 0 HYPRE_CUDA_ARCH)

cmake/ExternalMFEM.cmake

@@ -364,6 +364,7 @@ set(MFEM_PATCH_FILES
   "${CMAKE_SOURCE_DIR}/extern/patch/mfem/patch_par_tet_mesh_fix.diff"
   "${CMAKE_SOURCE_DIR}/extern/patch/mfem/patch_mesh_prism_vtu_fix.diff"
   "${CMAKE_SOURCE_DIR}/extern/patch/mfem/patch_workspace_vectors.diff"
+  "${CMAKE_SOURCE_DIR}/extern/patch/mfem/patch_workspace_vectors_2.diff"
   "${CMAKE_SOURCE_DIR}/extern/patch/mfem/patch_hypre_runtime_compute_policy.diff"
 )
 

extern/patch/mfem/patch_workspace_vectors_2.diff

@@ -0,0 +1,76 @@
+diff --git a/general/workspace.cpp b/general/workspace.cpp
+index ad8210a61..69622ec70 100644
+--- a/general/workspace.cpp
++++ b/general/workspace.cpp
+@@ -15,11 +15,11 @@ namespace mfem
+ {
+
+ WorkspaceVector::WorkspaceVector(
+-   internal::WorkspaceChunk &chunk_, int offset_, int n)
++   internal::WorkspaceChunk &chunk_, int offset_, int n, int padding)
+    : Vector(chunk_.GetData(), chunk_.GetOffset(), n),
+      chunk(chunk_),
+      offset(offset_),
+-     original_size(n)
++     original_size(n + padding)
+ {
+    UseDevice(true);
+ }
+@@ -47,9 +47,12 @@ WorkspaceChunk::WorkspaceChunk(int capacity)
+
+ WorkspaceVector WorkspaceChunk::NewVector(int n)
+ {
+-   MFEM_ASSERT(HasCapacityFor(n), "Requested vector is too large.");
+-   WorkspaceVector vector(*this, offset, n);
+-   offset += n;
++   constexpr int alignment = 16;
++   const int s = (n * sizeof(double)) % alignment;
++   const int padding = (s > 0) ? (alignment - s) / sizeof(double) : 0;
++   MFEM_ASSERT(HasCapacityFor(n + padding), "Requested vector is too large.");
++   WorkspaceVector vector(*this, offset, n, padding);
++   offset += n + padding;
+    vector_count += 1;
+    return vector;
+ }
+diff --git a/general/workspace.hpp b/general/workspace.hpp
+index b5c0ebae8..e41e71b11 100644
+--- a/general/workspace.hpp
++++ b/general/workspace.hpp
+@@ -50,7 +50,8 @@ class WorkspaceVector : public Vector
+    bool moved_from = false;
+
+    /// Private constructor, create with Workspace::NewVector() instead.
+-   WorkspaceVector(internal::WorkspaceChunk &chunk_, int offset_, int n);
++   WorkspaceVector(internal::WorkspaceChunk &chunk_, int offset_, int n,
++                   int padding = 0);
+
+ public:
+    /// @brief Move constructor. The moved-from WorkspaceVector has @a
+@@ -104,6 +105,9 @@ public:
+    /// Create a WorkspaceChunk with the given @a capacity.
+    WorkspaceChunk(int capacity);
+
++   /// Return the data offset.
++   int GetOffset() const { return offset; }
++
+    /// @brief Return the available capacity (i.e. the largest vector that will
+    /// fit in this chunk).
+    int GetAvailableCapacity() const { return data.Size() - offset; }
+@@ -115,15 +119,12 @@ public:
+    /// "original capacity" remains unchained.
+    int GetOriginalCapacity() const { return original_capacity; }
+
+-   /// Return the data offset.
+-   int GetOffset() const { return offset; }
++   /// Returns true if this chunk can fit a new vector of size @a n.
++   bool HasCapacityFor(int n) const { return n <= GetAvailableCapacity(); }
+
+    /// Sets whether the chunk is in the front of the list
+    void SetFront(bool front_) { front = front_; }
+
+-   /// Returns true if this chunk can fit a new vector of size @a n.
+-   bool HasCapacityFor(int n) const { return n <= GetAvailableCapacity(); }
+-
+    /// Returns true if this chunk is empty.
+    bool IsEmpty() const { return vector_count == 0; }
+

palace/linalg/vector.cpp

@@ -48,6 +48,7 @@ void ComplexVector::UseDevice(bool use_dev)
 
 void ComplexVector::SetSize(int size)
 {
+  MFEM_VERIFY(xr.OwnsData() && xi.OwnsData(), "Alias vectors should not be resized!");
   xr.SetSize(size);
   xi.SetSize(size);
 }

palace/linalg/vector.hpp

@@ -21,7 +21,7 @@ using Vector = mfem::Vector;
 // A complex-valued vector represented as two real vectors, one for each component.
 class ComplexVector
 {
-private:
+protected:
   Vector xr, xi;
 
 public:

palace/utils/workspace.hpp

@@ -18,12 +18,14 @@ namespace palace::workspace
 class ComplexWorkspaceVector : public ComplexVector
 {
 private:
-  mfem::WorkspaceVector data;
+  mfem::WorkspaceVector data_xr, data_xi;
 
 public:
-  ComplexWorkspaceVector(int size) : data(mfem::Workspace::NewVector(2 * size))
+  ComplexWorkspaceVector(int size)
+    : data_xr(mfem::Workspace::NewVector(size)), data_xi(mfem::Workspace::NewVector(size))
   {
-    MakeRef(data, 0, size);
+    xr.MakeRef(data_xr, 0, size);
+    xi.MakeRef(data_xi, 0, size);
   }
 
   // No copy constructor.