Opm::gpuistl::detail Namespace Reference

Namespaces

namespace  DILU
 
namespace  ILU0
 
namespace  JAC
 

Classes

class  CuBlasHandle
 The CuBlasHandle class provides a singleton for the simulator universal cuBlasHandle. More...
 
class  CuSparseHandle
 The CuSparseHandle class provides a singleton for the simulator universal cuSparseHandle. More...
 
class  CuSparseResource
 The CuSparseResource class wraps a CuSparse resource in a proper RAII pattern. More...
 
class  has_communication
 The has_communication class checks if the type has the member function getCommunication. More...
 
class  has_should_call_post
 The has_should_call_post class detects the presence of the method shouldCallPost. More...
 
class  has_should_call_pre
 The has_should_call_pre class detects the presence of the method shouldCallPre. More...
 
class  is_a_well_operator
 The is_a_well_operator class tries to guess if the operator is a well type operator. More...
 

Typedefs

using GpuSparseMatrixDescription = CuSparseResource< cusparseMatDescr_t >
 
using GpuSparseMatrixDescriptionPtr = std::shared_ptr< CuSparseResource< cusparseMatDescr_t > >
 

Functions

template<typename func >
int tuneThreadBlockSize (func &f, std::string descriptionOfFunction)
 Function that tests the best thread block size, assumes the provided function depends on threadblock-size. More...
 
std::vector< int > createReorderedToNatural (const Opm::SparseTable< size_t > &levelSets)
 
std::vector< int > createNaturalToReordered (const Opm::SparseTable< size_t > &levelSets)
 
template<class M , class field_type , class GPUM >
std::unique_ptr< GPUM > createReorderedMatrix (const M &naturalMatrix, const std::vector< int > &reorderedToNatural)
 
template<class M , class field_type , class GPUM >
std::tuple< std::unique_ptr< GPUM >, std::unique_ptr< GPUM > > extractLowerAndUpperMatrices (const M &naturalMatrix, const std::vector< int > &reorderedToNatural)
 
std::string getCublasErrorMessage (cublasStatus_t error, const std::string_view &expression, const std::string_view &filename, const std::string_view &functionName, size_t lineNumber)
 getCublasErrorMessage generates the error message to display for a given error. More...
 
void cublasSafeCall (cublasStatus_t error, const std::string_view &expression, const std::string_view &filename, const std::string_view &functionName, size_t lineNumber)
 cublasSafeCall checks the return type of the CUBLAS expression (function call) and throws an exception if it does not equal CUBLAS_STATUS_SUCCESS. More...
 
cublasStatus_t cublasWarnIfError (cublasStatus_t error, const std::string_view &expression, const std::string_view &filename, const std::string_view &functionName, size_t lineNumber)
 cublasWarnIfError checks the return type of the CUBLAS expression (function call) and issues a warning if it does not equal CUBLAS_STATUS_SUCCESS. More...
 
cublasStatus_t cublasScal (cublasHandle_t handle, int n, const double *alpha, double *x, int incx)
 
cublasStatus_t cublasScal (cublasHandle_t handle, int n, const float *alpha, float *x, int incx)
 
cublasStatus_t cublasScal (cublasHandle_t handle, int n, const int *alpha, int *x, int incx)
 
cublasStatus_t cublasAxpy (cublasHandle_t handle, int n, const double *alpha, const double *x, int incx, double *y, int incy)
 
cublasStatus_t cublasAxpy (cublasHandle_t handle, int n, const float *alpha, const float *x, int incx, float *y, int incy)
 
cublasStatus_t cublasAxpy (cublasHandle_t handle, int n, const int *alpha, const int *x, int incx, int *y, int incy)
 
cublasStatus_t cublasDot (cublasHandle_t handle, int n, const double *x, int incx, const double *y, int incy, double *result)
 
cublasStatus_t cublasDot (cublasHandle_t handle, int n, const float *x, int incx, const float *y, int incy, float *result)
 
cublasStatus_t cublasDot (cublasHandle_t handle, int n, const int *x, int incx, const int *y, int incy, int *result)
 
cublasStatus_t cublasNrm2 (cublasHandle_t handle, int n, const double *x, int incx, double *result)
 
cublasStatus_t cublasNrm2 (cublasHandle_t handle, int n, const float *x, int incx, float *result)
 
cublasStatus_t cublasNrm2 (cublasHandle_t handle, int n, const int *x, int incx, int *result)
 
GpuSparseMatrixDescriptionPtr createMatrixDescription ()
 createMatrixDescription creates a default matrix description More...
 
GpuSparseMatrixDescriptionPtr createLowerDiagonalDescription ()
 createLowerDiagonalDescription creates a lower diagonal matrix description More...
 
GpuSparseMatrixDescriptionPtr createUpperDiagonalDescription ()
 createUpperDiagonalDescription creates an upper diagonal matrix description More...
 
std::string getCusparseErrorCodeToString (int code)
 getCusparseErrorCodeToString Converts an error code returned from a cusparse function a human readable string. More...
 
std::string getCusparseErrorMessage (cusparseStatus_t error, const std::string_view &expression, const std::string_view &filename, const std::string_view &functionName, size_t lineNumber)
 getCusparseErrorMessage generates the error message to display for a given error. More...
 
void cusparseSafeCall (cusparseStatus_t error, const std::string_view &expression, const std::string_view &filename, const std::string_view &functionName, size_t lineNumber)
 cusparseSafeCall checks the return type of the CUSPARSE expression (function call) and throws an exception if it does not equal CUSPARSE_STATUS_SUCCESS. More...
 
cusparseStatus_t cusparseWarnIfError (cusparseStatus_t error, const std::string_view &expression, const std::string_view &filename, const std::string_view &functionName, size_t lineNumber)
 cusparseWarnIfError checks the return type of the CUSPARSE expression (function call) and issues a warning if it does not equal CUSPARSE_STATUS_SUCCESS. More...
 
cusparseStatus_t cusparseBsrilu02_analysis (cusparseHandle_t handle, cusparseDirection_t dirA, int mb, int nnzb, const cusparseMatDescr_t descrA, double *bsrSortedVal, const int *bsrSortedRowPtr, const int *bsrSortedColInd, int blockDim, bsrilu02Info_t info, cusparseSolvePolicy_t policy, void *pBuffer)
 
