From 10f7a8553285f8d70c56ff807f710bce6690b9a7 Mon Sep 17 00:00:00 2001 From: jakobtorben Date: Wed, 6 Nov 2024 16:00:19 +0100 Subject: [PATCH 01/13] Add Hypre preconditioner --- CMakeLists.txt | 15 ++ CMakeLists_files.cmake | 6 + opm-simulators-prereqs.cmake | 3 + opm/simulators/linalg/HyprePreconditioner.hpp | 208 ++++++++++++++++++ .../linalg/PreconditionerFactory_impl.hpp | 11 + 5 files changed, 243 insertions(+) create mode 100644 opm/simulators/linalg/HyprePreconditioner.hpp diff --git a/CMakeLists.txt b/CMakeLists.txt index 90af1219b..cc789f43f 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -35,6 +35,7 @@ option(USE_DAMARIS_LIB "Use the Damaris library for asynchronous I/O?" OFF) option(USE_GPU_BRIDGE "Enable the GPU bridge (GPU/AMGCL solvers)" ON) option(USE_TRACY_PROFILER "Enable tracy profiling" OFF) option(CONVERT_CUDA_TO_HIP "Convert CUDA code to HIP (to run on AMD cards)" OFF) +option(USE_HYPRE "Use the Hypre library for linear solvers?" OFF) set(OPM_COMPILE_COMPONENTS "2;3;4;5;6;7" CACHE STRING "The components to compile support for") # Ensure OPM_COMPILE_COMPONENTS is a list of at least one element @@ -715,3 +716,17 @@ endif() install(DIRECTORY doc/man1 DESTINATION ${CMAKE_INSTALL_MANDIR} FILES_MATCHING PATTERN "*.1") + +if(USE_HYPRE) + find_package(HYPRE) + if(HYPRE_FOUND) + set(HAVE_HYPRE 1) + else() + message(WARNING "Hypre requested but not found. Continuing without Hypre support.") + set(USE_HYPRE OFF) + endif() +endif() + +if(HYPRE_FOUND) + target_link_libraries(opmsimulators PUBLIC HYPRE::HYPRE) +endif() diff --git a/CMakeLists_files.cmake b/CMakeLists_files.cmake index 353567f16..e3ece638b 100644 --- a/CMakeLists_files.cmake +++ b/CMakeLists_files.cmake @@ -1156,3 +1156,9 @@ if(dune-alugrid_FOUND) examples/fracture_discretefracture.cpp ) endif() + +if(HYPRE_FOUND) + list(APPEND PUBLIC_HEADER_FILES + opm/simulators/linalg/HyprePreconditioner.hpp + ) +endif() diff --git a/opm-simulators-prereqs.cmake b/opm-simulators-prereqs.cmake index f2c6fb9be..2af64605f 100644 --- a/opm-simulators-prereqs.cmake +++ b/opm-simulators-prereqs.cmake @@ -23,6 +23,7 @@ set (opm-simulators_CONFIG_VAR HAVE_SUITESPARSE_UMFPACK HAVE_DAMARIS HAVE_HDF5 + HAVE_HYPRE USE_HIP USE_TRACY FLOW_INSTANTIATE_FLOAT @@ -56,6 +57,8 @@ set (opm-simulators_DEPS # packages from ROCm framework "rocblas" "rocsparse" + # HYPRE solver library + "HYPRE" # OPM dependency "opm-common REQUIRED" "opm-grid REQUIRED" diff --git a/opm/simulators/linalg/HyprePreconditioner.hpp b/opm/simulators/linalg/HyprePreconditioner.hpp new file mode 100644 index 000000000..a2e4d50b3 --- /dev/null +++ b/opm/simulators/linalg/HyprePreconditioner.hpp @@ -0,0 +1,208 @@ +/* + Copyright 2024 SINTEF AS + This file is part of the Open Porous Media project (OPM). + OPM is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + OPM is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + You should have received a copy of the GNU General Public License + along with OPM. If not, see . +*/ + +#ifndef OPM_HYPRE_PRECONDITIONER_HEADER_INCLUDED +#define OPM_HYPRE_PRECONDITIONER_HEADER_INCLUDED + +#include +#include +#include + +#include +#include + +#include +#include +#include + +#include +#include + +namespace Dune { + +/// Wrapper for Hypre's BoomerAMG preconditioner +template +class HyprePreconditioner : public PreconditionerWithUpdate { +public: + //! \brief The matrix type the preconditioner is for + using matrix_type = M; + //! \brief The domain type of the preconditioner + using domain_type = X; + //! \brief The range type of the preconditioner + using range_type = Y; + //! \brief The field type of the preconditioner + using field_type = typename X::field_type; + + // Constructor + HyprePreconditioner (const M& A) + : A_(A) + { + OPM_TIMEBLOCK(prec_construct); + + // Initialize Hypre + HYPRE_Init(); + + // Create the solver (BoomerAMG) + HYPRE_BoomerAMGCreate(&solver_); + + // Set some default parameters + HYPRE_BoomerAMGSetPrintLevel(solver_, 0); // Reduce output + HYPRE_BoomerAMGSetOldDefault(solver_); // Falgout coarsening with modified classical interpolation + HYPRE_BoomerAMGSetRelaxType(solver_, 3); // G-S/Jacobi hybrid relaxation + HYPRE_BoomerAMGSetRelaxOrder(solver_, 1); // Uses C/F relaxation + HYPRE_BoomerAMGSetNumSweeps(solver_, 1); // Sweeps on each level + HYPRE_BoomerAMGSetMaxLevels(solver_, 20); // Maximum number of levels + HYPRE_BoomerAMGSetTol(solver_, 1e-7); // Convergence tolerance + + + // Create the vectors + const int N = A_.N(); + HYPRE_IJVectorCreate(MPI_COMM_WORLD, 0, N-1, &x_hypre_); + HYPRE_IJVectorCreate(MPI_COMM_WORLD, 0, N-1, &b_hypre_); + HYPRE_IJVectorSetObjectType(x_hypre_, HYPRE_PARCSR); + HYPRE_IJVectorSetObjectType(b_hypre_, HYPRE_PARCSR); + HYPRE_IJVectorInitialize(x_hypre_); + HYPRE_IJVectorInitialize(b_hypre_); + + update(); + } + + // Destructor + ~HyprePreconditioner() { + if (solver_) { + HYPRE_BoomerAMGDestroy(solver_); + } + if (parcsr_A_) { + HYPRE_IJMatrixDestroy(A_hypre_); + } + HYPRE_Finalize(); + } + + void update() override { + OPM_TIMEBLOCK(prec_update); + copyMatrixToHypre(); + } + + void pre(X& x, Y& b) override { + DUNE_UNUSED_PARAMETER(x); + DUNE_UNUSED_PARAMETER(b); + } + + void apply(X& v, const Y& d) override { + OPM_TIMEBLOCK(prec_apply); + + // Copy vectors to Hypre format + copyVectorsToHypre(v, d); + + // Apply the preconditioner (one AMG V-cycle) + HYPRE_BoomerAMGSolve(solver_, parcsr_A_, par_b_, par_x_); + + // Copy result back + copyVectorFromHypre(v); + } + + void post(X& x) override { + DUNE_UNUSED_PARAMETER(x); + } + + SolverCategory::Category category() const override { + return SolverCategory::sequential; + } + + bool hasPerfectUpdate() const override + { + // The Hypre preconditioner can depend on the values of the matrix, so it must be recreated + return false; + } + +private: + void copyMatrixToHypre() { + const int N = A_.N(); + HYPRE_IJMatrixCreate(MPI_COMM_WORLD, 0, N-1, 0, N-1, &A_hypre_); + HYPRE_IJMatrixSetObjectType(A_hypre_, HYPRE_PARCSR); + HYPRE_IJMatrixInitialize(A_hypre_); + + // Copy values from Dune matrix to Hypre matrix + for (auto row = A_.begin(); row != A_.end(); ++row) { + int i = row.index(); + std::vector cols; + std::vector vals; + + for (auto col = row->begin(); col != row->end(); ++col) { + cols.push_back(col.index()); + vals.push_back(*col); + } + + int ncols = cols.size(); + HYPRE_IJMatrixSetValues(A_hypre_, 1, &ncols, &i, cols.data(), vals.data()); + } + + HYPRE_IJMatrixAssemble(A_hypre_); + HYPRE_IJMatrixGetObject(A_hypre_, (void**)&parcsr_A_); + + // Setup the solver with the new matrix + HYPRE_BoomerAMGSetup(solver_, parcsr_A_, par_b_, par_x_); + } + + void copyVectorsToHypre(const X& v, const Y& d) { + const int N = A_.N(); + + // Copy values + std::vector indices(N); + std::vector x_vals(N); + std::vector b_vals(N); + + for (int i = 0; i < N; ++i) { + indices[i] = i; + x_vals[i] = v[i]; + b_vals[i] = d[i]; + } + + HYPRE_IJVectorSetValues(x_hypre_, N, indices.data(), x_vals.data()); + HYPRE_IJVectorSetValues(b_hypre_, N, indices.data(), b_vals.data()); + + HYPRE_IJVectorGetObject(x_hypre_, (void**)&par_x_); + HYPRE_IJVectorGetObject(b_hypre_, (void**)&par_b_); + } + + void copyVectorFromHypre(X& v) { + const int N = A_.N(); + std::vector indices(N); + std::vector vals(N); + + for (int i = 0; i < N; ++i) { + indices[i] = i; + } + + HYPRE_IJVectorGetValues(x_hypre_, N, indices.data(), vals.data()); + + for (int i = 0; i < N; ++i) { + v[i] = vals[i]; + } + } + + const M& A_; + HYPRE_Solver solver_ = nullptr; + HYPRE_IJMatrix A_hypre_ = nullptr; + HYPRE_ParCSRMatrix parcsr_A_ = nullptr; + HYPRE_IJVector x_hypre_ = nullptr; + HYPRE_IJVector b_hypre_ = nullptr; + HYPRE_ParVector par_x_ = nullptr; + HYPRE_ParVector par_b_ = nullptr; +}; + +} // namespace Dune + +#endif diff --git a/opm/simulators/linalg/PreconditionerFactory_impl.hpp b/opm/simulators/linalg/PreconditionerFactory_impl.hpp index c53c069e7..7b4e85285 100644 --- a/opm/simulators/linalg/PreconditionerFactory_impl.hpp +++ b/opm/simulators/linalg/PreconditionerFactory_impl.hpp @@ -49,6 +49,10 @@ // Include all cuistl/GPU preconditioners inside of this headerfile #include +#if HAVE_HYPRE +#include +#endif + namespace Opm { @@ -543,6 +547,13 @@ struct StandardPreconditioners { return getRebuildOnUpdateWrapper>(op, crit, parms); } }); + // Only add Hypre for scalar matrices + if constexpr (M::block_type::rows == 1 && M::block_type::cols == 1) { + F::addCreator("HypreBoomerAMG", [](const O& op, const P& prm, const std::function&, std::size_t) { + DUNE_UNUSED_PARAMETER(prm); + return std::make_shared>(op.getmat()); + }); + } } if constexpr (std::is_same_v>) { F::addCreator( From 84193fa53d3abb0d363a6e019f5fd306fbc8f895 Mon Sep 17 00:00:00 2001 From: jakobtorben Date: Wed, 6 Nov 2024 16:16:54 +0100 Subject: [PATCH 02/13] Set Hypre AMG parameters to Jutul default --- opm/simulators/linalg/HyprePreconditioner.hpp | 15 ++++++++------- 1 file changed, 8 insertions(+), 7 deletions(-) diff --git a/opm/simulators/linalg/HyprePreconditioner.hpp b/opm/simulators/linalg/HyprePreconditioner.hpp index a2e4d50b3..3ec62c24e 100644 --- a/opm/simulators/linalg/HyprePreconditioner.hpp +++ b/opm/simulators/linalg/HyprePreconditioner.hpp @@ -58,13 +58,14 @@ public: HYPRE_BoomerAMGCreate(&solver_); // Set some default parameters - HYPRE_BoomerAMGSetPrintLevel(solver_, 0); // Reduce output - HYPRE_BoomerAMGSetOldDefault(solver_); // Falgout coarsening with modified classical interpolation - HYPRE_BoomerAMGSetRelaxType(solver_, 3); // G-S/Jacobi hybrid relaxation - HYPRE_BoomerAMGSetRelaxOrder(solver_, 1); // Uses C/F relaxation - HYPRE_BoomerAMGSetNumSweeps(solver_, 1); // Sweeps on each level - HYPRE_BoomerAMGSetMaxLevels(solver_, 20); // Maximum number of levels - HYPRE_BoomerAMGSetTol(solver_, 1e-7); // Convergence tolerance + HYPRE_BoomerAMGSetPrintLevel(solver_, 0); // Reduce output + HYPRE_BoomerAMGSetCoarsenType(solver_, 10); // HMIS coarsening + HYPRE_BoomerAMGSetStrongThreshold(solver_, 0.5); // Strength threshold for 3D + HYPRE_BoomerAMGSetAggNumLevels(solver_, 1); // Aggressive coarsening on first level + HYPRE_BoomerAMGSetAggTruncFactor(solver_, 0.3); // Remove weak connections + HYPRE_BoomerAMGSetInterpType(solver_, 6); // ext+i interpolation + HYPRE_BoomerAMGSetMaxLevels(solver_, 15); // Maximum number of levels + HYPRE_BoomerAMGSetTol(solver_, 0.1); // Convergence tolerance // Create the vectors From 6d22b4f72ae2718090b8efbb786c94793ac84803 Mon Sep 17 00:00:00 2001 From: jakobtorben Date: Wed, 6 Nov 2024 18:06:45 +0100 Subject: [PATCH 03/13] Optimise the transfer of memory between Dune and Hypre --- opm/simulators/linalg/HyprePreconditioner.hpp | 95 ++++++++++--------- 1 file changed, 51 insertions(+), 44 deletions(-) diff --git a/opm/simulators/linalg/HyprePreconditioner.hpp b/opm/simulators/linalg/HyprePreconditioner.hpp index 3ec62c24e..54d0911cb 100644 --- a/opm/simulators/linalg/HyprePreconditioner.hpp +++ b/opm/simulators/linalg/HyprePreconditioner.hpp @@ -29,6 +29,7 @@ #include #include +#include namespace Dune { @@ -68,7 +69,7 @@ public: HYPRE_BoomerAMGSetTol(solver_, 0.1); // Convergence tolerance - // Create the vectors + // Create Hypre vectors const int N = A_.N(); HYPRE_IJVectorCreate(MPI_COMM_WORLD, 0, N-1, &x_hypre_); HYPRE_IJVectorCreate(MPI_COMM_WORLD, 0, N-1, &b_hypre_); @@ -76,7 +77,16 @@ public: HYPRE_IJVectorSetObjectType(b_hypre_, HYPRE_PARCSR); HYPRE_IJVectorInitialize(x_hypre_); HYPRE_IJVectorInitialize(b_hypre_); + // Create indices vector + indices_.resize(A_.N()); + std::iota(indices_.begin(), indices_.end(), 0); + // Create Hypre matrix + HYPRE_IJMatrixCreate(MPI_COMM_WORLD, 0, N-1, 0, N-1, &A_hypre_); + HYPRE_IJMatrixSetObjectType(A_hypre_, HYPRE_PARCSR); + HYPRE_IJMatrixInitialize(A_hypre_); + + setupSparsityPattern(); update(); } @@ -94,6 +104,7 @@ public: void update() override { OPM_TIMEBLOCK(prec_update); copyMatrixToHypre(); + HYPRE_BoomerAMGSetup(solver_, parcsr_A_, par_b_, par_x_); } void pre(X& x, Y& b) override { @@ -129,69 +140,57 @@ public: } private: - void copyMatrixToHypre() { + void setupSparsityPattern() { const int N = A_.N(); - HYPRE_IJMatrixCreate(MPI_COMM_WORLD, 0, N-1, 0, N-1, &A_hypre_); - HYPRE_IJMatrixSetObjectType(A_hypre_, HYPRE_PARCSR); - HYPRE_IJMatrixInitialize(A_hypre_); + const int nnz = A_.nonzeroes(); - // Copy values from Dune matrix to Hypre matrix + // Allocate arrays required by Hypre + ncols_.resize(N); + rows_.resize(N); + cols_.resize(nnz); + + // Setup arrays and fill column indices + int pos = 0; for (auto row = A_.begin(); row != A_.end(); ++row) { - int i = row.index(); - std::vector cols; - std::vector vals; + const int rowIdx = row.index(); + rows_[rowIdx] = rowIdx; + ncols_[rowIdx] = row->size(); for (auto col = row->begin(); col != row->end(); ++col) { - cols.push_back(col.index()); - vals.push_back(*col); + cols_[pos++] = col.index(); } - - int ncols = cols.size(); - HYPRE_IJMatrixSetValues(A_hypre_, 1, &ncols, &i, cols.data(), vals.data()); } + } + + void copyMatrixToHypre() { + // Get pointer to matrix values array + const double* vals = &(A_[0][0][0][0]); // Indexing explanation: + // A_[row] - First row of the matrix + // [0] - First block in that row + // [0] - First row within the 1x1 block + // [0] - First column within the 1x1 block + + // Set all values at once using stored sparsity pattern + HYPRE_IJMatrixSetValues(A_hypre_, A_.N(), ncols_.data(), rows_.data(), cols_.data(), vals); HYPRE_IJMatrixAssemble(A_hypre_); HYPRE_IJMatrixGetObject(A_hypre_, (void**)&parcsr_A_); - - // Setup the solver with the new matrix - HYPRE_BoomerAMGSetup(solver_, parcsr_A_, par_b_, par_x_); } void copyVectorsToHypre(const X& v, const Y& d) { - const int N = A_.N(); + const double* x_vals = &(v[0][0]); + const double* b_vals = &(d[0][0]); - // Copy values - std::vector indices(N); - std::vector x_vals(N); - std::vector b_vals(N); - - for (int i = 0; i < N; ++i) { - indices[i] = i; - x_vals[i] = v[i]; - b_vals[i] = d[i]; - } - - HYPRE_IJVectorSetValues(x_hypre_, N, indices.data(), x_vals.data()); - HYPRE_IJVectorSetValues(b_hypre_, N, indices.data(), b_vals.data()); + HYPRE_IJVectorSetValues(x_hypre_, A_.N(), indices_.data(), x_vals); + HYPRE_IJVectorSetValues(b_hypre_, A_.N(), indices_.data(), b_vals); HYPRE_IJVectorGetObject(x_hypre_, (void**)&par_x_); HYPRE_IJVectorGetObject(b_hypre_, (void**)&par_b_); } void copyVectorFromHypre(X& v) { - const int N = A_.N(); - std::vector indices(N); - std::vector vals(N); - - for (int i = 0; i < N; ++i) { - indices[i] = i; - } - - HYPRE_IJVectorGetValues(x_hypre_, N, indices.data(), vals.data()); - - for (int i = 0; i < N; ++i) { - v[i] = vals[i]; - } + double* vals = &(v[0][0]); + HYPRE_IJVectorGetValues(x_hypre_, A_.N(), indices_.data(), vals); } const M& A_; @@ -202,6 +201,14 @@ private: HYPRE_IJVector b_hypre_ = nullptr; HYPRE_ParVector par_x_ = nullptr; HYPRE_ParVector par_b_ = nullptr; + + // Store sparsity pattern + std::vector ncols_; + std::vector rows_; + std::vector cols_; + + // Store indices vector + std::vector indices_; }; } // namespace Dune From 9d47aa605ee6dfca05fff2441c6803dbdc9ece7a Mon Sep 17 00:00:00 2001 From: jakobtorben Date: Thu, 7 Nov 2024 16:08:25 +0100 Subject: [PATCH 04/13] Add property tree to constructor --- opm/simulators/linalg/HyprePreconditioner.hpp | 14 ++++++++------ .../linalg/PreconditionerFactory_impl.hpp | 5 ++--- 2 files changed, 10 insertions(+), 9 deletions(-) diff --git a/opm/simulators/linalg/HyprePreconditioner.hpp b/opm/simulators/linalg/HyprePreconditioner.hpp index 54d0911cb..5af6c3a5b 100644 --- a/opm/simulators/linalg/HyprePreconditioner.hpp +++ b/opm/simulators/linalg/HyprePreconditioner.hpp @@ -31,11 +31,11 @@ #include #include -namespace Dune { +namespace Hypre { /// Wrapper for Hypre's BoomerAMG preconditioner template -class HyprePreconditioner : public PreconditionerWithUpdate { +class HyprePreconditioner : public Dune::PreconditionerWithUpdate { public: //! \brief The matrix type the preconditioner is for using matrix_type = M; @@ -47,8 +47,8 @@ public: using field_type = typename X::field_type; // Constructor - HyprePreconditioner (const M& A) - : A_(A) + HyprePreconditioner (const M& A, const Opm::PropertyTree& prm) + : A_(A), prm_(prm) { OPM_TIMEBLOCK(prec_construct); @@ -129,8 +129,8 @@ public: DUNE_UNUSED_PARAMETER(x); } - SolverCategory::Category category() const override { - return SolverCategory::sequential; + Dune::SolverCategory::Category category() const override { + return Dune::SolverCategory::sequential; } bool hasPerfectUpdate() const override @@ -194,6 +194,8 @@ private: } const M& A_; + const Opm::PropertyTree& prm_; + HYPRE_Solver solver_ = nullptr; HYPRE_IJMatrix A_hypre_ = nullptr; HYPRE_ParCSRMatrix parcsr_A_ = nullptr; diff --git a/opm/simulators/linalg/PreconditionerFactory_impl.hpp b/opm/simulators/linalg/PreconditionerFactory_impl.hpp index 7b4e85285..4c56a9404 100644 --- a/opm/simulators/linalg/PreconditionerFactory_impl.hpp +++ b/opm/simulators/linalg/PreconditionerFactory_impl.hpp @@ -549,9 +549,8 @@ struct StandardPreconditioners { }); // Only add Hypre for scalar matrices if constexpr (M::block_type::rows == 1 && M::block_type::cols == 1) { - F::addCreator("HypreBoomerAMG", [](const O& op, const P& prm, const std::function&, std::size_t) { - DUNE_UNUSED_PARAMETER(prm); - return std::make_shared>(op.getmat()); + F::addCreator("hypre", [](const O& op, const P& prm, const std::function&, std::size_t) { + return std::make_shared>(op.getmat(), prm); }); } } From 30e3dc358214c4c1cb358a74eccc9e36b410003f Mon Sep 17 00:00:00 2001 From: jakobtorben Date: Mon, 11 Nov 2024 17:10:47 +0100 Subject: [PATCH 05/13] Only apply Hypre AMG one V-cycle --- opm/simulators/linalg/HyprePreconditioner.hpp | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/opm/simulators/linalg/HyprePreconditioner.hpp b/opm/simulators/linalg/HyprePreconditioner.hpp index 5af6c3a5b..1eb64a08d 100644 --- a/opm/simulators/linalg/HyprePreconditioner.hpp +++ b/opm/simulators/linalg/HyprePreconditioner.hpp @@ -60,13 +60,14 @@ public: // Set some default parameters HYPRE_BoomerAMGSetPrintLevel(solver_, 0); // Reduce output + HYPRE_BoomerAMGSetMaxIter(solver_, 1); // Only one V-cycle, as used as preconditioner HYPRE_BoomerAMGSetCoarsenType(solver_, 10); // HMIS coarsening HYPRE_BoomerAMGSetStrongThreshold(solver_, 0.5); // Strength threshold for 3D HYPRE_BoomerAMGSetAggNumLevels(solver_, 1); // Aggressive coarsening on first level HYPRE_BoomerAMGSetAggTruncFactor(solver_, 0.3); // Remove weak connections HYPRE_BoomerAMGSetInterpType(solver_, 6); // ext+i interpolation HYPRE_BoomerAMGSetMaxLevels(solver_, 15); // Maximum number of levels - HYPRE_BoomerAMGSetTol(solver_, 0.1); // Convergence tolerance + HYPRE_BoomerAMGSetTol(solver_, 0); // Convergence tolerance, 0 as used as preconditioner with one V-cycle // Create Hypre vectors From 032a7f5ad8047b08649f07f682f78579e9f8374b Mon Sep 17 00:00:00 2001 From: jakobtorben Date: Mon, 11 Nov 2024 17:31:57 +0100 Subject: [PATCH 06/13] Create coarse matrix with contiguous memory layout --- opm/simulators/linalg/PressureBhpTransferPolicy.hpp | 2 +- opm/simulators/linalg/PressureTransferPolicy.hpp | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/opm/simulators/linalg/PressureBhpTransferPolicy.hpp b/opm/simulators/linalg/PressureBhpTransferPolicy.hpp index b80f3c749..f8fedcb25 100644 --- a/opm/simulators/linalg/PressureBhpTransferPolicy.hpp +++ b/opm/simulators/linalg/PressureBhpTransferPolicy.hpp @@ -139,7 +139,7 @@ namespace Opm } } else { coarseLevelMatrix_.reset( - new CoarseMatrix(fineLevelMatrix.N(), fineLevelMatrix.M(), CoarseMatrix::row_wise)); + new CoarseMatrix(fineLevelMatrix.N(), fineLevelMatrix.M(), fineLevelMatrix.nonzeroes(), CoarseMatrix::row_wise)); auto createIter = coarseLevelMatrix_->createbegin(); for (const auto& row : fineLevelMatrix) { for (auto col = row.begin(), cend = row.end(); col != cend; ++col) { diff --git a/opm/simulators/linalg/PressureTransferPolicy.hpp b/opm/simulators/linalg/PressureTransferPolicy.hpp index db52ba576..60d19904a 100644 --- a/opm/simulators/linalg/PressureTransferPolicy.hpp +++ b/opm/simulators/linalg/PressureTransferPolicy.hpp @@ -74,7 +74,7 @@ public: { using CoarseMatrix = typename CoarseOperator::matrix_type; const auto& fineLevelMatrix = fineOperator.getmat(); - coarseLevelMatrix_.reset(new CoarseMatrix(fineLevelMatrix.N(), fineLevelMatrix.M(), CoarseMatrix::row_wise)); + coarseLevelMatrix_.reset(new CoarseMatrix(fineLevelMatrix.N(), fineLevelMatrix.M(), fineLevelMatrix.nonzeroes(), CoarseMatrix::row_wise)); auto createIter = coarseLevelMatrix_->createbegin(); for (const auto& row : fineLevelMatrix) { From 1af2556385b4765ca5a2c18f21d5ad99f63d2114 Mon Sep 17 00:00:00 2001 From: jakobtorben Date: Wed, 20 Nov 2024 19:30:50 +0100 Subject: [PATCH 07/13] Add support for GPU with Hypre --- CMakeLists.txt | 11 +- opm/simulators/linalg/HyprePreconditioner.hpp | 179 ++++++++++++++---- 2 files changed, 148 insertions(+), 42 deletions(-) diff --git a/CMakeLists.txt b/CMakeLists.txt index cc789f43f..923ea2db3 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -649,8 +649,13 @@ add_custom_target(extra_test ${CMAKE_CTEST_COMMAND} -C ExtraTests) if(CUDA_FOUND) if (NOT USE_HIP) - target_link_libraries( opmsimulators PUBLIC ${CUDA_cusparse_LIBRARY} ) - target_link_libraries( opmsimulators PUBLIC ${CUDA_cublas_LIBRARY} ) + target_link_libraries(opmsimulators + PUBLIC + ${CUDA_cusparse_LIBRARY} + ${CUDA_cublas_LIBRARY} + ${CUDA_nvptxcompiler_static_LIBRARY} + ) + foreach(tgt test_gpu_safe_call test_cuda_check_last_error test_GpuVector) target_link_libraries(${tgt} CUDA::cudart) endforeach() @@ -729,4 +734,4 @@ endif() if(HYPRE_FOUND) target_link_libraries(opmsimulators PUBLIC HYPRE::HYPRE) -endif() +endif() \ No newline at end of file diff --git a/opm/simulators/linalg/HyprePreconditioner.hpp b/opm/simulators/linalg/HyprePreconditioner.hpp index 1eb64a08d..6435effa6 100644 --- a/opm/simulators/linalg/HyprePreconditioner.hpp +++ b/opm/simulators/linalg/HyprePreconditioner.hpp @@ -26,6 +26,7 @@ #include #include #include +#include <_hypre_utilities.h> #include #include @@ -52,38 +53,75 @@ public: { OPM_TIMEBLOCK(prec_construct); - // Initialize Hypre - HYPRE_Init(); + use_gpu_ = prm_.get("use_gpu", false); + + HYPRE_Initialize(); + + // Set memory location and execution policy + if (use_gpu_) { + HYPRE_SetMemoryLocation(HYPRE_MEMORY_DEVICE); + HYPRE_SetExecutionPolicy(HYPRE_EXEC_DEVICE); + // use hypre's SpGEMM instead of vendor implementation + HYPRE_SetSpGemmUseVendor(false); + // use cuRand for PMIS + HYPRE_SetUseGpuRand(1); + HYPRE_DeviceInitialize(); + HYPRE_PrintDeviceInfo(); + } + else { + HYPRE_SetMemoryLocation(HYPRE_MEMORY_HOST); + HYPRE_SetExecutionPolicy(HYPRE_EXEC_HOST); + } // Create the solver (BoomerAMG) HYPRE_BoomerAMGCreate(&solver_); - // Set some default parameters - HYPRE_BoomerAMGSetPrintLevel(solver_, 0); // Reduce output - HYPRE_BoomerAMGSetMaxIter(solver_, 1); // Only one V-cycle, as used as preconditioner - HYPRE_BoomerAMGSetCoarsenType(solver_, 10); // HMIS coarsening - HYPRE_BoomerAMGSetStrongThreshold(solver_, 0.5); // Strength threshold for 3D - HYPRE_BoomerAMGSetAggNumLevels(solver_, 1); // Aggressive coarsening on first level - HYPRE_BoomerAMGSetAggTruncFactor(solver_, 0.3); // Remove weak connections - HYPRE_BoomerAMGSetInterpType(solver_, 6); // ext+i interpolation - HYPRE_BoomerAMGSetMaxLevels(solver_, 15); // Maximum number of levels - HYPRE_BoomerAMGSetTol(solver_, 0); // Convergence tolerance, 0 as used as preconditioner with one V-cycle + // Set parameters from property tree with defaults + HYPRE_BoomerAMGSetPrintLevel(solver_, prm_.get("print_level", 0)); + HYPRE_BoomerAMGSetMaxIter(solver_, prm_.get("max_iter", 1)); + HYPRE_BoomerAMGSetStrongThreshold(solver_, prm_.get("strong_threshold", 0.5)); + HYPRE_BoomerAMGSetAggTruncFactor(solver_, prm_.get("agg_trunc_factor", 0.3)); + HYPRE_BoomerAMGSetInterpType(solver_, prm_.get("interp_type", 6)); + HYPRE_BoomerAMGSetMaxLevels(solver_, prm_.get("max_levels", 15)); + HYPRE_BoomerAMGSetTol(solver_, prm_.get("tolerance", 0.0)); + if (use_gpu_) { + HYPRE_BoomerAMGSetRelaxType(solver_, 16); + HYPRE_BoomerAMGSetCoarsenType(solver_, 8); + HYPRE_BoomerAMGSetAggNumLevels(solver_, 0); + HYPRE_BoomerAMGSetAggInterpType(solver_, 6); + // Keep transpose to avoid SpMTV + HYPRE_BoomerAMGSetKeepTranspose(solver_, true); + } + else { + HYPRE_BoomerAMGSetRelaxType(solver_, prm_.get("relax_type", 13)); + HYPRE_BoomerAMGSetCoarsenType(solver_, prm_.get("coarsen_type", 10)); + HYPRE_BoomerAMGSetAggNumLevels(solver_, prm_.get("agg_num_levels", 1)); + HYPRE_BoomerAMGSetAggInterpType(solver_, prm_.get("agg_interp_type", 4)); + } // Create Hypre vectors - const int N = A_.N(); - HYPRE_IJVectorCreate(MPI_COMM_WORLD, 0, N-1, &x_hypre_); - HYPRE_IJVectorCreate(MPI_COMM_WORLD, 0, N-1, &b_hypre_); + N_ = A_.N(); + nnz_ = A_.nonzeroes(); + HYPRE_IJVectorCreate(MPI_COMM_WORLD, 0, N_-1, &x_hypre_); + HYPRE_IJVectorCreate(MPI_COMM_WORLD, 0, N_-1, &b_hypre_); HYPRE_IJVectorSetObjectType(x_hypre_, HYPRE_PARCSR); HYPRE_IJVectorSetObjectType(b_hypre_, HYPRE_PARCSR); HYPRE_IJVectorInitialize(x_hypre_); HYPRE_IJVectorInitialize(b_hypre_); // Create indices vector - indices_.resize(A_.N()); + indices_.resize(N_); std::iota(indices_.begin(), indices_.end(), 0); + if (use_gpu_) { + indices_device_ = hypre_CTAlloc(HYPRE_BigInt, N_, HYPRE_MEMORY_DEVICE); + hypre_TMemcpy(indices_device_, indices_.data(), HYPRE_BigInt, N_, HYPRE_MEMORY_DEVICE, HYPRE_MEMORY_HOST); + // Allocate device vectors + x_values_device_ = hypre_CTAlloc(HYPRE_Real, N_, HYPRE_MEMORY_DEVICE); + b_values_device_ = hypre_CTAlloc(HYPRE_Real, N_, HYPRE_MEMORY_DEVICE); + } // Create Hypre matrix - HYPRE_IJMatrixCreate(MPI_COMM_WORLD, 0, N-1, 0, N-1, &A_hypre_); + HYPRE_IJMatrixCreate(MPI_COMM_WORLD, 0, N_-1, 0, N_-1, &A_hypre_); HYPRE_IJMatrixSetObjectType(A_hypre_, HYPRE_PARCSR); HYPRE_IJMatrixInitialize(A_hypre_); @@ -96,9 +134,27 @@ public: if (solver_) { HYPRE_BoomerAMGDestroy(solver_); } - if (parcsr_A_) { + if (A_hypre_) { HYPRE_IJMatrixDestroy(A_hypre_); } + if (x_hypre_) { + HYPRE_IJVectorDestroy(x_hypre_); + } + if (b_hypre_) { + HYPRE_IJVectorDestroy(b_hypre_); + } + if (values_device_) { + hypre_TFree(values_device_, HYPRE_MEMORY_DEVICE); + } + if (x_values_device_) { + hypre_TFree(x_values_device_, HYPRE_MEMORY_DEVICE); + } + if (b_values_device_) { + hypre_TFree(b_values_device_, HYPRE_MEMORY_DEVICE); + } + if (indices_device_) { + hypre_TFree(indices_device_, HYPRE_MEMORY_DEVICE); + } HYPRE_Finalize(); } @@ -142,13 +198,10 @@ public: private: void setupSparsityPattern() { - const int N = A_.N(); - const int nnz = A_.nonzeroes(); - // Allocate arrays required by Hypre - ncols_.resize(N); - rows_.resize(N); - cols_.resize(nnz); + ncols_.resize(N_); + rows_.resize(N_); + cols_.resize(nnz_); // Setup arrays and fill column indices int pos = 0; @@ -161,41 +214,80 @@ private: cols_[pos++] = col.index(); } } + if (use_gpu_) { + // Allocate device arrays + ncols_device_ = hypre_CTAlloc(HYPRE_Int, N_, HYPRE_MEMORY_DEVICE); + rows_device_ = hypre_CTAlloc(HYPRE_BigInt, N_, HYPRE_MEMORY_DEVICE); + cols_device_ = hypre_CTAlloc(HYPRE_BigInt, nnz_, HYPRE_MEMORY_DEVICE); + values_device_ = hypre_CTAlloc(HYPRE_Real, nnz_, HYPRE_MEMORY_DEVICE); + + // Copy to device + hypre_TMemcpy(ncols_device_, ncols_.data(), HYPRE_Int, N_, HYPRE_MEMORY_DEVICE, HYPRE_MEMORY_HOST); + hypre_TMemcpy(rows_device_, rows_.data(), HYPRE_BigInt, N_, HYPRE_MEMORY_DEVICE, HYPRE_MEMORY_HOST); + hypre_TMemcpy(cols_device_, cols_.data(), HYPRE_BigInt, nnz_, HYPRE_MEMORY_DEVICE, HYPRE_MEMORY_HOST); + } } void copyMatrixToHypre() { + OPM_TIMEBLOCK(prec_copy_matrix); // Get pointer to matrix values array - const double* vals = &(A_[0][0][0][0]); // Indexing explanation: - // A_[row] - First row of the matrix - // [0] - First block in that row - // [0] - First row within the 1x1 block - // [0] - First column within the 1x1 block + const HYPRE_Real* values = &(A_[0][0][0][0]); + // Indexing explanation: + // A_[row] - First row of the matrix + // [0] - First block in that row + // [0] - First row within the 1x1 block + // [0] - First column within the 1x1 block - // Set all values at once using stored sparsity pattern - HYPRE_IJMatrixSetValues(A_hypre_, A_.N(), ncols_.data(), rows_.data(), cols_.data(), vals); + if (use_gpu_) { + hypre_TMemcpy(values_device_, values, HYPRE_Real, nnz_, HYPRE_MEMORY_DEVICE, HYPRE_MEMORY_HOST); + HYPRE_IJMatrixSetValues(A_hypre_, N_, ncols_device_, rows_device_, cols_device_, values_device_); + } + else { + HYPRE_IJMatrixSetValues(A_hypre_, N_, ncols_.data(), rows_.data(), cols_.data(), values); + } HYPRE_IJMatrixAssemble(A_hypre_); HYPRE_IJMatrixGetObject(A_hypre_, (void**)&parcsr_A_); } void copyVectorsToHypre(const X& v, const Y& d) { - const double* x_vals = &(v[0][0]); - const double* b_vals = &(d[0][0]); + OPM_TIMEBLOCK(prec_copy_vectors_to_hypre); + const HYPRE_Real* x_vals = &(v[0][0]); + const HYPRE_Real* b_vals = &(d[0][0]); - HYPRE_IJVectorSetValues(x_hypre_, A_.N(), indices_.data(), x_vals); - HYPRE_IJVectorSetValues(b_hypre_, A_.N(), indices_.data(), b_vals); + if (use_gpu_) { + hypre_TMemcpy(x_values_device_, x_vals, HYPRE_Real, N_, HYPRE_MEMORY_DEVICE, HYPRE_MEMORY_HOST); + hypre_TMemcpy(b_values_device_, b_vals, HYPRE_Real, N_, HYPRE_MEMORY_DEVICE, HYPRE_MEMORY_HOST); + HYPRE_IJVectorSetValues(x_hypre_, N_, indices_device_, x_values_device_); + HYPRE_IJVectorSetValues(b_hypre_, N_, indices_device_, b_values_device_); + } + else { + HYPRE_IJVectorSetValues(x_hypre_, N_, indices_.data(), x_vals); + HYPRE_IJVectorSetValues(b_hypre_, N_, indices_.data(), b_vals); + } + + HYPRE_IJVectorAssemble(x_hypre_); + HYPRE_IJVectorAssemble(b_hypre_); HYPRE_IJVectorGetObject(x_hypre_, (void**)&par_x_); HYPRE_IJVectorGetObject(b_hypre_, (void**)&par_b_); } void copyVectorFromHypre(X& v) { - double* vals = &(v[0][0]); - HYPRE_IJVectorGetValues(x_hypre_, A_.N(), indices_.data(), vals); + OPM_TIMEBLOCK(prec_copy_vector_from_hypre); + HYPRE_Real* values = &(v[0][0]); + if (use_gpu_) { + HYPRE_IJVectorGetValues(x_hypre_, N_, indices_device_, x_values_device_); + hypre_TMemcpy(values, x_values_device_, HYPRE_Real, N_, HYPRE_MEMORY_HOST, HYPRE_MEMORY_DEVICE); + } + else { + HYPRE_IJVectorGetValues(x_hypre_, N_, indices_.data(), values); + } } const M& A_; const Opm::PropertyTree& prm_; + bool use_gpu_ = false; HYPRE_Solver solver_ = nullptr; HYPRE_IJMatrix A_hypre_ = nullptr; @@ -209,11 +301,20 @@ private: std::vector ncols_; std::vector rows_; std::vector cols_; - + HYPRE_Int* ncols_device_ = nullptr; + HYPRE_BigInt* rows_device_ = nullptr; + HYPRE_BigInt* cols_device_ = nullptr; + HYPRE_Real* values_device_ = nullptr; // Store indices vector std::vector indices_; + HYPRE_BigInt* indices_device_ = nullptr; + HYPRE_Int N_ = -1; + HYPRE_Int nnz_ = -1; + + HYPRE_Real* x_values_device_ = nullptr; + HYPRE_Real* b_values_device_ = nullptr; }; -} // namespace Dune +} // namespace Hypre #endif From 56897e6d11b72c8451269130caff24603e215dc7 Mon Sep 17 00:00:00 2001 From: jakobtorben Date: Mon, 25 Nov 2024 16:56:59 +0100 Subject: [PATCH 08/13] Guard Hypre preconditioner creation with HAVE_HYPRE macro --- opm/simulators/linalg/PreconditionerFactory_impl.hpp | 2 ++ 1 file changed, 2 insertions(+) diff --git a/opm/simulators/linalg/PreconditionerFactory_impl.hpp b/opm/simulators/linalg/PreconditionerFactory_impl.hpp index 4c56a9404..2d7607a3a 100644 --- a/opm/simulators/linalg/PreconditionerFactory_impl.hpp +++ b/opm/simulators/linalg/PreconditionerFactory_impl.hpp @@ -547,12 +547,14 @@ struct StandardPreconditioners { return getRebuildOnUpdateWrapper>(op, crit, parms); } }); +#if HAVE_HYPRE // Only add Hypre for scalar matrices if constexpr (M::block_type::rows == 1 && M::block_type::cols == 1) { F::addCreator("hypre", [](const O& op, const P& prm, const std::function&, std::size_t) { return std::make_shared>(op.getmat(), prm); }); } +#endif } if constexpr (std::is_same_v>) { F::addCreator( From 6fa9c25ba53f69f6c23585c4c764a0e7c1a5ff9b Mon Sep 17 00:00:00 2001 From: jakobtorben Date: Tue, 26 Nov 2024 14:55:27 +0100 Subject: [PATCH 09/13] Move Hypre initialization and finalization out to Main --- opm/simulators/flow/Main.cpp | 12 ++++++++++++ opm/simulators/linalg/HyprePreconditioner.hpp | 3 --- 2 files changed, 12 insertions(+), 3 deletions(-) diff --git a/opm/simulators/flow/Main.cpp b/opm/simulators/flow/Main.cpp index 1d40593ed..963fee78b 100644 --- a/opm/simulators/flow/Main.cpp +++ b/opm/simulators/flow/Main.cpp @@ -35,6 +35,10 @@ #include #endif +#if HAVE_HYPRE +#include +#endif + namespace Opm { Main::Main(int argc, char** argv, bool ownMPI) @@ -88,6 +92,10 @@ Main::~Main() } #endif // HAVE_MPI +#if HAVE_HYPRE + HYPRE_Finalize(); +#endif + if (ownMPI_) { FlowGenericVanguard::setCommunication(nullptr); } @@ -163,6 +171,10 @@ void Main::initMPI() #endif #endif // HAVE_MPI + +#if HAVE_HYPRE + HYPRE_Initialize(); +#endif } void Main::handleVersionCmdLine_(int argc, char** argv, diff --git a/opm/simulators/linalg/HyprePreconditioner.hpp b/opm/simulators/linalg/HyprePreconditioner.hpp index 6435effa6..c773fd6ab 100644 --- a/opm/simulators/linalg/HyprePreconditioner.hpp +++ b/opm/simulators/linalg/HyprePreconditioner.hpp @@ -55,8 +55,6 @@ public: use_gpu_ = prm_.get("use_gpu", false); - HYPRE_Initialize(); - // Set memory location and execution policy if (use_gpu_) { HYPRE_SetMemoryLocation(HYPRE_MEMORY_DEVICE); @@ -155,7 +153,6 @@ public: if (indices_device_) { hypre_TFree(indices_device_, HYPRE_MEMORY_DEVICE); } - HYPRE_Finalize(); } void update() override { From 3b67d6dc54cb195b3a1395157fa6b77e1a3475d2 Mon Sep 17 00:00:00 2001 From: Arne Morten Kvarving Date: Wed, 4 Dec 2024 11:31:04 +0100 Subject: [PATCH 10/13] Small fixes for hypre integration (#1) * use appropriate hypre init function * properly handle hypre without device support * don't add hypre preconditioner to factory if scalar type does not match --- opm/simulators/flow/Main.cpp | 5 +++++ opm/simulators/linalg/HyprePreconditioner.hpp | 5 ++++- opm/simulators/linalg/PreconditionerFactory_impl.hpp | 3 ++- 3 files changed, 11 insertions(+), 2 deletions(-) diff --git a/opm/simulators/flow/Main.cpp b/opm/simulators/flow/Main.cpp index 963fee78b..91864fba9 100644 --- a/opm/simulators/flow/Main.cpp +++ b/opm/simulators/flow/Main.cpp @@ -36,6 +36,7 @@ #endif #if HAVE_HYPRE +#include #include #endif @@ -173,7 +174,11 @@ void Main::initMPI() #endif // HAVE_MPI #if HAVE_HYPRE +#if HYPRE_RELEASE_NUMBER >= 22900 HYPRE_Initialize(); +#else + HYPRE_Init(); +#endif #endif } diff --git a/opm/simulators/linalg/HyprePreconditioner.hpp b/opm/simulators/linalg/HyprePreconditioner.hpp index c773fd6ab..9539967c7 100644 --- a/opm/simulators/linalg/HyprePreconditioner.hpp +++ b/opm/simulators/linalg/HyprePreconditioner.hpp @@ -56,6 +56,7 @@ public: use_gpu_ = prm_.get("use_gpu", false); // Set memory location and execution policy +#if HYPRE_USING_CUDA || HYPRE_USING_HIP if (use_gpu_) { HYPRE_SetMemoryLocation(HYPRE_MEMORY_DEVICE); HYPRE_SetExecutionPolicy(HYPRE_EXEC_DEVICE); @@ -66,7 +67,9 @@ public: HYPRE_DeviceInitialize(); HYPRE_PrintDeviceInfo(); } - else { + else +#endif + { HYPRE_SetMemoryLocation(HYPRE_MEMORY_HOST); HYPRE_SetExecutionPolicy(HYPRE_EXEC_HOST); } diff --git a/opm/simulators/linalg/PreconditionerFactory_impl.hpp b/opm/simulators/linalg/PreconditionerFactory_impl.hpp index 2d7607a3a..04bedf0e2 100644 --- a/opm/simulators/linalg/PreconditionerFactory_impl.hpp +++ b/opm/simulators/linalg/PreconditionerFactory_impl.hpp @@ -549,7 +549,8 @@ struct StandardPreconditioners { }); #if HAVE_HYPRE // Only add Hypre for scalar matrices - if constexpr (M::block_type::rows == 1 && M::block_type::cols == 1) { + if constexpr (M::block_type::rows == 1 && M::block_type::cols == 1 && + std::is_same_v) { F::addCreator("hypre", [](const O& op, const P& prm, const std::function&, std::size_t) { return std::make_shared>(op.getmat(), prm); }); From e60aa7da135c5303cd380861f34a10042efdf3f8 Mon Sep 17 00:00:00 2001 From: jakobtorben Date: Wed, 4 Dec 2024 15:22:14 +0100 Subject: [PATCH 11/13] Propagate need for recreating solver properly for two level methods --- opm/simulators/linalg/HyprePreconditioner.hpp | 6 ++++-- opm/simulators/linalg/OwningTwoLevelPreconditioner.hpp | 2 +- opm/simulators/linalg/PressureSolverPolicy.hpp | 5 +++++ opm/simulators/linalg/twolevelmethodcpr.hh | 6 ++++++ 4 files changed, 16 insertions(+), 3 deletions(-) diff --git a/opm/simulators/linalg/HyprePreconditioner.hpp b/opm/simulators/linalg/HyprePreconditioner.hpp index 9539967c7..17bbfec09 100644 --- a/opm/simulators/linalg/HyprePreconditioner.hpp +++ b/opm/simulators/linalg/HyprePreconditioner.hpp @@ -192,8 +192,10 @@ public: bool hasPerfectUpdate() const override { - // The Hypre preconditioner can depend on the values of the matrix, so it must be recreated - return false; + // The Hypre preconditioner can depend on the values of the matrix so it does not have perfect update. + // However, copying the matrix to Hypre requires to setup the solver again, so this is handled internally. + // So for ISTLSolver, we can return true. + return true; } private: diff --git a/opm/simulators/linalg/OwningTwoLevelPreconditioner.hpp b/opm/simulators/linalg/OwningTwoLevelPreconditioner.hpp index 5f32cccc6..98ab1cac1 100644 --- a/opm/simulators/linalg/OwningTwoLevelPreconditioner.hpp +++ b/opm/simulators/linalg/OwningTwoLevelPreconditioner.hpp @@ -174,7 +174,7 @@ public: } virtual bool hasPerfectUpdate() const override { - return false; + return twolevel_method_.hasPerfectUpdate(); } private: diff --git a/opm/simulators/linalg/PressureSolverPolicy.hpp b/opm/simulators/linalg/PressureSolverPolicy.hpp index 684ed1c77..0669a725a 100644 --- a/opm/simulators/linalg/PressureSolverPolicy.hpp +++ b/opm/simulators/linalg/PressureSolverPolicy.hpp @@ -88,6 +88,11 @@ namespace Amg return linsolver_->category(); } + bool hasPerfectUpdate() const + { + return linsolver_->preconditioner().hasPerfectUpdate(); + } + void apply(X& x, X& b, double reduction, Dune::InverseOperatorResult& res) override { linsolver_->apply(x, b, reduction, res); diff --git a/opm/simulators/linalg/twolevelmethodcpr.hh b/opm/simulators/linalg/twolevelmethodcpr.hh index 97435c82c..fbf8638ed 100644 --- a/opm/simulators/linalg/twolevelmethodcpr.hh +++ b/opm/simulators/linalg/twolevelmethodcpr.hh @@ -495,6 +495,12 @@ public: { } + bool hasPerfectUpdate() const + { + // The two-level method has perfect update if both the finesmoother and coarse solver do. + return smoother_->hasPerfectUpdate() && coarseSolver_->hasPerfectUpdate(); + } + void apply(FineDomainType& v, const FineRangeType& d) { FineDomainType u(v); From d663f16bea69482f8de76ae26fb3d7f8c58f8f84 Mon Sep 17 00:00:00 2001 From: jakobtorben Date: Wed, 4 Dec 2024 19:49:01 +0100 Subject: [PATCH 12/13] Add tests for Hypre Preconditioner --- CMakeLists.txt | 10 +- CMakeLists_files.cmake | 7 + opm/simulators/linalg/HyprePreconditioner.hpp | 28 ++-- tests/HyprePreconditionerTestHelper.hpp | 123 ++++++++++++++++++ tests/test_HyprePreconditionerCPU.cpp | 68 ++++++++++ tests/test_HyprePreconditionerGPU.cpp | 60 +++++++++ 6 files changed, 282 insertions(+), 14 deletions(-) create mode 100644 tests/HyprePreconditionerTestHelper.hpp create mode 100644 tests/test_HyprePreconditionerCPU.cpp create mode 100644 tests/test_HyprePreconditionerGPU.cpp diff --git a/CMakeLists.txt b/CMakeLists.txt index 923ea2db3..30cfe278f 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -726,12 +726,14 @@ if(USE_HYPRE) find_package(HYPRE) if(HYPRE_FOUND) set(HAVE_HYPRE 1) + target_link_libraries(opmsimulators PUBLIC HYPRE::HYPRE) + if (HYPRE_USING_CUDA) + set_tests_properties(HyprePreconditionerGPU PROPERTIES LABELS gpu_cuda) + elseif (HYPRE_USING_HIP) + set_tests_properties(HyprePreconditionerGPU PROPERTIES LABELS gpu_hip) + endif() else() message(WARNING "Hypre requested but not found. Continuing without Hypre support.") set(USE_HYPRE OFF) endif() endif() - -if(HYPRE_FOUND) - target_link_libraries(opmsimulators PUBLIC HYPRE::HYPRE) -endif() \ No newline at end of file diff --git a/CMakeLists_files.cmake b/CMakeLists_files.cmake index e3ece638b..e67d7667e 100644 --- a/CMakeLists_files.cmake +++ b/CMakeLists_files.cmake @@ -492,6 +492,13 @@ if(HDF5_FOUND) list(APPEND TEST_SOURCE_FILES tests/test_HDF5Serializer.cpp) endif() +if(HYPRE_FOUND) + list(APPEND TEST_SOURCE_FILES tests/test_HyprePreconditionerCPU.cpp) + if(HYPRE_USING_CUDA OR HYPRE_USING_HIP) + list(APPEND TEST_SOURCE_FILES tests/test_HyprePreconditionerGPU.cpp) + endif() +endif() + list (APPEND TEST_DATA_FILES tests/equil_base.DATA tests/equil_capillary.DATA diff --git a/opm/simulators/linalg/HyprePreconditioner.hpp b/opm/simulators/linalg/HyprePreconditioner.hpp index 17bbfec09..1f1919ca2 100644 --- a/opm/simulators/linalg/HyprePreconditioner.hpp +++ b/opm/simulators/linalg/HyprePreconditioner.hpp @@ -1,14 +1,19 @@ /* Copyright 2024 SINTEF AS + Copyright 2024 Equinor ASA + This file is part of the Open Porous Media project (OPM). + OPM is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version. + OPM is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details. + You should have received a copy of the GNU General Public License along with OPM. If not, see . */ @@ -19,6 +24,7 @@ #include #include #include +#include #include #include @@ -28,7 +34,6 @@ #include #include <_hypre_utilities.h> -#include #include #include @@ -48,11 +53,17 @@ public: using field_type = typename X::field_type; // Constructor - HyprePreconditioner (const M& A, const Opm::PropertyTree& prm) + HyprePreconditioner (const M& A, const Opm::PropertyTree prm) : A_(A), prm_(prm) { OPM_TIMEBLOCK(prec_construct); + int size; + MPI_Comm_size(MPI_COMM_WORLD, &size); + if (size > 1) { + OPM_THROW(std::runtime_error, "HyprePreconditioner is currently only implemented for sequential runs"); + } + use_gpu_ = prm_.get("use_gpu", false); // Set memory location and execution policy @@ -104,8 +115,8 @@ public: // Create Hypre vectors N_ = A_.N(); nnz_ = A_.nonzeroes(); - HYPRE_IJVectorCreate(MPI_COMM_WORLD, 0, N_-1, &x_hypre_); - HYPRE_IJVectorCreate(MPI_COMM_WORLD, 0, N_-1, &b_hypre_); + HYPRE_IJVectorCreate(MPI_COMM_SELF, 0, N_-1, &x_hypre_); + HYPRE_IJVectorCreate(MPI_COMM_SELF, 0, N_-1, &b_hypre_); HYPRE_IJVectorSetObjectType(x_hypre_, HYPRE_PARCSR); HYPRE_IJVectorSetObjectType(b_hypre_, HYPRE_PARCSR); HYPRE_IJVectorInitialize(x_hypre_); @@ -122,7 +133,7 @@ public: } // Create Hypre matrix - HYPRE_IJMatrixCreate(MPI_COMM_WORLD, 0, N_-1, 0, N_-1, &A_hypre_); + HYPRE_IJMatrixCreate(MPI_COMM_SELF, 0, N_-1, 0, N_-1, &A_hypre_); HYPRE_IJMatrixSetObjectType(A_hypre_, HYPRE_PARCSR); HYPRE_IJMatrixInitialize(A_hypre_); @@ -164,9 +175,7 @@ public: HYPRE_BoomerAMGSetup(solver_, parcsr_A_, par_b_, par_x_); } - void pre(X& x, Y& b) override { - DUNE_UNUSED_PARAMETER(x); - DUNE_UNUSED_PARAMETER(b); + void pre(X& /*x*/, Y& /*b*/) override { } void apply(X& v, const Y& d) override { @@ -182,8 +191,7 @@ public: copyVectorFromHypre(v); } - void post(X& x) override { - DUNE_UNUSED_PARAMETER(x); + void post(X& /*x*/) override { } Dune::SolverCategory::Category category() const override { diff --git a/tests/HyprePreconditionerTestHelper.hpp b/tests/HyprePreconditionerTestHelper.hpp new file mode 100644 index 000000000..af817b4bf --- /dev/null +++ b/tests/HyprePreconditionerTestHelper.hpp @@ -0,0 +1,123 @@ +/* + Copyright 2024 SINTEF AS + Copyright 2024 Equinor ASA + + This file is part of the Open Porous Media project (OPM). + + OPM is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + OPM is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with OPM. If not, see . +*/ + +#ifndef TEST_HYPREPRECONDITIONER_HELPER_HPP +#define TEST_HYPREPRECONDITIONER_HELPER_HPP + +#include + +#include +#include +#include +#include +#include + +#include +#include + + +template +void setupLaplace2d(int N, Matrix& mat) +{ + const int nonZeroes = N*N * 5; // max 5 entries per row (diagonal + 4 neighbors) + mat.setBuildMode(Matrix::row_wise); + mat.setSize(N*N, N*N, nonZeroes); + + // Set up sparsity pattern + for (auto row = mat.createbegin(); row != mat.createend(); ++row) { + const int i = row.index(); + int x = i % N; + int y = i / N; + + row.insert(i); // diagonal + if (x > 0) // left neighbor + row.insert(i-1); + if (x < N-1) // right neighbor + row.insert(i+1); + if (y > 0) // upper neighbor + row.insert(i-N); + if (y < N-1) // lower neighbor + row.insert(i+N); + } + + // Fill the matrix with values + for (auto row = mat.begin(); row != mat.end(); ++row) { + const int i = row.index(); + int x = i % N; + int y = i / N; + + // Set diagonal + (*row)[i] = 4.0; + + // Set off-diagonal entries + if (x > 0) // left neighbor + (*row)[i-1] = -1.0; + if (x < N-1) // right neighbor + (*row)[i+1] = -1.0; + if (y > 0) // upper neighbor + (*row)[i-N] = -1.0; + if (y < N-1) // lower neighbor + (*row)[i+N] = -1.0; + } +} + +inline void testHyprePreconditioner(bool use_gpu) +{ + constexpr int N = 100; // 100x100 grid + using Matrix = Dune::BCRSMatrix>; + using Vector = Dune::BlockVector>; + + // Create matrix + Matrix matrix; + setupLaplace2d(N, matrix); + + // Create vectors + Vector x(N*N), b(N*N); + x = 100.0; // Initial guess + b = 0.0; // RHS + + // Create operator + using Operator = Dune::MatrixAdapter; + Operator op(matrix); + + // Set up HYPRE parameters + Opm::PropertyTree prm; + prm.put("use_gpu", use_gpu ? 1 : 0); + + // Create preconditioner + auto prec = std::make_shared>(matrix, prm); + + // Create solver + double reduction = 1e-8; + int maxit = 300; + int verbosity = 0; + Dune::LoopSolver solver(op, *prec, reduction, maxit, verbosity); + + // Solve + Dune::InverseOperatorResult res; + solver.apply(x, b, res); + + // Check convergence + BOOST_CHECK(res.converged); + BOOST_CHECK_LT(res.reduction, 1e-8); +} + + +#endif // TEST_HYPREPRECONDITIONER_HELPER_HPP diff --git a/tests/test_HyprePreconditionerCPU.cpp b/tests/test_HyprePreconditionerCPU.cpp new file mode 100644 index 000000000..0a029da9b --- /dev/null +++ b/tests/test_HyprePreconditionerCPU.cpp @@ -0,0 +1,68 @@ +/* + Copyright 2024 SINTEF AS + Copyright 2024 Equinor ASA + + This file is part of the Open Porous Media project (OPM). + + OPM is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + OPM is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with OPM. If not, see . +*/ + +#include + +#define BOOST_TEST_MODULE TestHyprePreconditionerCPU +#define BOOST_TEST_NO_MAIN + +#include +#include "MpiFixture.hpp" +#include "HyprePreconditionerTestHelper.hpp" + +#include + +#include +#include +#include +#include +#include + +#include +#include + +BOOST_GLOBAL_FIXTURE(MPIFixture); + +BOOST_AUTO_TEST_CASE(TestHyprePreconditionerCPU) +{ + testHyprePreconditioner(false); +} + +bool init_unit_test_func() +{ + return true; +} + +int main(int argc, char** argv) +{ + Dune::MPIHelper::instance(argc, argv); + +#if HYPRE_RELEASE_NUMBER >= 22900 + HYPRE_Initialize(); +#else + HYPRE_Init(); +#endif + + int result = boost::unit_test::unit_test_main(&init_unit_test_func, argc, argv); + + HYPRE_Finalize(); + + return result; +} diff --git a/tests/test_HyprePreconditionerGPU.cpp b/tests/test_HyprePreconditionerGPU.cpp new file mode 100644 index 000000000..e6630f86d --- /dev/null +++ b/tests/test_HyprePreconditionerGPU.cpp @@ -0,0 +1,60 @@ +/* + Copyright 2024 SINTEF AS + Copyright 2024 Equinor ASA + + This file is part of the Open Porous Media project (OPM). + + OPM is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + OPM is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with OPM. If not, see . +*/ + +#include + +#define BOOST_TEST_MODULE TestHyprePreconditionerGPU +#define BOOST_TEST_NO_MAIN +#include + +#include +#include + +#include "MpiFixture.hpp" +#include "HyprePreconditionerTestHelper.hpp" + +BOOST_GLOBAL_FIXTURE(MPIFixture); + +BOOST_AUTO_TEST_CASE(TestHyprePreconditionerGPU) +{ + testHyprePreconditioner(true); +} + +bool init_unit_test_func() +{ + return true; +} + +int main(int argc, char** argv) +{ + Dune::MPIHelper::instance(argc, argv); + +#if HYPRE_RELEASE_NUMBER >= 22900 + HYPRE_Initialize(); +#else + HYPRE_Init(); +#endif + + int result = boost::unit_test::unit_test_main(&init_unit_test_func, argc, argv); + + HYPRE_Finalize(); + + return result; +} From 1536de038b63c069da7687fde162de53079b63c4 Mon Sep 17 00:00:00 2001 From: jakobtorben Date: Thu, 5 Dec 2024 12:56:58 +0100 Subject: [PATCH 13/13] Document the HyprePreconditioner class --- opm/simulators/linalg/HyprePreconditioner.hpp | 181 +++++++++++++----- 1 file changed, 129 insertions(+), 52 deletions(-) diff --git a/opm/simulators/linalg/HyprePreconditioner.hpp b/opm/simulators/linalg/HyprePreconditioner.hpp index 1f1919ca2..067fe6e66 100644 --- a/opm/simulators/linalg/HyprePreconditioner.hpp +++ b/opm/simulators/linalg/HyprePreconditioner.hpp @@ -39,22 +39,30 @@ namespace Hypre { -/// Wrapper for Hypre's BoomerAMG preconditioner +/** + * @brief Wrapper for Hypre's BoomerAMG preconditioner. + * + * This class provides an interface to the BoomerAMG preconditioner from the Hypre library. + * It is designed to work with matrices, update vectors, and defect vectors specified by the template parameters. + * + * @tparam M The matrix type the preconditioner is for. + * @tparam X The type of the update vector. + * @tparam Y The type of the defect vector. + */ template class HyprePreconditioner : public Dune::PreconditionerWithUpdate { public: - //! \brief The matrix type the preconditioner is for - using matrix_type = M; - //! \brief The domain type of the preconditioner - using domain_type = X; - //! \brief The range type of the preconditioner - using range_type = Y; - //! \brief The field type of the preconditioner - using field_type = typename X::field_type; - // Constructor + /** + * @brief Constructor for the HyprePreconditioner class. + * + * Initializes the preconditioner with the given matrix and property tree. + * + * @param A The matrix for which the preconditioner is constructed. + * @param prm The property tree containing configuration parameters. + */ HyprePreconditioner (const M& A, const Opm::PropertyTree prm) - : A_(A), prm_(prm) + : A_(A) { OPM_TIMEBLOCK(prec_construct); @@ -64,7 +72,7 @@ public: OPM_THROW(std::runtime_error, "HyprePreconditioner is currently only implemented for sequential runs"); } - use_gpu_ = prm_.get("use_gpu", false); + use_gpu_ = prm.get("use_gpu", false); // Set memory location and execution policy #if HYPRE_USING_CUDA || HYPRE_USING_HIP @@ -89,13 +97,13 @@ public: HYPRE_BoomerAMGCreate(&solver_); // Set parameters from property tree with defaults - HYPRE_BoomerAMGSetPrintLevel(solver_, prm_.get("print_level", 0)); - HYPRE_BoomerAMGSetMaxIter(solver_, prm_.get("max_iter", 1)); - HYPRE_BoomerAMGSetStrongThreshold(solver_, prm_.get("strong_threshold", 0.5)); - HYPRE_BoomerAMGSetAggTruncFactor(solver_, prm_.get("agg_trunc_factor", 0.3)); - HYPRE_BoomerAMGSetInterpType(solver_, prm_.get("interp_type", 6)); - HYPRE_BoomerAMGSetMaxLevels(solver_, prm_.get("max_levels", 15)); - HYPRE_BoomerAMGSetTol(solver_, prm_.get("tolerance", 0.0)); + HYPRE_BoomerAMGSetPrintLevel(solver_, prm.get("print_level", 0)); + HYPRE_BoomerAMGSetMaxIter(solver_, prm.get("max_iter", 1)); + HYPRE_BoomerAMGSetStrongThreshold(solver_, prm.get("strong_threshold", 0.5)); + HYPRE_BoomerAMGSetAggTruncFactor(solver_, prm.get("agg_trunc_factor", 0.3)); + HYPRE_BoomerAMGSetInterpType(solver_, prm.get("interp_type", 6)); + HYPRE_BoomerAMGSetMaxLevels(solver_, prm.get("max_levels", 15)); + HYPRE_BoomerAMGSetTol(solver_, prm.get("tolerance", 0.0)); if (use_gpu_) { HYPRE_BoomerAMGSetRelaxType(solver_, 16); @@ -106,10 +114,10 @@ public: HYPRE_BoomerAMGSetKeepTranspose(solver_, true); } else { - HYPRE_BoomerAMGSetRelaxType(solver_, prm_.get("relax_type", 13)); - HYPRE_BoomerAMGSetCoarsenType(solver_, prm_.get("coarsen_type", 10)); - HYPRE_BoomerAMGSetAggNumLevels(solver_, prm_.get("agg_num_levels", 1)); - HYPRE_BoomerAMGSetAggInterpType(solver_, prm_.get("agg_interp_type", 4)); + HYPRE_BoomerAMGSetRelaxType(solver_, prm.get("relax_type", 13)); + HYPRE_BoomerAMGSetCoarsenType(solver_, prm.get("coarsen_type", 10)); + HYPRE_BoomerAMGSetAggNumLevels(solver_, prm.get("agg_num_levels", 1)); + HYPRE_BoomerAMGSetAggInterpType(solver_, prm.get("agg_interp_type", 4)); } // Create Hypre vectors @@ -141,7 +149,11 @@ public: update(); } - // Destructor + /** + * @brief Destructor for the HyprePreconditioner class. + * + * Cleans up resources allocated by the preconditioner. + */ ~HyprePreconditioner() { if (solver_) { HYPRE_BoomerAMGDestroy(solver_); @@ -169,15 +181,36 @@ public: } } + /** + * @brief Updates the preconditioner with the current matrix values. + * + * This method should be called whenever the matrix values change. + */ void update() override { OPM_TIMEBLOCK(prec_update); copyMatrixToHypre(); HYPRE_BoomerAMGSetup(solver_, parcsr_A_, par_b_, par_x_); } - void pre(X& /*x*/, Y& /*b*/) override { + /** + * @brief Pre-processing step before applying the preconditioner. + * + * This method is currently a no-op. + * + * @param v The update vector. + * @param d The defect vector. + */ + void pre(X& /*v*/, Y& /*d*/) override { } + /** + * @brief Applies the preconditioner to a vector. + * + * Performs one AMG V-cycle to solve the system. + * + * @param v The update vector. + * @param d The defect vector. + */ void apply(X& v, const Y& d) override { OPM_TIMEBLOCK(prec_apply); @@ -191,13 +224,30 @@ public: copyVectorFromHypre(v); } - void post(X& /*x*/) override { + /** + * @brief Post-processing step after applying the preconditioner. + * + * This method is currently a no-op. + * + * @param v The update vector. + */ + void post(X& /*v*/) override { } + /** + * @brief Returns the solver category. + * + * @return The solver category, which is sequential. + */ Dune::SolverCategory::Category category() const override { return Dune::SolverCategory::sequential; } + /** + * @brief Checks if the preconditioner has a perfect update. + * + * @return True, indicating that the preconditioner can be perfectly updated. + */ bool hasPerfectUpdate() const override { // The Hypre preconditioner can depend on the values of the matrix so it does not have perfect update. @@ -207,6 +257,11 @@ public: } private: + /** + * @brief Sets up the sparsity pattern for the Hypre matrix. + * + * Allocates and initializes arrays required by Hypre. + */ void setupSparsityPattern() { // Allocate arrays required by Hypre ncols_.resize(N_); @@ -238,12 +293,19 @@ private: } } + /** + * @brief Copies the matrix values to the Hypre matrix. + * + * This method transfers the matrix data from the host to the Hypre matrix. + * It assumes that the values of the matrix are stored in a contiguous array. + * If GPU is used, the data is transferred to the device. + */ void copyMatrixToHypre() { OPM_TIMEBLOCK(prec_copy_matrix); // Get pointer to matrix values array const HYPRE_Real* values = &(A_[0][0][0][0]); // Indexing explanation: - // A_[row] - First row of the matrix + // A_[0] - First row of the matrix // [0] - First block in that row // [0] - First row within the 1x1 block // [0] - First column within the 1x1 block @@ -260,6 +322,15 @@ private: HYPRE_IJMatrixGetObject(A_hypre_, (void**)&parcsr_A_); } + /** + * @brief Copies vectors to the Hypre format. + * + * Transfers the update and defect vectors to Hypre. + * If GPU is used, the data is transferred from the host to the device. + * + * @param v The update vector. + * @param d The defect vector. + */ void copyVectorsToHypre(const X& v, const Y& d) { OPM_TIMEBLOCK(prec_copy_vectors_to_hypre); const HYPRE_Real* x_vals = &(v[0][0]); @@ -283,6 +354,14 @@ private: HYPRE_IJVectorGetObject(b_hypre_, (void**)&par_b_); } + /** + * @brief Copies the solution vector from Hypre. + * + * Transfers the solution vector from Hypre back to the host. + * If GPU is used, the data is transferred from the device to the host. + * + * @param v The update vector. + */ void copyVectorFromHypre(X& v) { OPM_TIMEBLOCK(prec_copy_vector_from_hypre); HYPRE_Real* values = &(v[0][0]); @@ -295,34 +374,32 @@ private: } } - const M& A_; - const Opm::PropertyTree& prm_; - bool use_gpu_ = false; + const M& A_; //!< The matrix for which the preconditioner is constructed. + bool use_gpu_ = false; //!< Flag indicating whether to use GPU acceleration. - HYPRE_Solver solver_ = nullptr; - HYPRE_IJMatrix A_hypre_ = nullptr; - HYPRE_ParCSRMatrix parcsr_A_ = nullptr; - HYPRE_IJVector x_hypre_ = nullptr; - HYPRE_IJVector b_hypre_ = nullptr; - HYPRE_ParVector par_x_ = nullptr; - HYPRE_ParVector par_b_ = nullptr; + HYPRE_Solver solver_ = nullptr; //!< The Hypre solver object. + HYPRE_IJMatrix A_hypre_ = nullptr; //!< The Hypre matrix object. + HYPRE_ParCSRMatrix parcsr_A_ = nullptr; //!< The parallel CSR matrix object. + HYPRE_IJVector x_hypre_ = nullptr; //!< The Hypre solution vector. + HYPRE_IJVector b_hypre_ = nullptr; //!< The Hypre right-hand side vector. + HYPRE_ParVector par_x_ = nullptr; //!< The parallel solution vector. + HYPRE_ParVector par_b_ = nullptr; //!< The parallel right-hand side vector. - // Store sparsity pattern - std::vector ncols_; - std::vector rows_; - std::vector cols_; - HYPRE_Int* ncols_device_ = nullptr; - HYPRE_BigInt* rows_device_ = nullptr; - HYPRE_BigInt* cols_device_ = nullptr; - HYPRE_Real* values_device_ = nullptr; - // Store indices vector - std::vector indices_; - HYPRE_BigInt* indices_device_ = nullptr; - HYPRE_Int N_ = -1; - HYPRE_Int nnz_ = -1; + std::vector ncols_; //!< Number of columns per row. + std::vector rows_; //!< Row indices. + std::vector cols_; //!< Column indices. + HYPRE_Int* ncols_device_ = nullptr; //!< Device array for number of columns per row. + HYPRE_BigInt* rows_device_ = nullptr; //!< Device array for row indices. + HYPRE_BigInt* cols_device_ = nullptr; //!< Device array for column indices. + HYPRE_Real* values_device_ = nullptr; //!< Device array for matrix values. - HYPRE_Real* x_values_device_ = nullptr; - HYPRE_Real* b_values_device_ = nullptr; + std::vector indices_; //!< Indices vector for copying vectors to/from Hypre. + HYPRE_BigInt* indices_device_ = nullptr; //!< Device array for indices. + HYPRE_Int N_ = -1; //!< Number of rows in the matrix. + HYPRE_Int nnz_ = -1; //!< Number of non-zero elements in the matrix. + + HYPRE_Real* x_values_device_ = nullptr; //!< Device array for solution vector values. + HYPRE_Real* b_values_device_ = nullptr; //!< Device array for right-hand side vector values. }; } // namespace Hypre