Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
Public Types | Private Attributes | List of all members
CUDAWrappers::SparseMatrix< Number > Class Template Reference

#include <deal.II/lac/cuda_sparse_matrix.h>

Inheritance diagram for CUDAWrappers::SparseMatrix< Number >:
[legend]

Public Types

using size_type = int
 
using value_type = Number
 
using real_type = Number
 

Public Member Functions

Constructors and initialization
 SparseMatrix ()
 
 SparseMatrix (Utilities::CUDA::Handle &handle, const ::SparseMatrix< Number > &sparse_matrix_host)
 
 SparseMatrix (CUDAWrappers::SparseMatrix< Number > &&)
 
 SparseMatrix (const CUDAWrappers::SparseMatrix< Number > &)=delete
 
 ~SparseMatrix ()
 
SparseMatrixoperator= (CUDAWrappers::SparseMatrix< Number > &&)
 
SparseMatrixoperator= (const CUDAWrappers::SparseMatrix< Number > &)=delete
 
void reinit (Utilities::CUDA::Handle &handle, const ::SparseMatrix< Number > &sparse_matrix_host)
 
Information on the matrix
size_type m () const
 
size_type n () const
 
std::size_t n_nonzero_elements () const
 
template<class StreamType >
void print (StreamType &out, const bool across=false, const bool diagonal_first=true) const
 
void print_formatted (std::ostream &out, const unsigned int precision=3, const bool scientific=true, const unsigned int width=0, const char *zero_string=" ", const double denominator=1.) const
 
Modifying entries
SparseMatrixoperator*= (const Number factor)
 
SparseMatrixoperator/= (const Number factor)
 
Multiplications
void vmult (LinearAlgebra::CUDAWrappers::Vector< Number > &dst, const LinearAlgebra::CUDAWrappers::Vector< Number > &src) const
 
void Tvmult (LinearAlgebra::CUDAWrappers::Vector< Number > &dst, const LinearAlgebra::CUDAWrappers::Vector< Number > &src) const
 
void vmult_add (LinearAlgebra::CUDAWrappers::Vector< Number > &dst, const LinearAlgebra::CUDAWrappers::Vector< Number > &src) const
 
void Tvmult_add (LinearAlgebra::CUDAWrappers::Vector< Number > &dst, const LinearAlgebra::CUDAWrappers::Vector< Number > &src) const
 
Number matrix_norm_square (const LinearAlgebra::CUDAWrappers::Vector< Number > &v) const
 
Number matrix_scalar_product (const LinearAlgebra::CUDAWrappers::Vector< Number > &u, const LinearAlgebra::CUDAWrappers::Vector< Number > &v) const
 
Number residual (LinearAlgebra::CUDAWrappers::Vector< Number > &dst, const LinearAlgebra::CUDAWrappers::Vector< Number > &x, const LinearAlgebra::CUDAWrappers::Vector< Number > &b) const
 
Matrix norms
Number l1_norm () const
 
Number linfty_norm () const
 
Number frobenius_norm () const
 
Access to underlying CUDA data
std::tuple< Number *, int *, int *, cusparseMatDescr_t > get_cusparse_matrix () const
 
- Public Member Functions inherited from Subscriptor
 Subscriptor ()
 
 Subscriptor (const Subscriptor &)
 
 Subscriptor (Subscriptor &&) noexcept
 
virtual ~Subscriptor ()
 
Subscriptoroperator= (const Subscriptor &)
 
Subscriptoroperator= (Subscriptor &&) noexcept
 
void subscribe (std::atomic< bool > *const validity, const std::string &identifier="") const
 
void unsubscribe (std::atomic< bool > *const validity, const std::string &identifier="") const
 
unsigned int n_subscriptions () const
 
template<typename StreamType >
void list_subscribers (StreamType &stream) const
 
void list_subscribers () const
 
template<class Archive >
void serialize (Archive &ar, const unsigned int version)
 

Private Attributes

cusparseHandle_t cusparse_handle
 
int nnz
 
int n_rows
 
int n_cols
 
std::unique_ptr< Number[], void(*)(Number *)> val_dev
 
std::unique_ptr< int[], void(*)(int *)> column_index_dev
 
std::unique_ptr< int[], void(*)(int *)> row_ptr_dev
 
cusparseMatDescr_t descr
 

Additional Inherited Members

- Static Public Member Functions inherited from Subscriptor
static ::ExceptionBaseExcInUse (int arg1, std::string arg2, std::string arg3)
 