cusparseStatus_t cusparseBsrsv2_analysis (cusparseHandle_t handle, cusparseDirection_t dirA, cusparseOperation_t transA, int mb, int nnzb, const cusparseMatDescr_t descrA, const double *bsrSortedValA, const int *bsrSortedRowPtrA, const int *bsrSortedColIndA, int blockDim, bsrsv2Info_t info, cusparseSolvePolicy_t policy, void *pBuffer)
 
cusparseStatus_t cusparseBsrsv2_analysis (cusparseHandle_t handle, cusparseDirection_t dirA, cusparseOperation_t transA, int mb, int nnzb, const cusparseMatDescr_t descrA, const float *bsrSortedValA, const int *bsrSortedRowPtrA, const int *bsrSortedColIndA, int blockDim, bsrsv2Info_t info, cusparseSolvePolicy_t policy, void *pBuffer)
 
cusparseStatus_t cusparseBsrilu02_analysis (cusparseHandle_t handle, cusparseDirection_t dirA, int mb, int nnzb, const cusparseMatDescr_t descrA, float *bsrSortedVal, const int *bsrSortedRowPtr, const int *bsrSortedColInd, int blockDim, bsrilu02Info_t info, cusparseSolvePolicy_t policy, void *pBuffer)
 
cusparseStatus_t cusparseBsrsv2_solve (cusparseHandle_t handle, cusparseDirection_t dirA, cusparseOperation_t transA, int mb, int nnzb, const double *alpha, const cusparseMatDescr_t descrA, const double *bsrSortedValA, const int *bsrSortedRowPtrA, const int *bsrSortedColIndA, int blockDim, bsrsv2Info_t info, const double *f, double *x, cusparseSolvePolicy_t policy, void *pBuffer)
 
cusparseStatus_t cusparseBsrsv2_solve (cusparseHandle_t handle, cusparseDirection_t dirA, cusparseOperation_t transA, int mb, int nnzb, const float *alpha, const cusparseMatDescr_t descrA, const float *bsrSortedValA, const int *bsrSortedRowPtrA, const int *bsrSortedColIndA, int blockDim, bsrsv2Info_t info, const float *f, float *x, cusparseSolvePolicy_t policy, void *pBuffer)
 
cusparseStatus_t cusparseBsrilu02_bufferSize (cusparseHandle_t handle, cusparseDirection_t dirA, int mb, int nnzb, const cusparseMatDescr_t descrA, double *bsrSortedVal, const int *bsrSortedRowPtr, const int *bsrSortedColInd, int blockDim, bsrilu02Info_t info, int *pBufferSizeInBytes)
 
cusparseStatus_t cusparseBsrilu02_bufferSize (cusparseHandle_t handle, cusparseDirection_t dirA, int mb, int nnzb, const cusparseMatDescr_t descrA, float *bsrSortedVal, const int *bsrSortedRowPtr, const int *bsrSortedColInd, int blockDim, bsrilu02Info_t info, int *pBufferSizeInBytes)
 
cusparseStatus_t cusparseBsrsv2_bufferSize (cusparseHandle_t handle, cusparseDirection_t dirA, cusparseOperation_t transA, int mb, int nnzb, const cusparseMatDescr_t descrA, double *bsrSortedValA, const int *bsrSortedRowPtrA, const int *bsrSortedColIndA, int blockDim, bsrsv2Info_t info, int *pBufferSizeInBytes)
 
cusparseStatus_t cusparseBsrsv2_bufferSize (cusparseHandle_t handle, cusparseDirection_t dirA, cusparseOperation_t transA, int mb, int nnzb, const cusparseMatDescr_t descrA, float *bsrSortedValA, const int *bsrSortedRowPtrA, const int *bsrSortedColIndA, int blockDim, bsrsv2Info_t info, int *pBufferSizeInBytes)
 
cusparseStatus_t cusparseBsrilu02 (cusparseHandle_t handle, cusparseDirection_t dirA, int mb, int nnzb, const cusparseMatDescr_t descrA, double *bsrSortedVal, const int *bsrSortedRowPtr, const int *bsrSortedColInd, int blockDim, bsrilu02Info_t info, cusparseSolvePolicy_t policy, void *pBuffer)
 
cusparseStatus_t cusparseBsrilu02 (cusparseHandle_t handle, cusparseDirection_t dirA, int mb, int nnzb, const cusparseMatDescr_t descrA, float *bsrSortedVal, const int *bsrSortedRowPtr, const int *bsrSortedColInd, int blockDim, bsrilu02Info_t info, cusparseSolvePolicy_t policy, void *pBuffer)
 
cusparseStatus_t cusparseBsrmv (cusparseHandle_t handle, cusparseDirection_t dirA, cusparseOperation_t transA, int mb, int nb, int nnzb, const double *alpha, const cusparseMatDescr_t descrA, const double *bsrSortedValA, const int *bsrSortedRowPtrA, const int *bsrSortedColIndA, int blockDim, const double *x, const double *beta, double *y)
 
cusparseStatus_t cusparseBsrmv (cusparseHandle_t handle, cusparseDirection_t dirA, cusparseOperation_t transA, int mb, int nb, int nnzb, const float *alpha, const cusparseMatDescr_t descrA, const float *bsrSortedValA, const int *bsrSortedRowPtrA, const int *bsrSortedColIndA, int blockDim, const float *x, const float *beta, float *y)
 
template<class Matrix >
const Matrix makeMatrixWithNonzeroDiagonal (const Matrix &matrix, const typename Matrix::field_type replacementValue=std::numeric_limits< typename Matrix::field_type >::epsilon())
 makeMatrixWithNonzeroDiagonal creates a new matrix with the zero diagonal elements (when viewed as a matrix of scalrars) set to replacementValue More...
 
std::string getCudaErrorMessage (cudaError_t error, const std::string_view &expression, const std::string_view &filename, const std::string_view &functionName, size_t lineNumber)
 getCudaErrorMessage generates the error message to display for a given error. More...
 
void cudaSafeCall (cudaError_t error, const std::string_view &expression, const std::string_view &filename, const std::string_view &functionName, size_t lineNumber)
 cudaSafeCall checks the return type of the GPU expression (function call) and throws an exception if it does not equal cudaSuccess. More...
 
