cudss.hpp

Go to the documentation of this file.

// Copyright (c) 2010-2026, Lawrence Livermore National Security, LLC. Produced
// at the Lawrence Livermore National Laboratory. All Rights reserved. See files
// LICENSE and NOTICE for details. LLNL-CODE-806117.
//
// This file is part of the MFEM library. For more information and source code
// availability visit https://mfem.org.
//
// MFEM is free software; you can redistribute it and/or modify it under the
// terms of the BSD-3 license. We welcome feedback and contributions, see file
// CONTRIBUTING.md for details.
 
#ifndef MFEM_CUDSS
#define MFEM_CUDSS
 
#include "../config/config.hpp"
 
#ifdef MFEM_USE_CUDSS
 
#include "cudss.h"
#include <memory>
 
#ifdef MFEM_USE_MPI
#include <mpi.h>
#include "hypre.hpp"
#else
#include "operator.hpp"
#include "sparsemat.hpp"
#endif
 
namespace mfem
{
/**
 * @brief cuDSS: A high-performance CUDA Library for Direct Sparse Solvers
 *
 * Interface for the distributed cuDSS solver
 */
class CuDSSSolver : public Solver
{
public:
   /// Specify the type of matrix we are applying the solver to
   enum MatType
   {
      /// CUDSS_MTYPE_GENERAL: General matrix [default].
      NONSYMMETRIC = 0,
      /// CUDSS_MTYPE_SYMMETRIC: Real symmetric matrix.
      SYMMETRIC_INDEFINITE = 1,
      /// CUDSS_MTYPE_SPD: Symmetric positive-definite matrix.
      SYMMETRIC_POSITIVE_DEFINITE = 2,
   };
 
   /// Specify the view type of matrix we are applying the solver to
   enum MatViewType
   {
      /// CUDSS_MVIEW_FULL: Full matrix [default]
      FULL = 0,
      /// CUDSS_MVIEW_LOWER: Lower-triangular matrix (including the diagonal).
      LOWER = 1,
      /// CUDSS_MVIEW_UPPER: Upper-triangular matrix (including the diagonal).
      UPPER = 2,
   };
 
   /**
    * @brief Constructor.
    */
   CuDSSSolver();
 
#ifdef MFEM_USE_MPI
   /**
    * @brief Constructor with MPI_Comm parameter.
    */
   CuDSSSolver(MPI_Comm comm);
#endif
 
   // Note: CuDSSSolver disables the move copy constructor and move assignment
   // operator
   CuDSSSolver(CuDSSSolver &&) = delete;
   CuDSSSolver &operator=(CuDSSSolver &&) = delete;
 
   /**
    * @brief Set the matrix type
    *
    * Supported matrix types:
    *          CuDSSSolver::NONSYMMETRIC,
    *          CuDSSSolver::SYMMETRIC_INDEFINITE,
    *      and CuDSSSolver::SYMMETRIC_POSITIVE_DEFINITE
    *
    * @param mtype_ Matrix type
    *
    * @note This method has to be called before SetOperator
    */
   void SetMatrixSymType(MatType mtype_);
 
   /**
    * @brief Set the matrix view type
    *
    * Supported matrix types:
    *          CuDSSSolver::FULL,
    *          CuDSSSolver::LOWER,
    *      and CuDSSSolver::UPPER
    *
    * @param mvtype Matrix view type
    *
    * @note This method has to be called before SetOperator
    */
   void SetMatrixViewType(MatViewType mvtype);
 
   /**
    * @brief Set the flag controlling reuse of the symbolic factorization
    * for multiple operators
    *
    * @param reuse Flag to reuse symbolic factorization
    *
    * @note This method has to be called before repeated calls to SetOperator
    */
   void SetReorderingReuse(bool reuse);
 
   void SetOperator(const Operator &op) override;
 
   /**
    * @brief Solve $ y = Op^{-1} x $
    *
    * @param x RHS vector
    * @param y Solution vector
    */
   void Mult(const Vector &x, Vector &y) const override;
 
   /**
    * @brief Solve $ Y_i = Op^{-1} X_i $
    *
    * @param X Array of RHS vectors
    * @param Y Array of Solution vectors
    */
   void ArrayMult(const Array<const Vector *> &X,
                  Array<Vector *> &Y) const override;
 
   ~CuDSSSolver();
 
private:
#ifdef MFEM_USE_MPI
   // MPI_Comm
   MPI_Comm mpi_comm = MPI_COMM_NULL;
 
   int row_start = 0;  // the first row index in CSR matrix operator
   int row_end = 0;    // the end row index in CSR matrix operator
#endif
 
   // Parameter controlling whether or not to reuse the symbolic factorization
   // for multiple calls to SetOperator
   bool reorder_reuse = false;
 
   // Parameter controlling the matrix type
   cudssMatrixType_t mat_type = CUDSS_MTYPE_GENERAL;
 
   int n_global = 0;      // global number of rows
   int n_loc = 0;         // the number of the rows in CSR matrix operator
 
   mutable int nrhs = 0;  // the number of the RHSs
   int nnz = 0;        // the number of non zeros
 
   // copy and keep the I and J arrays in device memory
   void *csr_offsets_d = NULL; // copy and keep I in device
   void *csr_columns_d = NULL; // copy and keep J in device
   void *csr_values_d = NULL;  // copy and keep csr data in device
 
   // cuDSS object specifies available matrix types for sparse matrices
   cudssMatrixViewType_t mview = CUDSS_MVIEW_FULL;
 
   // cuDSS objects storage for sparse matrix Ac, RHS yc and solution xc
   std::unique_ptr<cudssMatrix_t> Ac;
   mutable cudssMatrix_t xc, yc;
 
   // common for all cuDSS solver instances.
   // cuDSS object holds the cuDSS library context
   cudssHandle_t handle;
 
   // cuDSS object stores configuration settings for the solver
   mutable cudssConfig_t solverConfig;
   // cuDSS object holds internal data
   mutable cudssData_t solverData;
 
   /// Method for configuring storage for distributed/centralized RHS and
   /// solution
   void SetNumRHS(int nrhs_) const;
 
#ifdef MFEM_USE_MPI
   /**
    * @brief Set the HypreParMatrix object
    *
    * @param op HypreParMatrix object
    *
    * @note This method is called inside SetOperator
    */
   void SetMatrix(const HypreParMatrix &op);
#endif
 
   /**
    * @brief Set the SparseMatrix object
    *
    * @param op SparseMatrix object
    *
    * @note This method is called inside SetOperator
    */
   void SetMatrix(const SparseMatrix &op);
 
   /**
    * @brief Set the matrix values for cuDSS
    *
    * @param csr_offsets Row offsets of the CSR matrix
    * @param csr_columns Column indices of the CSR matrix
    * @param csr_values Non-zero values of the CSR matrix
    *
    * @note This method is called inside SetMatrix.
   */
   void SetMatrixCuDSS(int* csr_offsets, int* csr_columns, real_t* csr_values);
 
   /// Method for initializing the cuDSS library and creating the cuDSS handle
   void InitCuDSS();
};
 
} // namespace mfem
 
#endif // MFEM_USE_CUDSS
#endif // MFEM_CUDSS