static ::ExceptionBaseExcNoSubscriber (std::string arg1, std::string arg2)
 

Detailed Description

template<typename Number>
class CUDAWrappers::SparseMatrix< Number >

This class is a wrapper around cuSPARSE csr sparse matrix. Unlike deal.II's own SparseMatrix all elements within each row are stored in increasing column index order.

Note
Instantiations for this template are provided for <float> and <double>.
Author
Bruno Turcksin
Date
2018

Definition at line 44 of file cuda_precondition.h.

Member Typedef Documentation

◆ size_type

template<typename Number>
using CUDAWrappers::SparseMatrix< Number >::size_type = int

Declare type for container size.

Definition at line 56 of file cuda_sparse_matrix.h.

◆ value_type

template<typename Number>
using CUDAWrappers::SparseMatrix< Number >::value_type = Number

Type of the matrix entries.

Definition at line 61 of file cuda_sparse_matrix.h.

◆ real_type

template<typename Number>
using CUDAWrappers::SparseMatrix< Number >::real_type = Number

Declare a type that holds real-valued numbers with the same precision as the template argument to this class.

Definition at line 67 of file cuda_sparse_matrix.h.

Constructor & Destructor Documentation

◆ SparseMatrix() [1/4]

template<typename Number>
CUDAWrappers::SparseMatrix< Number >::SparseMatrix ( )

Constructor. Initialize the matrix to be empty, without any structure, i.e., the matrix is not usable at all. This constructor is therefore only useful for matrices which are members of a class.

You have to initialize the matrix before usage with reinit.

◆ SparseMatrix() [2/4]

template<typename Number>
CUDAWrappers::SparseMatrix< Number >::SparseMatrix ( Utilities::CUDA::Handle handle,
const ::SparseMatrix< Number > &  sparse_matrix_host 
)

Constructor. Takes a Utilities::CUDA::Handle and a sparse matrix on the host. The sparse matrix on the host is copied on the device and the elements are reordered according to the format supported by cuSPARSE.

◆ SparseMatrix() [3/4]

template<typename Number>
CUDAWrappers::SparseMatrix< Number >::SparseMatrix ( CUDAWrappers::SparseMatrix< Number > &&  )

Move constructor. Create a new SparseMatrix by stealing the internal data.

◆ SparseMatrix() [4/4]

template<typename Number>
CUDAWrappers::SparseMatrix< Number >::SparseMatrix ( const CUDAWrappers::SparseMatrix< Number > &  )
delete

Copy constructor is deleted.

◆ ~SparseMatrix()

template<typename Number>
CUDAWrappers::SparseMatrix< Number >::~SparseMatrix ( )

Destructor. Free all memory.

Member Function Documentation

◆ operator=() [1/2]

template<typename Number>
SparseMatrix& CUDAWrappers::SparseMatrix< Number >::operator= ( CUDAWrappers::SparseMatrix< Number > &&  )

Move assignment operator.

◆ operator=() [2/2]

template<typename Number>
SparseMatrix& CUDAWrappers::SparseMatrix< Number >::operator= ( const CUDAWrappers::SparseMatrix< Number > &  )
delete

Copy assignment is deleted.

◆ reinit()

template<typename Number>
void CUDAWrappers::SparseMatrix< Number >::reinit ( Utilities::CUDA::Handle handle,
const ::SparseMatrix< Number > &  sparse_matrix_host 
)

Reinitialize the sparse matrix. The sparse matrix on the host is copied to the device and the elementes are reordered according to the format supported by cuSPARSE.

◆ m()

template<typename Number >
SparseMatrix< Number >::size_type SparseMatrix< Number >::m ( ) const
inline

Return the dimension of the codomain (or range) space. Note that the matrix is of dimension \(m \times n\).

Definition at line 378 of file cuda_sparse_matrix.h.

◆ n()

template<typename Number >
SparseMatrix< Number >::size_type SparseMatrix< Number >::n ( ) const
inline

Return the dimension of the domain space. Note that the matrix is of dimension \(m \times n\).

Definition at line 387 of file cuda_sparse_matrix.h.

◆ n_nonzero_elements()

template<typename Number >
std::size_t SparseMatrix< Number >::n_nonzero_elements ( ) const
inline

Return the number of nonzero elements of this matrix. Actually, it returns the number of entries in the sparsity pattern; if any of the entries should happen to be zero, it is counted anyway.

Definition at line 396 of file cuda_sparse_matrix.h.