void cudaWarnIfError (cudaError_t error, const std::string_view &expression, const std::string_view &filename, const std::string_view &functionName, size_t lineNumber)
 cudaWarnIfError checks the return type of the GPU expression (function call) and issues a warning if it does not equal cudaSuccess. More...
 
template<class T , int blocksize>
void copyMatDataToReordered (T *srcMatrix, int *srcRowIndices, T *dstMatrix, int *dstRowIndices, int *naturalToReordered, size_t numberOfRows, int threadBlockSize)
 Reorders the elements of a matrix by copying them from one matrix to another using a permutation list. More...
 
template<class T , int blocksize>
void copyMatDataToReorderedSplit (T *srcMatrix, int *srcRowIndices, int *srcColumnIndices, T *dstLowerMatrix, int *dstLowerRowIndices, T *dstUpperMatrix, int *dstUpperRowIndices, T *dstDiag, int *naturalToReordered, size_t numberOfRows, int threadBlockSize)
 Reorders the elements of a matrix by copying them from one matrix to a split matrix using a permutation list. More...
 
constexpr size_t getThreads (size_t numberOfRows)
 
size_t getBlocks (size_t numberOfRows)
 
template<class Kernel >
int getCudaRecomendedThreadBlockSize (Kernel k, int suggestedThrBlockSize=-1)
 
int getNumberOfBlocks (int wantedThreads, int threadBlockSize)
 
template<class PreconditionerType >
constexpr bool shouldCallPreconditionerPre ()
 Tests (compile time) if the preconditioner type needs to call pre() before a call to apply() More...
 
template<class PreconditionerType >
constexpr bool shouldCallPreconditionerPost ()
 Tests (compile time) if the preconditioner type needs to call post() after a call to apply(...) More...
 
int to_int (std::size_t s)
 to_int converts a (on most relevant platforms) 64 bits unsigned size_t to a signed 32 bits signed int More...
 
__host__ __device__ std::size_t to_size_t (int i)
 to_size_t converts a (on most relevant platforms) a 32 bit signed int to a 64 bits unsigned int More...
 
template<class T >
void setVectorValue (T *deviceData, size_t numberOfElements, const T &value)
 setVectorValue sets every element of deviceData to value More...
 
template<class T >
void setZeroAtIndexSet (T *deviceData, size_t numberOfElements, const int *indices)
 setZeroAtIndexSet sets deviceData to zero in the indices of contained in indices More...
 
template<class T >
innerProductAtIndices (cublasHandle_t cublasHandle, const T *deviceA, const T *deviceB, T *buffer, size_t numberOfElements, const int *indices)
 innerProductAtIndices computes the inner product between deviceA[indices] and deviceB[indices] More...
 
template<class T >
void prepareSendBuf (const T *deviceA, T *buffer, size_t numberOfElements, const int *indices)
 
template<class T >
void syncFromRecvBuf (T *deviceA, T *buffer, size_t numberOfElements, const int *indices)
 
template<class T >
void weightedDiagMV (const T *squareBlockVector, const size_t numberOfRows, const size_t blocksize, T relaxationFactor, const T *srcVec, T *dstVec)
 Compue the weighted matrix vector product where the matrix is diagonal, the diagonal is a vector, meaning we compute the Hadamard product. More...
 

Variables

const constexpr auto CUSPARSE_MATRIX_ORDER = CUSPARSE_DIRECTION_ROW
 

Detailed Description

Contains wrappers to make the CuBLAS library behave as a modern C++ library with function overlading.

In simple terms, this allows one to call say cublasScal on both double and single precisision, instead of calling cublasDscal and cublasSscal respectively.

Contains wrappers to make the CuSPARSE library behave as a modern C++ library with function overlading.

In simple terms, this allows one to call say cusparseBsrilu02_analysis on both double and single precisision, instead of calling cusparseDbsrilu02_analysis and cusparseDbsrilu02_analysis respectively.

Simple utility structs to test for the existence of functions in types.

Note that there are alternatives to this, see for instance https://stackoverflow.com/questions/257288/templated-check-for-the-existence-of-a-class-member-function , however, this is by far the cleanest approach for where this is going to be used for now.

