HypreParMatrix and BlockOperator cleanup

src/redecomp/partition/PartitionElements.cpp

@@ -49,6 +49,9 @@ std::vector<CoordList<NDIMS>> PartitionElements<NDIMS>::EntityCoordinates(
       //    b) Interpolate to Quadrature points (Q-vector)
       qi->Values( e_vec, mesh_centroids );
 
+      // NOTE: mesh_centroids will be fed into the RCB PartitionMethod and CoordList which are host only. make sure host
+      // is updated.
+      mesh_centroids.HostRead();
       return mesh_centroids;
     };
 

src/redecomp/transfer/MatrixTransfer.cpp

@@ -10,6 +10,7 @@
 #include "mfem/general/forall.hpp"
 
 #include "redecomp/RedecompMesh.hpp"
+#include "shared/math/ParSparseMat.hpp"
 
 namespace redecomp {
 
@@ -37,13 +38,14 @@ MatrixTransfer::MatrixTransfer( const mfem::ParFiniteElementSpace& parent_test_f
   test_r2p_elem_rank_ = buildRedecomp2ParentElemRank( *test_redecomp, true );
 }
 
-std::unique_ptr<mfem::HypreParMatrix> MatrixTransfer::TransferToParallel(
-    const axom::Array<int>& test_elem_idx, const axom::Array<int>& trial_elem_idx,
-    const axom::Array<mfem::DenseMatrix>& src_elem_mat, bool parallel_assemble ) const
+shared::ParSparseMat MatrixTransfer::TransferToParallel( const axom::Array<int>& test_elem_idx,
+                                                         const axom::Array<int>& trial_elem_idx,
+                                                         const axom::Array<mfem::DenseMatrix>& src_elem_mat,
+                                                         bool parallel_assemble ) const
 {
   auto J_sparse = TransferToParallelSparse( test_elem_idx, trial_elem_idx, src_elem_mat );
   J_sparse.Finalize();
-  return ConvertToHypreParMatrix( J_sparse, parallel_assemble );
+  return ConvertToParSparseMat( std::move( J_sparse ), parallel_assemble );
 }
 
 mfem::SparseMatrix MatrixTransfer::TransferToParallelSparse( const axom::Array<int>& test_elem_idx,
@@ -142,8 +144,7 @@ mfem::SparseMatrix MatrixTransfer::TransferToParallelSparse( const axom::Array<i
   return parentJ;
 }
 
-std::unique_ptr<mfem::HypreParMatrix> MatrixTransfer::ConvertToHypreParMatrix( mfem::SparseMatrix& sparse,
-                                                                               bool parallel_assemble ) const
+shared::ParSparseMat MatrixTransfer::ConvertToParSparseMat( mfem::SparseMatrix&& sparse, bool parallel_assemble ) const
 {
   SLIC_ERROR_IF( sparse.Height() != parent_test_fes_.GetVSize(),
                  "Height of sparse must match number of test ParFiniteElementSpace L-dofs." );
@@ -161,15 +162,17 @@ std::unique_ptr<mfem::HypreParMatrix> MatrixTransfer::ConvertToHypreParMatrix( m
   // update the host pointer
   J_bigint.HostRead();
 
-  auto J_full = std::make_unique<mfem::HypreParMatrix>(
-      getMPIUtility().MPIComm(), parent_test_fes_.GetVSize(), parent_test_fes_.GlobalVSize(),
-      parent_trial_fes_.GlobalVSize(), sparse.GetI(), J_bigint.GetData(), sparse.GetData(),
-      parent_test_fes_.GetDofOffsets(), parent_trial_fes_.GetDofOffsets() );
+  shared::ParSparseMat J_full( getMPIUtility().MPIComm(), parent_test_fes_.GetVSize(), parent_test_fes_.GlobalVSize(),
+                               parent_trial_fes_.GlobalVSize(), sparse.GetI(), J_bigint.GetData(), sparse.GetData(),
+                               parent_test_fes_.GetDofOffsets(), parent_trial_fes_.GetDofOffsets() );
   if ( !parallel_assemble ) {
     return J_full;
   } else {
-    auto J_true = std::unique_ptr<mfem::HypreParMatrix>(
-        mfem::RAP( parent_test_fes_.Dof_TrueDof_Matrix(), J_full.get(), parent_trial_fes_.Dof_TrueDof_Matrix() ) );
+    auto P_test = parent_test_fes_.Dof_TrueDof_Matrix();
+    P_test->HostRead();
+    auto P_trial = parent_trial_fes_.Dof_TrueDof_Matrix();
+    P_trial->HostRead();
+    auto J_true = shared::ParSparseMat::RAP( P_test, J_full, P_trial );
     return J_true;
   }
 }

src/redecomp/transfer/MatrixTransfer.hpp

@@ -8,7 +8,11 @@
 
 #include "mfem.hpp"
 
-#include "redecomp/common/TypeDefs.hpp"
+#include "axom/core.hpp"
+
+#include "shared/math/ParSparseMat.hpp"
+
+#include "redecomp/utils/MPIArray.hpp"
 
 namespace redecomp {
 
@@ -26,7 +30,7 @@ class RedecompMesh;
  * two-stage transfer process is available. First, TransferToParallel() creates
  * an un-finalized mfem::SparseMatrix with ldofs on the rows and global ldofs on
  * the columns. This mfem::SparseMatrix is designed to be passed to a
- * mfem::HypreParMatrix constructor (done through the ConvertToHypreParMatrix()
+ * mfem::HypreParMatrix constructor (done through the ConvertToParSparseMat()
  * method).  The two-stage process allows easier manipulation of matrix
  * contributions before the matrix is finalized.  Both square and rectangular
  * matrices are supported, necessitating test and trial finite element spaces in
@@ -61,10 +65,10 @@ class MatrixTransfer {
    * HypreParMatrix returned to transform it to the t-dofs if parallel_assemble
    * is false.
    */
-  std::unique_ptr<mfem::HypreParMatrix> TransferToParallel( const axom::Array<int>& test_elem_idx,
-                                                            const axom::Array<int>& trial_elem_idx,
-                                                            const axom::Array<mfem::DenseMatrix>& src_elem_mat,
-                                                            bool parallel_assemble = true ) const;
+  shared::ParSparseMat TransferToParallel( const axom::Array<int>& test_elem_idx,
+                                           const axom::Array<int>& trial_elem_idx,
+                                           const axom::Array<mfem::DenseMatrix>& src_elem_mat,
+                                           bool parallel_assemble = true ) const;
 
   /**
    * @brief Transfers element RedecompMesh matrices to parent mfem::ParMesh
@@ -94,8 +98,7 @@ class MatrixTransfer {
    * HypreParMatrix returned to transform it to the t-dofs if parallel_assemble
    * is false.
    */
-  std::unique_ptr<mfem::HypreParMatrix> ConvertToHypreParMatrix( mfem::SparseMatrix& sparse,
-                                                                 bool parallel_assemble = true ) const;
+  shared::ParSparseMat ConvertToParSparseMat( mfem::SparseMatrix&& sparse, bool parallel_assemble = true ) const;
 
  private:
   /**

src/shared/CMakeLists.txt

@@ -6,18 +6,24 @@
 ## list of headers
 set(shared_headers
 
+    common/BasicTypes.hpp
+    common/ExecModel.hpp
     infrastructure/Profiling.hpp
+    math/ParSparseMat.hpp
+    math/ParVector.hpp
     mesh/MeshBuilder.hpp
     )
 
 ## list of sources
 set(shared_sources
 
+    math/ParSparseMat.cpp
+    math/ParVector.cpp
     mesh/MeshBuilder.cpp
     )
 
 ## setup the dependency list for tribol shared library
-set(shared_depends mfem)
+set(shared_depends mfem axom::core)
 blt_list_append(TO shared_depends ELEMENTS blt::mpi IF TRIBOL_USE_MPI )
 blt_list_append(TO shared_depends ELEMENTS caliper IF TRIBOL_USE_CALIPER )
 message(STATUS "Tribol shared library dependencies: ${shared_depends}")
@@ -31,7 +37,7 @@ blt_add_library(
     FOLDER shared
     )
 
-# Add separable compilation flag
+# Add debug flag
 if(ENABLE_HIP)
     target_compile_options(tribol_shared PRIVATE
         $<$<CONFIG:Debug>:

src/shared/common/BasicTypes.hpp

@@ -0,0 +1,108 @@
+// Copyright (c) 2017-2025, Lawrence Livermore National Security, LLC and
+// other Tribol Project Developers. See the top-level LICENSE file for details.
+//
+// SPDX-License-Identifier: (MIT)
+
+#ifndef SRC_SHARED_COMMON_BASICTYPES_HPP_
+#define SRC_SHARED_COMMON_BASICTYPES_HPP_
+
+// Tribol config include
+#include "shared/config.hpp"
+
+// C includes
+#include <cstddef>
+
+// C++ includes
+#include <type_traits>
+
+// MPI includes
+#ifdef TRIBOL_USE_MPI
+#include <mpi.h>
+#endif
+
+// Axom includes
+#include "axom/core/Types.hpp"
+
+// MFEM includes
+#include "mfem.hpp"
+
+namespace shared {
+
+#ifdef TRIBOL_USE_MPI
+
+using CommT = MPI_Comm;
+#define TRIBOL_COMM_WORLD MPI_COMM_WORLD
+#define TRIBOL_COMM_NULL MPI_COMM_NULL
+
+#else
+
+using CommT = int;
+#define TRIBOL_COMM_WORLD 0
+#define TRIBOL_COMM_NULL -1
+
+#endif
+
+// match index type used in axom (since data is held in axom data structures)
+using IndexT = axom::IndexType;
+
+// size type matching size of addressable memory
+using SizeT = size_t;
+
+#ifdef TRIBOL_USE_SINGLE_PRECISION
+
+#error "Tribol does not support single precision."
+using RealT = float;
+
+#else
+
+using RealT = double;
+
+#endif
+
+// mfem's real_t should match ours
+static_assert( std::is_same_v<RealT, mfem::real_t>, "tribol::RealT and mfem::real_t are required to match" );
+
+#define TRIBOL_UNUSED_VAR AXOM_UNUSED_VAR
+#define TRIBOL_UNUSED_PARAM AXOM_UNUSED_PARAM
+
+// Execution space specifiers
+#if defined( TRIBOL_USE_CUDA ) || defined( TRIBOL_USE_HIP )
+#ifndef __device__
+#error "TRIBOL_USE_CUDA or TRIBOL_USE_HIP but __device__ is undefined.  Check include files"
+#endif
+#define TRIBOL_DEVICE __device__
+#define TRIBOL_HOST_DEVICE __host__ __device__
+#else
+#define TRIBOL_DEVICE
+#define TRIBOL_HOST_DEVICE
+#endif
+
+// Execution space identifier for defaulted constructors and destructors
+#ifdef TRIBOL_USE_HIP
+#define TRIBOL_DEFAULT_DEVICE __device__
+#define TRIBOL_DEFAULT_HOST_DEVICE __host__ __device__
+#else
+#define TRIBOL_DEFAULT_DEVICE
+#define TRIBOL_DEFAULT_HOST_DEVICE
+#endif
+
+// Defined when Tribol doesn't have a device available
+#if !( defined( TRIBOL_USE_CUDA ) || defined( TRIBOL_USE_HIP ) )
+#define TRIBOL_USE_HOST
+#endif
+
+// Define variable when in device code
+#if defined( __CUDA_ARCH__ ) || defined( __HIP_DEVICE_COMPILE__ )
+#define TRIBOL_DEVICE_CODE
+#endif
+
+// Ignore host code in __host__ __device__ code warning on NVCC
+#ifdef TRIBOL_USE_CUDA
+#define TRIBOL_NVCC_EXEC_CHECK_DISABLE #pragma nv_exec_check_disable
+#else
+#define TRIBOL_NVCC_EXEC_CHECK_DISABLE
+#endif
+
+}  // namespace shared
+
+#endif /* SRC_SHARED_COMMON_BASICTYPES_HPP_ */