◆ print()

template<typename Number >
template<class StreamType >
void SparseMatrix< Number >::print ( StreamType &  out,
const bool  across = false,
const bool  diagonal_first = true 
) const
inline

Print the matrix to the given stream, using the format (row,column) value, i.e. one nonzero entry of the matrix per line. If across is true, print all entries on a single line, using the format row,column:value.

If the argument diagonal_first is true, diagonal elements of quadratic matrices are printed first in their row. If it is false, the elements in a row are written in ascending column order.

Definition at line 406 of file cuda_sparse_matrix.h.

◆ print_formatted()

template<typename Number >
void SparseMatrix< Number >::print_formatted ( std::ostream &  out,
const unsigned int  precision = 3,
const bool  scientific = true,
const unsigned int  width = 0,
const char *  zero_string = " ",
const double  denominator = 1. 
) const

Print the matrix in the usual format, i.e. as a matrix and not as a list of nonzero elements. For better readability, elements not in the matrix are displayed as empty space, while matrix elements which are explicitly set to zero are displayed as such.

The parameters allow for a flexible setting of the output format: precision and scientific are used to determine the number format, where scientific = false means fixed point notation. A zero entry for width makes the function compute a width, but it may be changed to a positive value, if output is crude.

Additionally, a character for an empty value may be specified.

Finally, the whole matrix can be multiplied with a common denominator to produce more readable output, even integers.

Attention
This function may produce large amounts of output if applied to a large matrix!

Definition at line 462 of file cuda_sparse_matrix.h.

◆ operator*=()

template<typename Number>
SparseMatrix& CUDAWrappers::SparseMatrix< Number >::operator*= ( const Number  factor)

Multiply the entire matrix by a fixed factor.

◆ operator/=()

template<typename Number>
SparseMatrix& CUDAWrappers::SparseMatrix< Number >::operator/= ( const Number  factor)

Divide the entire matrix by a fixed factor.

◆ vmult()

template<typename Number>
void CUDAWrappers::SparseMatrix< Number >::vmult ( LinearAlgebra::CUDAWrappers::Vector< Number > &  dst,
const LinearAlgebra::CUDAWrappers::Vector< Number > &  src 
) const

Matrix-vector multiplication: let \(dst = M \cdot src\) with \(M\) being this matrix.

◆ Tvmult()

template<typename Number>
void CUDAWrappers::SparseMatrix< Number >::Tvmult ( LinearAlgebra::CUDAWrappers::Vector< Number > &  dst,
const LinearAlgebra::CUDAWrappers::Vector< Number > &  src 
) const

Matrix-vector multiplication: let \(dst = M^T \cdot src\) with \(M\) being this matrix. This function does the same as vmult() but takes this transposed matrix.

◆ vmult_add()

template<typename Number>
void CUDAWrappers::SparseMatrix< Number >::vmult_add ( LinearAlgebra::CUDAWrappers::Vector< Number > &  dst,
const LinearAlgebra::CUDAWrappers::Vector< Number > &  src 
) const

Adding matrix-vector multiplication. Add \(M \cdot src\) on \(dst\) with \(M\) being this matrix.

◆ Tvmult_add()

template<typename Number>
void CUDAWrappers::SparseMatrix< Number >::Tvmult_add ( LinearAlgebra::CUDAWrappers::Vector< Number > &  dst,
const LinearAlgebra::CUDAWrappers::Vector< Number > &  src 
) const

Adding matrix-vector multiplication. Add \(M^T \cdot src\) to \(dst\) with \(M\) being this matrix. This function foes the same as vmult_add() but takes the transposed matrix.

◆ matrix_norm_square()

template<typename Number>
Number CUDAWrappers::SparseMatrix< Number >::matrix_norm_square ( const LinearAlgebra::CUDAWrappers::Vector< Number > &  v) const

Return the square of the norm of the vector \(v\) with respect to the norm induced by this matrix, i.e., \(\left(v,Mv\right)\). This is useful, e.g., in the finite context, where the \(L_2\) norm of a function equals the matrix norm with respect to the mass matrix of the vector representing the nodal values of the finite element function.

Obviously, the matrix needs to be quadratic for this operation.

◆ matrix_scalar_product()

template<typename Number>
Number CUDAWrappers::SparseMatrix< Number >::matrix_scalar_product ( const LinearAlgebra::CUDAWrappers::Vector< Number > &  u,
const LinearAlgebra::CUDAWrappers::Vector< Number > &  v 
) const