TODO: Use the requires-keyword once C++20 becomes availble (https://en.cppreference.com/w/cpp/language/constraints ). With C++20 this file can be removed.

Provides various utilities for doing signed to unsigned conversion, unsigned to signed, 32 bits to 64 bits and 64 bits to 32 bits.

The main use case within cuistl is that the cusparse library requires signed int for all its size parameters, while Dune::BlockVector (and relatives) use unsigned size_t.

Typedef Documentation

◆ GpuSparseMatrixDescription

GpuSparseMatrixDescription holder. This is internal information needed for most calls to the CuSparse API.

◆ GpuSparseMatrixDescriptionPtr

using Opm::gpuistl::detail::GpuSparseMatrixDescriptionPtr = typedef std::shared_ptr<CuSparseResource<cusparseMatDescr_t> >

Pointer to GpuSparseMatrixDescription holder. This is internal information needed for most calls to the CuSparse API.

Function Documentation

◆ copyMatDataToReordered()

template<class T , int blocksize>
void Opm::gpuistl::detail::copyMatDataToReordered ( T *  srcMatrix,
int *  srcRowIndices,
T *  dstMatrix,
int *  dstRowIndices,
int *  naturalToReordered,
size_t  numberOfRows,
int  threadBlockSize 
)

Reorders the elements of a matrix by copying them from one matrix to another using a permutation list.

Parameters
srcMatrixThe source matrix we will copy data from
srcRowIndicesPointer to vector on GPU containing row indices for the source matrix compliant wiht bsr format
[out]dstMatrixThe destination matrix that we copy data to
dstRowIndicesPointer to vector on GPU containing riw indices for the destination matrix compliant wiht bsr format
naturalToReorderedPermuation list that converts indices in the src matrix to the indices in the dst matrix
numberOfRowsThe number of rows in the matrices

◆ copyMatDataToReorderedSplit()

template<class T , int blocksize>
void Opm::gpuistl::detail::copyMatDataToReorderedSplit ( T *  srcMatrix,
int *  srcRowIndices,
int *  srcColumnIndices,
T *  dstLowerMatrix,
int *  dstLowerRowIndices,
T *  dstUpperMatrix,
int *  dstUpperRowIndices,
T *  dstDiag,
int *  naturalToReordered,
size_t  numberOfRows,
int  threadBlockSize 
)

Reorders the elements of a matrix by copying them from one matrix to a split matrix using a permutation list.

Parameters
srcMatrixThe source matrix we will copy data from
srcRowIndicesPointer to vector on GPU containing row indices for the source matrix compliant wiht bsr format
[out]dstLowerMatrixThe destination of entries that originates from the strictly lower triangular matrix
dstRowIndicesPointer to vector on GPU containing rww indices for the destination lower matrix compliant wiht bsr format
[out]dstUpperMatrixThe destination of entries that originates from the strictly upper triangular matrix
dstRowIndicesPointer to vector on GPU containing riw indices for the destination upper matrix compliant wiht bsr format
[out]dstDiagThe destination buffer for the diagonal part of the matrix
naturalToReorderedPermuation list that converts indices in the src matrix to the indices in the dst matrix
numberOfRowsThe number of rows in the matrices

◆ createLowerDiagonalDescription()

GpuSparseMatrixDescriptionPtr Opm::gpuistl::detail::createLowerDiagonalDescription ( )
inline

createLowerDiagonalDescription creates a lower diagonal matrix description

Returns
a lower diagonal matrix description overlapped with options from Opm::gpuistl::detail::createMatrixDescription()
Note
This will assume it has a unit diagonal

References createMatrixDescription(), and OPM_CUSPARSE_SAFE_CALL.

◆ createMatrixDescription()

GpuSparseMatrixDescriptionPtr Opm::gpuistl::detail::createMatrixDescription ( )
inline

createMatrixDescription creates a default matrix description

Returns
a matrix description to a general sparse matrix with zero based indexing.

References OPM_CUSPARSE_SAFE_CALL.

Referenced by createLowerDiagonalDescription(), and createUpperDiagonalDescription().

◆ createNaturalToReordered()

std::vector< int > Opm::gpuistl::detail::createNaturalToReordered ( const Opm::SparseTable< size_t > &  levelSets)
inline

References to_size_t().

◆ createReorderedMatrix()

template<class M , class field_type , class GPUM >
std::unique_ptr< GPUM > Opm::gpuistl::detail::createReorderedMatrix ( const M &  naturalMatrix,
const std::vector< int > &  reorderedToNatural 
)
inline

◆ createReorderedToNatural()

std::vector< int > Opm::gpuistl::detail::createReorderedToNatural ( const Opm::SparseTable< size_t > &  levelSets)
inline

References to_size_t().

◆ createUpperDiagonalDescription()

GpuSparseMatrixDescriptionPtr Opm::gpuistl::detail::createUpperDiagonalDescription ( )
inline

createUpperDiagonalDescription creates an upper diagonal matrix description

Returns
an upper diagonal matrix description overlapped with options from Opm::gpuistl::detail::createMatrixDescription()
Note
This will assume it has a non-unit diagonal.

References createMatrixDescription(), and OPM_CUSPARSE_SAFE_CALL.

◆ cublasAxpy() [1/3]

cublasStatus_t Opm::gpuistl::detail::cublasAxpy ( cublasHandle_t  handle,
int  n,
const double *  alpha,
const double *  x,
int  incx,
double *  y,
int  incy 
)
inline

◆ cublasAxpy() [2/3]

cublasStatus_t Opm::gpuistl::detail::cublasAxpy ( cublasHandle_t  handle,
int  n,
const float *  alpha,
const float *  x,
int  incx,
float *  y,
int  incy 
)
inline

◆ cublasAxpy() [3/3]

cublasStatus_t Opm::gpuistl::detail::cublasAxpy ( cublasHandle_t  handle,
int  n,
const int *  alpha,
const int *  x,
int  incx,
int *  y,
int  incy 
)
inline

◆ cublasDot() [1/3]

cublasStatus_t Opm::gpuistl::detail::cublasDot ( cublasHandle_t  handle,
int  n,
const double *  x,
int  incx,
const double *  y,
int  incy,
double *  result 
)
inline

◆ cublasDot() [2/3]

cublasStatus_t Opm::gpuistl::detail::cublasDot ( cublasHandle_t  handle,
int  n,
const float *  x,
int  incx,
const float *  y,
int  incy,
float *  result 
)
inline

◆ cublasDot() [3/3]

cublasStatus_t Opm::gpuistl::detail::cublasDot ( cublasHandle_t  handle,
int  n,
const int *  x,
int  incx,
const int *  y,
int  incy,
int *  result 
)
inline

◆ cublasNrm2() [1/3]

cublasStatus_t Opm::gpuistl::detail::cublasNrm2 ( cublasHandle_t  handle,
int  n,
const double *  x,
int  incx,
double *  result 
)
inline

◆ cublasNrm2() [2/3]

cublasStatus_t Opm::gpuistl::detail::cublasNrm2 ( cublasHandle_t  handle,
int  n,
const float *  x,
int  incx,
float *  result 
)
inline

◆ cublasNrm2() [3/3]

cublasStatus_t Opm::gpuistl::detail::cublasNrm2 ( cublasHandle_t  handle,
int  n,
const int *  x,
int  incx,
int *  result 
)
inline

◆ cublasSafeCall()

void Opm::gpuistl::detail::cublasSafeCall ( cublasStatus_t  error,
const std::string_view &  expression,
const std::string_view &  filename,
const std::string_view &  functionName,
size_t  lineNumber 
)
inline

cublasSafeCall checks the return type of the CUBLAS expression (function call) and throws an exception if it does not equal CUBLAS_STATUS_SUCCESS.

Parameters
errorthe error code from cublas
expressionthe expresison (say "cublasCreate(&handle)")
filenamethe code file the error occured in (typically FILE)
functionNamename of the function the error occured in (typically func)
lineNumberthe line number the error occured in (typically LINE)

Example usage:

#include <cublas_v2.h>
void some_function() {
cublasHandle_t cublasHandle;
cudaSafeCall(cublasCreate(&cublasHandle), "cublasCreate(&cublasHandle)", __FILE__, __func__, __LINE__);
}
void cudaSafeCall(cudaError_t error, const std::string_view &expression, const std::string_view &filename, const std::string_view &functionName, size_t lineNumber)
cudaSafeCall checks the return type of the GPU expression (function call) and throws an exception if ...
Definition: gpu_safe_call.hpp:82
Note
It is probably easier to use the macro OPM_CUBLAS_SAFE_CALL

References getCublasErrorMessage().

◆ cublasScal() [1/3]

cublasStatus_t Opm::gpuistl::detail::cublasScal ( cublasHandle_t  handle,
int  n,
const double *  alpha,
double *  x,
int  incx 
)
inline

◆ cublasScal() [2/3]

cublasStatus_t Opm::gpuistl::detail::cublasScal ( cublasHandle_t  handle,
int  n,
const float *  alpha,
float *  x,
int  incx 
)
inline

◆ cublasScal() [3/3]

cublasStatus_t Opm::gpuistl::detail::cublasScal ( cublasHandle_t  handle,
int  n,
const int *  alpha,
int *  x,
int  incx 
)
inline

◆ cublasWarnIfError()

cublasStatus_t Opm::gpuistl::detail::cublasWarnIfError ( cublasStatus_t  error,
const std::string_view &  expression,
const std::string_view &  filename,
const std::string_view &  functionName,
size_t  lineNumber 
)
inline

cublasWarnIfError checks the return type of the CUBLAS expression (function call) and issues a warning if it does not equal CUBLAS_STATUS_SUCCESS.

Parameters
errorthe error code from cublas
expressionthe expresison (say "cublasCreate(&handle)")
filenamethe code file the error occured in (typically FILE)
functionNamename of the function the error occured in (typically func)
lineNumberthe line number the error occured in (typically LINE)
Returns
the error sent in (for convenience).

Example usage:

#include <cublas_v2.h>
void some_function() {
cublasHandle_t cublasHandle;
cublasWarnIfError(cublasCreate(&cublasHandle), "cublasCreate(&cublasHandle)", __FILE__, __func__, __LINE__);
}
cublasStatus_t cublasWarnIfError(cublasStatus_t error, const std::string_view &expression, const std::string_view &filename, const std::string_view &functionName, size_t lineNumber)
cublasWarnIfError checks the return type of the CUBLAS expression (function call) and issues a warnin...
Definition: cublas_safe_call.hpp:165
Note
It is probably easier to use the macro OPM_CUBLAS_WARN_IF_ERROR
Prefer the cublasSafeCall/OPM_CUBLAS_SAFE_CALL counterpart unless you really don't want to throw an exception.

References getCublasErrorMessage().

◆ cudaSafeCall()

void Opm::gpuistl::detail::cudaSafeCall ( cudaError_t  error,
const std::string_view &  expression,
const std::string_view &  filename,
const std::string_view &  functionName,
size_t  lineNumber 
)
inline

cudaSafeCall checks the return type of the GPU expression (function call) and throws an exception if it does not equal cudaSuccess.

Example usage:

#include <cuda_runtime.h>
void some_function() {
void* somePointer;
cudaSafeCall(cudaMalloc(&somePointer, 1), "cudaMalloc(&somePointer, 1)", __FILE__, __func__, __LINE__);
}
Note
It is probably easier to use the macro OPM_GPU_SAFE_CALL

References getCudaErrorMessage().

◆ cudaWarnIfError()

void Opm::gpuistl::detail::cudaWarnIfError ( cudaError_t  error,
const std::string_view &  expression,
const std::string_view &  filename,
const std::string_view &  functionName,
size_t  lineNumber 
)
inline

cudaWarnIfError checks the return type of the GPU expression (function call) and issues a warning if it does not equal cudaSuccess.

Parameters
errorthe error code from cublas
expressionthe expresison (say "cudaMalloc(&pointer, 1)")
filenamethe code file the error occured in (typically FILE)
functionNamename of the function the error occured in (typically func)
lineNumberthe line number the error occured in (typically LINE)

Example usage:

#include <cuda_runtime.h>
void some_function() {
void* somePointer;
cudaWarnIfError(cudaMalloc(&somePointer, 1), "cudaMalloc(&somePointer, 1)", __FILE__, __func__, __LINE__);
}
void cudaWarnIfError(cudaError_t error, const std::string_view &expression, const std::string_view &filename, const std::string_view &functionName, size_t lineNumber)
cudaWarnIfError checks the return type of the GPU expression (function call) and issues a warning if ...
Definition: gpu_safe_call.hpp:121
Note
It is probably easier to use the macro OPM_GPU_WARN_IF_ERROR
Prefer the cudaSafeCall/OPM_GPU_SAFE_CALL counterpart unless you really don't want to throw an exception.

References getCudaErrorMessage().

◆ cusparseBsrilu02() [1/2]

cusparseStatus_t Opm::gpuistl::detail::cusparseBsrilu02 ( cusparseHandle_t  handle,
cusparseDirection_t  dirA,
int  mb,
int  nnzb,
const cusparseMatDescr_t  descrA,
double *  bsrSortedVal,
const int *  bsrSortedRowPtr,
const int *  bsrSortedColInd,
int  blockDim,
bsrilu02Info_t  info,
cusparseSolvePolicy_t  policy,
void *  pBuffer 
)
inline

◆ cusparseBsrilu02() [2/2]

cusparseStatus_t Opm::gpuistl::detail::cusparseBsrilu02 ( cusparseHandle_t  handle,
cusparseDirection_t  dirA,
int  mb,
int  nnzb,
const cusparseMatDescr_t  descrA,
float *  bsrSortedVal,
const int *  bsrSortedRowPtr,
const int *  bsrSortedColInd,
int  blockDim,
bsrilu02Info_t  info,
cusparseSolvePolicy_t  policy,
void *  pBuffer 
)
inline

◆ cusparseBsrilu02_analysis() [1/2]

cusparseStatus_t Opm::gpuistl::detail::cusparseBsrilu02_analysis ( cusparseHandle_t  handle,
cusparseDirection_t  dirA,
int  mb,
int  nnzb,
const cusparseMatDescr_t  descrA,
double *  bsrSortedVal,
const int *  bsrSortedRowPtr,
const int *  bsrSortedColInd,
int  blockDim,
bsrilu02Info_t  info,
cusparseSolvePolicy_t  policy,
void *  pBuffer 
)
inline

◆ cusparseBsrilu02_analysis() [2/2]

cusparseStatus_t Opm::gpuistl::detail::cusparseBsrilu02_analysis ( cusparseHandle_t  handle,
cusparseDirection_t  dirA,
int  mb,
int  nnzb,
const cusparseMatDescr_t  descrA,
float *  bsrSortedVal,
const int *  bsrSortedRowPtr,
const int *  bsrSortedColInd,
int  blockDim,
bsrilu02Info_t  info,
cusparseSolvePolicy_t  policy,
void *  pBuffer 
)
inline

◆ cusparseBsrilu02_bufferSize() [1/2]

cusparseStatus_t Opm::gpuistl::detail::cusparseBsrilu02_bufferSize ( cusparseHandle_t  handle,
cusparseDirection_t  dirA,
int  mb,
int  nnzb,
const cusparseMatDescr_t  descrA,
double *  bsrSortedVal,
const int *  bsrSortedRowPtr,
const int *  bsrSortedColInd,
int  blockDim,
bsrilu02Info_t  info,
int *  pBufferSizeInBytes 
)
inline

◆ cusparseBsrilu02_bufferSize() [2/2]

cusparseStatus_t Opm::gpuistl::detail::cusparseBsrilu02_bufferSize ( cusparseHandle_t  handle,
cusparseDirection_t  dirA,
int  mb,
int  nnzb,
const cusparseMatDescr_t  descrA,
float *  bsrSortedVal,
const int *  bsrSortedRowPtr,
const int *  bsrSortedColInd,
int  blockDim,
bsrilu02Info_t  info,
int *  pBufferSizeInBytes 
)
inline

◆ cusparseBsrmv() [1/2]

cusparseStatus_t Opm::gpuistl::detail::cusparseBsrmv ( cusparseHandle_t  handle,
cusparseDirection_t  dirA,
cusparseOperation_t  transA,
int  mb,
int  nb,
int  nnzb,
const double *  alpha,
const cusparseMatDescr_t  descrA,
const double *  bsrSortedValA,
const int *  bsrSortedRowPtrA,
const int *  bsrSortedColIndA,
int  blockDim,
const double *  x,
const double *  beta,
double *  y 
)
inline

◆ cusparseBsrmv() [2/2]

cusparseStatus_t Opm::gpuistl::detail::cusparseBsrmv ( cusparseHandle_t  handle,
cusparseDirection_t  dirA,
cusparseOperation_t  transA,
int  mb,
int  nb,
int  nnzb,
const float *  alpha,
const cusparseMatDescr_t  descrA,
const float *  bsrSortedValA,
const int *  bsrSortedRowPtrA,
const int *  bsrSortedColIndA,
int  blockDim,
const float *  x,
const float *  beta,
float *  y 
)
inline

◆ cusparseBsrsv2_analysis() [1/2]

cusparseStatus_t Opm::gpuistl::detail::cusparseBsrsv2_analysis ( cusparseHandle_t  handle,
cusparseDirection_t  dirA,
cusparseOperation_t  transA,
int  mb,
int  nnzb,
const cusparseMatDescr_t  descrA,
const double *  bsrSortedValA,
const int *  bsrSortedRowPtrA,
const int *  bsrSortedColIndA,
int  blockDim,
bsrsv2Info_t  info,
cusparseSolvePolicy_t  policy,
void *  pBuffer 
)
inline

◆ cusparseBsrsv2_analysis() [2/2]

cusparseStatus_t Opm::gpuistl::detail::cusparseBsrsv2_analysis ( cusparseHandle_t  handle,
cusparseDirection_t  dirA,
cusparseOperation_t  transA,
int  mb,
int  nnzb,
const cusparseMatDescr_t  descrA,
const float *  bsrSortedValA,
const int *  bsrSortedRowPtrA,
const int *  bsrSortedColIndA,
int  blockDim,
bsrsv2Info_t  info,
cusparseSolvePolicy_t  policy,
void *  pBuffer 
)
inline

◆ cusparseBsrsv2_bufferSize() [1/2]

cusparseStatus_t Opm::gpuistl::detail::cusparseBsrsv2_bufferSize ( cusparseHandle_t  handle,
cusparseDirection_t  dirA,
cusparseOperation_t  transA,
int  mb,
int  nnzb,
const cusparseMatDescr_t  descrA,
double *  bsrSortedValA,
const int *  bsrSortedRowPtrA,
const int *  bsrSortedColIndA,
int  blockDim,
bsrsv2Info_t  info,
int *  pBufferSizeInBytes 
)
inline

◆ cusparseBsrsv2_bufferSize() [2/2]

cusparseStatus_t Opm::gpuistl::detail::cusparseBsrsv2_bufferSize ( cusparseHandle_t  handle,
cusparseDirection_t  dirA,
cusparseOperation_t  transA,
int  mb,
int  nnzb,
const cusparseMatDescr_t  descrA,
float *  bsrSortedValA,
const int *  bsrSortedRowPtrA,
const int *  bsrSortedColIndA,
int  blockDim,
bsrsv2Info_t  info,
int *  pBufferSizeInBytes 
)
inline

◆ cusparseBsrsv2_solve() [1/2]

cusparseStatus_t Opm::gpuistl::detail::cusparseBsrsv2_solve ( cusparseHandle_t  handle,
cusparseDirection_t  dirA,
cusparseOperation_t  transA,
int  mb,
int  nnzb,
const double *  alpha,
const cusparseMatDescr_t  descrA,
const double *  bsrSortedValA,
const int *  bsrSortedRowPtrA,
const int *  bsrSortedColIndA,
int  blockDim,
bsrsv2Info_t  info,
const double *  f,
double *  x,
cusparseSolvePolicy_t  policy,
void *  pBuffer 
)
inline

◆ cusparseBsrsv2_solve() [2/2]

cusparseStatus_t Opm::gpuistl::detail::cusparseBsrsv2_solve ( cusparseHandle_t  handle,
cusparseDirection_t  dirA,
cusparseOperation_t  transA,
int  mb,
int  nnzb,
const float *  alpha,
const cusparseMatDescr_t  descrA,
const float *  bsrSortedValA,
const int *  bsrSortedRowPtrA,
const int *  bsrSortedColIndA,
int  blockDim,
bsrsv2Info_t  info,
const float *  f,
float *  x,
cusparseSolvePolicy_t  policy,
void *  pBuffer 
)
inline

◆ cusparseSafeCall()

void Opm::gpuistl::detail::cusparseSafeCall ( cusparseStatus_t  error,
const std::string_view &  expression,
const std::string_view &  filename,
const std::string_view &  functionName,
size_t  lineNumber 
)
inline

cusparseSafeCall checks the return type of the CUSPARSE expression (function call) and throws an exception if it does not equal CUSPARSE_STATUS_SUCCESS.

Example usage:

#include <cublas_v2.h>
void some_function() {
cusparseHandle_t cusparseHandle;
cusparseSafeCall(cusparseCreate(&cusparseHandle), "cusparseCreate(&cusparseHandle)", __FILE__, __func__,
__LINE__);
}
void cusparseSafeCall(cusparseStatus_t error, const std::string_view &expression, const std::string_view &filename, const std::string_view &functionName, size_t lineNumber)
cusparseSafeCall checks the return type of the CUSPARSE expression (function call) and throws an exce...
Definition: cusparse_safe_call.hpp:111
Note
It is probably easier to use the macro OPM_CUBLAS_SAFE_CALL

References getCusparseErrorMessage().

◆ cusparseWarnIfError()

cusparseStatus_t Opm::gpuistl::detail::cusparseWarnIfError ( cusparseStatus_t  error,
const std::string_view &  expression,
const std::string_view &  filename,
const std::string_view &  functionName,
size_t  lineNumber 
)
inline

cusparseWarnIfError checks the return type of the CUSPARSE expression (function call) and issues a warning if it does not equal CUSPARSE_STATUS_SUCCESS.

Parameters
errorthe error code from cublas
expressionthe expresison (say "cublasCreate(&handle)")
filenamethe code file the error occured in (typically FILE)
functionNamename of the function the error occured in (typically func)
lineNumberthe line number the error occured in (typically LINE)
Returns
the error sent in (for convenience).

Example usage:

#include <cublas_v2.h>
void some_function() {
cusparseHandle_t cusparseHandle;
cusparseWarnIfError(cusparseCreate(&cusparseHandle), "cusparseCreate(&cusparseHandle)", __FILE__, __func__,
__LINE__);
}
cusparseStatus_t cusparseWarnIfError(cusparseStatus_t error, const std::string_view &expression, const std::string_view &filename, const std::string_view &functionName, size_t lineNumber)
cusparseWarnIfError checks the return type of the CUSPARSE expression (function call) and issues a wa...
Definition: cusparse_safe_call.hpp:152
Note
It is probably easier to use the macro OPM_CUSPARSE_WARN_IF_ERROR
Prefer the cusparseSafeCall/OPM_CUSPARSE_SAFE_CALL counterpart unless you really don't want to throw an exception.

References getCusparseErrorMessage().

◆ extractLowerAndUpperMatrices()

template<class M , class field_type , class GPUM >
std::tuple< std::unique_ptr< GPUM >, std::unique_ptr< GPUM > > Opm::gpuistl::detail::extractLowerAndUpperMatrices ( const M &  naturalMatrix,
const std::vector< int > &  reorderedToNatural 
)
inline

◆ getBlocks()

size_t Opm::gpuistl::detail::getBlocks ( size_t  numberOfRows)
inline

References getThreads().

◆ getCublasErrorMessage()

std::string Opm::gpuistl::detail::getCublasErrorMessage ( cublasStatus_t  error,
const std::string_view &  expression,
const std::string_view &  filename,
const std::string_view &  functionName,
size_t  lineNumber 
)
inline

getCublasErrorMessage generates the error message to display for a given error.

Parameters
errorthe error code from cublas
expressionthe expresison (say "cublasCreate(&handle)")
filenamethe code file the error occured in (typically FILE)
functionNamename of the function the error occured in (typically func)
lineNumberthe line number the error occured in (typically LINE)
Returns
An error message to be displayed.
Note
This function is mostly for internal use.

Referenced by cublasSafeCall(), and cublasWarnIfError().

◆ getCudaErrorMessage()

std::string Opm::gpuistl::detail::getCudaErrorMessage ( cudaError_t  error,
const std::string_view &  expression,
const std::string_view &  filename,
const std::string_view &  functionName,
size_t  lineNumber 
)
inline

getCudaErrorMessage generates the error message to display for a given error.

Parameters
errorthe error code from cublas
expressionthe expresison (say "cudaMalloc(&pointer, 1)")
filenamethe code file the error occured in (typically FILE)
functionNamename of the function the error occured in (typically func)
lineNumberthe line number the error occured in (typically LINE)
Returns
An error message to be displayed.
Note
This function is mostly for internal use.

Referenced by cudaSafeCall(), and cudaWarnIfError().

◆ getCudaRecomendedThreadBlockSize()

template<class Kernel >
int Opm::gpuistl::detail::getCudaRecomendedThreadBlockSize ( Kernel  k,
int  suggestedThrBlockSize = -1 
)
inline

References OPM_GPU_SAFE_CALL.

◆ getCusparseErrorCodeToString()

std::string Opm::gpuistl::detail::getCusparseErrorCodeToString ( int  code)
inline

getCusparseErrorCodeToString Converts an error code returned from a cusparse function a human readable string.

Parameters
codean error code from a cusparse routine
Returns
a human readable string.

References CHECK_CUSPARSE_ERROR_TYPE.

Referenced by getCusparseErrorMessage().

◆ getCusparseErrorMessage()

std::string Opm::gpuistl::detail::getCusparseErrorMessage ( cusparseStatus_t  error,
const std::string_view &  expression,
const std::string_view &  filename,
const std::string_view &  functionName,
size_t  lineNumber 
)
inline

getCusparseErrorMessage generates the error message to display for a given error.

Parameters
errorthe error code from cublas
expressionthe expresison (say "cusparseCreate(&handle)")
filenamethe code file the error occured in (typically FILE)
functionNamename of the function the error occured in (typically func)
lineNumberthe line number the error occured in (typically LINE)
Returns
An error message to be displayed.
Note
This function is mostly for internal use.

References getCusparseErrorCodeToString().

Referenced by cusparseSafeCall(), and cusparseWarnIfError().

◆ getNumberOfBlocks()

int Opm::gpuistl::detail::getNumberOfBlocks ( int  wantedThreads,
int  threadBlockSize 
)
inline

◆ getThreads()

constexpr size_t Opm::gpuistl::detail::getThreads ( size_t  numberOfRows)
inlineconstexpr

Referenced by getBlocks().

◆ innerProductAtIndices()

template<class T >
T Opm::gpuistl::detail::innerProductAtIndices ( cublasHandle_t  cublasHandle,
const T *  deviceA,
const T *  deviceB,
T *  buffer,
size_t  numberOfElements,
const int *  indices 
)

innerProductAtIndices computes the inner product between deviceA[indices] and deviceB[indices]

Parameters
cublasHandlea valid (initialized) cublas handle
deviceAdata A (device memory)
deviceBdata B (device memory)
buffera buffer with number of elements equal to numberOfElements (device memory)
numberOfElementsnumber of indices
indicesthe indices to compute the inner product over (device memory)
Returns
the result of the inner product
Note
This is equivalent to projecting the vectors to the indices contained in indices, then doing the inner product of those projected vectors.

◆ makeMatrixWithNonzeroDiagonal()

template<class Matrix >
const Matrix Opm::gpuistl::detail::makeMatrixWithNonzeroDiagonal ( const Matrix &  matrix,
const typename Matrix::field_type  replacementValue = std::numeric_limits<typename Matrix::field_type>::epsilon() 
)

makeMatrixWithNonzeroDiagonal creates a new matrix with the zero diagonal elements (when viewed as a matrix of scalrars) set to replacementValue

Parameters
matrixthe matrix to replace
replacementValuethe value to set in the diagonal elements that are zero
Returns
a new matrix with non non-zero diagonal elements.
Note
This modification is needed for the CuSparse implementation of ILU0. While the the algorithm operates on block matrices, it still requires that the scalar matrix has no zero diagonal elements.

◆ prepareSendBuf()

template<class T >
void Opm::gpuistl::detail::prepareSendBuf ( const T *  deviceA,
T *  buffer,
size_t  numberOfElements,
const int *  indices 
)

◆ setVectorValue()

template<class T >
void Opm::gpuistl::detail::setVectorValue ( T *  deviceData,
size_t  numberOfElements,
const T &  value 
)

setVectorValue sets every element of deviceData to value

Parameters
deviceDatapointer to GPU memory
numberOfElementsnumber of elements to set to value
valuethe value to use

◆ setZeroAtIndexSet()

template<class T >
void Opm::gpuistl::detail::setZeroAtIndexSet ( T *  deviceData,
size_t  numberOfElements,
const int *  indices 
)

setZeroAtIndexSet sets deviceData to zero in the indices of contained in indices

Parameters
deviceDatathe data to operate on (device memory)
numberOfElementsnumber of indices
indicesthe indices to use (device memory)

◆ shouldCallPreconditionerPost()

template<class PreconditionerType >
constexpr bool Opm::gpuistl::detail::shouldCallPreconditionerPost ( )
constexpr

Tests (compile time) if the preconditioner type needs to call post() after a call to apply(...)

Note
This is mostly used to avoid unneeded copying back and front to the GPU, as well as avoiding communication.

◆ shouldCallPreconditionerPre()

template<class PreconditionerType >
constexpr bool Opm::gpuistl::detail::shouldCallPreconditionerPre ( )
constexpr

Tests (compile time) if the preconditioner type needs to call pre() before a call to apply()

Note
This is mostly used to avoid unneeded copying back and front to the GPU, as well as avoiding communication.

◆ syncFromRecvBuf()

template<class T >
void Opm::gpuistl::detail::syncFromRecvBuf ( T *  deviceA,
T *  buffer,
size_t  numberOfElements,
const int *  indices 
)

◆ to_int()

int Opm::gpuistl::detail::to_int ( std::size_t  s)
inline

to_int converts a (on most relevant platforms) 64 bits unsigned size_t to a signed 32 bits signed int

Parameters
sthe unsigned integer
Exceptions
std::invalid_argumentexception if s is out of range for an int
Returns
converted s to int if s is within the range of int

Referenced by Opm::gpuistl::GPUAwareMPISender< field_type, block_size, OwnerOverlapCopyCommunicationType >::copyOwnerToAll().

◆ to_size_t()

__host__ __device__ std::size_t Opm::gpuistl::detail::to_size_t ( int  i)
inline

to_size_t converts a (on most relevant platforms) a 32 bit signed int to a 64 bits unsigned int

Parameters
ithe signed integer
Returns
converted i to size_t if it is a non-negative integer.
Exceptions
std::invalid_argumentif i is negative.

Referenced by Opm::gpuistl::GpuSparseMatrix< T >::blockSize(), createNaturalToReordered(), createReorderedToNatural(), Opm::gpuistl::GpuSparseMatrix< T >::dim(), Opm::gpuistl::GpuSparseMatrix< T >::N(), and Opm::gpuistl::GpuSparseMatrix< T >::nonzeroes().

◆ tuneThreadBlockSize()

template<typename func >
int Opm::gpuistl::detail::tuneThreadBlockSize ( func &  f,
std::string  descriptionOfFunction 
)

Function that tests the best thread block size, assumes the provided function depends on threadblock-size.

Template Parameters
Thetype of the function to tune
Parameters
fthe function to tune, which takes the thread block size as the input

References OPM_GPU_SAFE_CALL.

◆ weightedDiagMV()

template<class T >
void Opm::gpuistl::detail::weightedDiagMV ( const T *  squareBlockVector,
const size_t  numberOfRows,
const size_t  blocksize,
relaxationFactor,
const T *  srcVec,
T *  dstVec 
)

Compue the weighted matrix vector product where the matrix is diagonal, the diagonal is a vector, meaning we compute the Hadamard product.

Parameters
squareBlockVectorA GpuVector whose elements are NxN matrix blocks
numberOfRowsThe number of rows in the vector
blocksizeThe sidelength of the square block elements in the vector
src_vecA pointer to the data of the GpuVector we multiply the blockvector with
[out]dst_vecA pointer to the data of the GpuVector we store the result in
Note
This is implemented as a faster way to multiply a diagonal matrix with a blockvector. We need only store the diagonal of the matrix and use this product.

Variable Documentation

◆ CUSPARSE_MATRIX_ORDER

const constexpr auto Opm::gpuistl::detail::CUSPARSE_MATRIX_ORDER = CUSPARSE_DIRECTION_ROW
constexpr