Compute the matrix scalar product \(\left(u,Mv\right)\).

◆ residual()

template<typename Number>
Number CUDAWrappers::SparseMatrix< Number >::residual ( LinearAlgebra::CUDAWrappers::Vector< Number > &  dst,
const LinearAlgebra::CUDAWrappers::Vector< Number > &  x,
const LinearAlgebra::CUDAWrappers::Vector< Number > &  b 
) const

Compute the residual of an equation \(M \cdot x=b\), where the residual is defined to be \(r=b-M \cdot x\). Write the residual into \(dst\). The \(l_2\) norm of the residual vector is returned.

Source \(x\) and destination \(dst\) must not be the same vector.

◆ l1_norm()

template<typename Number>
Number CUDAWrappers::SparseMatrix< Number >::l1_norm ( ) const

Return the \(l_1\)-norm of the matrix, that is \(|M|_1=\max_{\mathrm{all\ columns\ }j}\sum_{\mathrm{all\ rows\ }i} |M_{ij}|\), (max. sum of columns). This is the natural matrix norm that is compatible to the \(l_1\)-norm for vectors, i.e., \(|Mv|_1\leq |M|_1 |v|_1\).

◆ linfty_norm()

template<typename Number>
Number CUDAWrappers::SparseMatrix< Number >::linfty_norm ( ) const

Return the \(l_\infty\)-norm of the matrix, that is \(|M|_\infty=\max_{\mathrm{all\ rows\ }i}\sum_{\mathrm{all\ columns\ }j} |M_{ij}|\), (max. sum of rows). This is the natural norm that is compatible to the \(l_\infty\)-norm of vectors, i.e., \(|Mv|_\infty \leq |M|_\infty |v|_\infty\).

◆ frobenius_norm()

template<typename Number>
Number CUDAWrappers::SparseMatrix< Number >::frobenius_norm ( ) const

Return the frobenius norm of the matrix, i.e., the square root of the sum of squares of all entries in the matrix.

◆ get_cusparse_matrix()

template<typename Number>
std::tuple<Number *, int *, int *, cusparseMatDescr_t> CUDAWrappers::SparseMatrix< Number >::get_cusparse_matrix ( ) const

Return a tuple containing the pointer to the values of matrix, the pointer to the columns indices, the pointer to the rows pointer, and the cuSPARSE matrix description.

Member Data Documentation

◆ cusparse_handle

template<typename Number>
cusparseHandle_t CUDAWrappers::SparseMatrix< Number >::cusparse_handle
private

cuSPARSE handle used to call cuSPARSE functions.

Definition at line 336 of file cuda_sparse_matrix.h.

◆ nnz

template<typename Number>
int CUDAWrappers::SparseMatrix< Number >::nnz
private

Number of non-zero elements in the sparse matrix.

Definition at line 341 of file cuda_sparse_matrix.h.

◆ n_rows

template<typename Number>
int CUDAWrappers::SparseMatrix< Number >::n_rows
private

Number of rows of the sparse matrix.

Definition at line 346 of file cuda_sparse_matrix.h.

◆ n_cols

template<typename Number>
int CUDAWrappers::SparseMatrix< Number >::n_cols
private

Number of columns of the sparse matrix.

Definition at line 351 of file cuda_sparse_matrix.h.

◆ val_dev

template<typename Number>
std::unique_ptr<Number[], void (*)(Number *)> CUDAWrappers::SparseMatrix< Number >::val_dev
private

Pointer to the values (on the device) of the sparse matrix.

Definition at line 356 of file cuda_sparse_matrix.h.

◆ column_index_dev

template<typename Number>
std::unique_ptr<int[], void (*)(int *)> CUDAWrappers::SparseMatrix< Number >::column_index_dev
private

Pointer to the column indices (on the device) of the sparse matrix.

Definition at line 361 of file cuda_sparse_matrix.h.

◆ row_ptr_dev

template<typename Number>
std::unique_ptr<int[], void (*)(int *)> CUDAWrappers::SparseMatrix< Number >::row_ptr_dev
private

Pointer to the row pointer (on the device) of the sparse matrix.

Definition at line 366 of file cuda_sparse_matrix.h.

◆ descr

template<typename Number>
cusparseMatDescr_t CUDAWrappers::SparseMatrix< Number >::descr
private

cuSPARSE description of the sparse matrix.

Definition at line 371 of file cuda_sparse_matrix.h.


The documentation for this class was generated from the following files: