Namespaces | Data Structures | Functions

viennacl::linalg Namespace Reference

Namespaces

namespace  kernels

Data Structures

struct  lower_tag
 A tag class representing a lower triangular matrix. More...
struct  upper_tag
 A tag class representing an upper triangular matrix. More...
struct  unit_lower_tag
 A tag class representing a lower triangular matrix with unit diagonal. More...
struct  unit_upper_tag
 A tag class representing an upper triangular matrix with unit diagonal. More...
class  no_precond
 A tag class representing the use of no preconditioner. More...
class  bicgstab_tag
 A tag for the stabilized Bi-conjugate gradient solver. Used for supplying solver parameters and for dispatching the solve() function. More...
class  cg_tag
 A tag for the conjugate gradient Used for supplying solver parameters and for dispatching the solve() function. More...
class  gmres_tag
 A tag for the solver GMRES. Used for supplying solver parameters and for dispatching the solve() function. More...
class  ilut_tag
 A tag for incomplete LU factorization with threshold (ILUT). More...
class  ilut_precond
 ILUT preconditioner class, can be supplied to solve()-routines. More...
class  ilut_precond< compressed_matrix< ScalarType, MAT_ALIGNMENT > >
 ILUT preconditioner class, can be supplied to solve()-routines. More...
class  jacobi_tag
 A tag for a jacobi preconditioner. More...
class  jacobi_precond
 Jacobi preconditioner class, can be supplied to solve()-routines. More...
class  jacobi_precond< compressed_matrix< ScalarType, MAT_ALIGNMENT > >
 Jacobi preconditioner class, can be supplied to solve()-routines. More...
class  row_scaling_tag
 A tag for a row preconditioner. More...
class  row_scaling
 Jacobi preconditioner class, can be supplied to solve()-routines. More...
class  row_scaling< compressed_matrix< ScalarType, MAT_ALIGNMENT > >
 Jacobi preconditioner class, can be supplied to solve()-routines. More...

Functions

template<class SCALARTYPE , unsigned int ALIGNMENT>
void norm_1_impl (const viennacl::vector< SCALARTYPE, ALIGNMENT > &vcl_vec, scalar< SCALARTYPE > &result)
 Computes the l^1-norm of a vector.
template<class SCALARTYPE , unsigned int ALIGNMENT>
void norm_2_impl (const viennacl::vector< SCALARTYPE, ALIGNMENT > &vcl_vec, scalar< SCALARTYPE > &result)
 Computes the l^2-norm of a vector - implementation.
template<class SCALARTYPE , unsigned int ALIGNMENT>
void norm_inf_impl (const viennacl::vector< SCALARTYPE, ALIGNMENT > &vcl_vec, scalar< SCALARTYPE > &result)
 Computes the supremum-norm of a vector.
template<class SCALARTYPE , typename F , unsigned int ALIGNMENT, unsigned int VECTOR_ALIGNMENT>
viennacl::vector_expression
< const viennacl::matrix
< SCALARTYPE, F, ALIGNMENT >
, const viennacl::vector
< SCALARTYPE, VECTOR_ALIGNMENT >
, op_prod > 
prod_impl (const viennacl::matrix< SCALARTYPE, F, ALIGNMENT > &mat, const viennacl::vector< SCALARTYPE, VECTOR_ALIGNMENT > &vec)
 Returns a proxy class that represents matrix-vector multiplication.
template<class SCALARTYPE , unsigned int ALIGNMENT, unsigned int VECTOR_ALIGNMENT>
viennacl::vector_expression
< const
viennacl::compressed_matrix
< SCALARTYPE, ALIGNMENT >
, const viennacl::vector
< SCALARTYPE, VECTOR_ALIGNMENT >
, op_prod > 
prod_impl (const compressed_matrix< SCALARTYPE, ALIGNMENT > &mat, const vector< SCALARTYPE, VECTOR_ALIGNMENT > &vec)
 Returns a proxy class that represents matrix-vector multiplication with a compressed_matrix.
template<class SCALARTYPE , unsigned int ALIGNMENT, unsigned int VECTOR_ALIGNMENT>
viennacl::vector_expression
< const
viennacl::coordinate_matrix
< SCALARTYPE, ALIGNMENT >
, const viennacl::vector
< SCALARTYPE, VECTOR_ALIGNMENT >
, op_prod > 
prod_impl (const coordinate_matrix< SCALARTYPE, ALIGNMENT > &mat, const vector< SCALARTYPE, VECTOR_ALIGNMENT > &vec)
 Returns a proxy class that represents matrix-vector multiplication with a compressed_matrix.
template<class SCALARTYPE , unsigned int ALIGNMENT1, unsigned int ALIGNMENT2>
viennacl::scalar_expression
< const viennacl::vector
< SCALARTYPE, ALIGNMENT1 >
, const viennacl::vector
< SCALARTYPE, ALIGNMENT2 >
, viennacl::op_inner_prod > 
inner_prod_impl (const viennacl::vector< SCALARTYPE, ALIGNMENT1 > &vec1, const viennacl::vector< SCALARTYPE, ALIGNMENT2 > &vec2)
 Computes the inner product of two vectors.
template<class SCALARTYPE , unsigned int ALIGNMENT>
void inner_prod_impl (const viennacl::vector< SCALARTYPE, ALIGNMENT > &vec1, const viennacl::vector< SCALARTYPE, ALIGNMENT > &vec2, scalar< SCALARTYPE > &result)
 Computes the inner product of two vectors - implementation. Library users should call inner_prod(vec1, vec2).
template<typename MatrixType , typename VectorType >
VectorType solve (const MatrixType &matrix, VectorType const &rhs, bicgstab_tag const &tag)
 Implementation of the stabilized Bi-conjugate gradient solver.
template<typename MatrixType , typename VectorType >
VectorType solve (const MatrixType &matrix, VectorType const &rhs, bicgstab_tag const &tag, viennacl::linalg::no_precond)
template<typename MatrixType , typename VectorType , typename PreconditionerType >
VectorType solve (const MatrixType &matrix, VectorType const &rhs, bicgstab_tag const &tag, PreconditionerType const &precond)
 Implementation of the preconditioned stabilized Bi-conjugate gradient solver.
template<typename MatrixType , typename VectorType >
VectorType solve (const MatrixType &matrix, VectorType const &rhs, cg_tag const &tag)
 Implementation of the conjugate gradient solver without preconditioner.
template<typename MatrixType , typename VectorType >
VectorType solve (const MatrixType &matrix, VectorType const &rhs, cg_tag const &tag, viennacl::linalg::no_precond)
template<typename MatrixType , typename VectorType , typename PreconditionerType >
VectorType solve (const MatrixType &matrix, VectorType const &rhs, cg_tag const &tag, PreconditionerType const &precond)
 Implementation of the preconditioned conjugate gradient solver.
template<class TYPE , unsigned int ALIGNMENT, unsigned int VECTOR_ALIGNMENT>
void prod_impl (const viennacl::compressed_matrix< TYPE, ALIGNMENT > &mat, const viennacl::vector< TYPE, VECTOR_ALIGNMENT > &vec, viennacl::vector< TYPE, VECTOR_ALIGNMENT > &result, size_t NUM_THREADS=0)
 Carries out matrix-vector multiplication with a compressed_matrix.
template<typename SCALARTYPE , unsigned int MAT_ALIGNMENT, unsigned int VEC_ALIGNMENT>
void inplace_solve (compressed_matrix< SCALARTYPE, MAT_ALIGNMENT > const &L, vector< SCALARTYPE, VEC_ALIGNMENT > &vec, viennacl::linalg::unit_lower_tag)
 Inplace solution of a lower triangular compressed_matrix with unit diagonal. Typically used for LU substitutions.
template<typename SCALARTYPE , unsigned int MAT_ALIGNMENT, unsigned int VEC_ALIGNMENT, typename TAG >
vector< SCALARTYPE, VEC_ALIGNMENT > solve (compressed_matrix< SCALARTYPE, MAT_ALIGNMENT > const &L, const vector< SCALARTYPE, VEC_ALIGNMENT > &vec, const viennacl::linalg::unit_lower_tag &tag)
 Convenience functions for result = solve(trans(mat), vec, unit_lower_tag()); Creates a temporary result vector and forwards the request to inplace_solve().
template<typename SCALARTYPE , unsigned int MAT_ALIGNMENT, unsigned int VEC_ALIGNMENT>
void inplace_solve (compressed_matrix< SCALARTYPE, MAT_ALIGNMENT > const &U, vector< SCALARTYPE, VEC_ALIGNMENT > &vec, viennacl::linalg::upper_tag)
 Inplace solution of a upper triangular compressed_matrix. Typically used for LU substitutions.
template<typename SCALARTYPE , unsigned int MAT_ALIGNMENT, unsigned int VEC_ALIGNMENT, typename TAG >
vector< SCALARTYPE, VEC_ALIGNMENT > solve (compressed_matrix< SCALARTYPE, MAT_ALIGNMENT > const &L, const vector< SCALARTYPE, VEC_ALIGNMENT > &vec, viennacl::linalg::upper_tag const &tag)
 Convenience functions for result = solve(trans(mat), vec, unit_lower_tag()); Creates a temporary result vector and forwards the request to inplace_solve().
template<class SCALARTYPE , unsigned int ALIGNMENT, unsigned int VECTOR_ALIGNMENT>
viennacl::vector_expression
< const
viennacl::coordinate_matrix
< SCALARTYPE, ALIGNMENT >
, const viennacl::vector
< SCALARTYPE, VECTOR_ALIGNMENT >
, viennacl::op_prod > 
prod_impl (const viennacl::coordinate_matrix< SCALARTYPE, ALIGNMENT > &mat, const viennacl::vector< SCALARTYPE, VECTOR_ALIGNMENT > &vec, size_t NUM_THREADS)
 Returns a proxy class that represents matrix-vector multiplication with a coordinate_matrix.
template<class TYPE , unsigned int ALIGNMENT, unsigned int VECTOR_ALIGNMENT>
void prod_impl (const viennacl::coordinate_matrix< TYPE, ALIGNMENT > &mat, const viennacl::vector< TYPE, VECTOR_ALIGNMENT > &vec, viennacl::vector< TYPE, VECTOR_ALIGNMENT > &result)
 Carries out matrix-vector multiplication with a coordinate_matrix.
template<typename SCALARTYPE , typename F1 , typename F2 , unsigned int A1, unsigned int A2, typename SOLVERTAG >
void inplace_solve (const matrix< SCALARTYPE, F1, A1 > &mat, matrix< SCALARTYPE, F2, A2 > &B, SOLVERTAG)
 Direct inplace solver for dense upper triangular systems.
template<typename SCALARTYPE , typename F1 , typename F2 , unsigned int A1, unsigned int A2, typename SOLVERTAG >
void inplace_solve (const matrix< SCALARTYPE, F1, A1 > &mat, const matrix_expression< const matrix< SCALARTYPE, F2, A2 >, const matrix< SCALARTYPE, F2, A2 >, op_trans > &B, SOLVERTAG)
 Direct inplace solver for dense upper triangular systems.
template<typename SCALARTYPE , typename F1 , typename F2 , unsigned int A1, unsigned int A2, typename SOLVERTAG >
void inplace_solve (const matrix_expression< const matrix< SCALARTYPE, F1, A1 >, const matrix< SCALARTYPE, F1, A1 >, op_trans > &proxy, matrix< SCALARTYPE, F2, A2 > &B, SOLVERTAG)
 Direct inplace solver for dense upper triangular systems that stem from transposed lower triangular systems.
template<typename SCALARTYPE , typename F1 , typename F2 , unsigned int A1, unsigned int A2, typename SOLVERTAG >
void inplace_solve (const matrix_expression< const matrix< SCALARTYPE, F1, A1 >, const matrix< SCALARTYPE, F1, A1 >, op_trans > &proxy, const matrix_expression< const matrix< SCALARTYPE, F2, A2 >, const matrix< SCALARTYPE, F2, A2 >, op_trans > &B, SOLVERTAG)
 Direct inplace solver for dense upper triangular systems that stem from transposed lower triangular systems.
template<typename SCALARTYPE , typename F , unsigned int ALIGNMENT, unsigned int VEC_ALIGNMENT, typename SOLVERTAG >
void inplace_solve (const matrix< SCALARTYPE, F, ALIGNMENT > &mat, vector< SCALARTYPE, VEC_ALIGNMENT > &vec, SOLVERTAG)
template<typename SCALARTYPE , typename F , unsigned int ALIGNMENT, unsigned int VEC_ALIGNMENT, typename SOLVERTAG >
void inplace_solve (const matrix_expression< const matrix< SCALARTYPE, F, ALIGNMENT >, const matrix< SCALARTYPE, F, ALIGNMENT >, op_trans > &proxy, vector< SCALARTYPE, VEC_ALIGNMENT > &vec, SOLVERTAG)
 Direct inplace solver for dense upper triangular systems that stem from transposed lower triangular systems.
template<typename SCALARTYPE , typename F1 , typename F2 , unsigned int ALIGNMENT_A, unsigned int ALIGNMENT_B, typename TAG >
matrix< SCALARTYPE, F2,
ALIGNMENT_B > 
solve (const matrix< SCALARTYPE, F1, ALIGNMENT_A > &A, const matrix< SCALARTYPE, F2, ALIGNMENT_B > &B, TAG const &tag)
 Convenience functions for C = solve(A, B, some_tag()); Creates a temporary result matrix and forwards the request to inplace_solve().
template<typename SCALARTYPE , typename F1 , typename F2 , unsigned int ALIGNMENT_A, unsigned int ALIGNMENT_B, typename TAG >
matrix< SCALARTYPE, F2,
ALIGNMENT_B > 
solve (const matrix< SCALARTYPE, F1, ALIGNMENT_A > &A, const matrix_expression< const matrix< SCALARTYPE, F2, ALIGNMENT_B >, const matrix< SCALARTYPE, F2, ALIGNMENT_B >, op_trans > &proxy, TAG const &tag)
 Convenience functions for C = solve(A, B^T, some_tag()); Creates a temporary result matrix and forwards the request to inplace_solve().
template<typename SCALARTYPE , typename F , unsigned int ALIGNMENT, unsigned int VEC_ALIGNMENT, typename TAG >
vector< SCALARTYPE, VEC_ALIGNMENT > solve (const matrix< SCALARTYPE, F, ALIGNMENT > &mat, const vector< SCALARTYPE, VEC_ALIGNMENT > &vec, TAG const &tag)
 Convenience functions for result = solve(mat, vec, some_tag()); Creates a temporary result vector and forwards the request to inplace_solve().
template<typename SCALARTYPE , typename F1 , typename F2 , unsigned int ALIGNMENT_A, unsigned int ALIGNMENT_B, typename TAG >
matrix< SCALARTYPE, F2,
ALIGNMENT_B > 
solve (const matrix_expression< const matrix< SCALARTYPE, F1, ALIGNMENT_A >, const matrix< SCALARTYPE, F1, ALIGNMENT_A >, op_trans > &proxy, const matrix< SCALARTYPE, F2, ALIGNMENT_B > &B, TAG const &tag)
 Convenience functions for result = solve(trans(mat), B, some_tag()); Creates a temporary result matrix and forwards the request to inplace_solve().
template<typename SCALARTYPE , typename F1 , typename F2 , unsigned int ALIGNMENT_A, unsigned int ALIGNMENT_B, typename TAG >
matrix< SCALARTYPE, F2,
ALIGNMENT_B > 
solve (const matrix_expression< const matrix< SCALARTYPE, F1, ALIGNMENT_A >, const matrix< SCALARTYPE, F1, ALIGNMENT_A >, op_trans > &proxy_A, const matrix_expression< const matrix< SCALARTYPE, F2, ALIGNMENT_B >, const matrix< SCALARTYPE, F2, ALIGNMENT_B >, op_trans > &proxy_B, TAG const &tag)
 Convenience functions for result = solve(trans(mat), vec, some_tag()); Creates a temporary result vector and forwards the request to inplace_solve().
template<typename SCALARTYPE , typename F , unsigned int ALIGNMENT, unsigned int VEC_ALIGNMENT, typename TAG >
vector< SCALARTYPE, VEC_ALIGNMENT > solve (const matrix_expression< const matrix< SCALARTYPE, F, ALIGNMENT >, const matrix< SCALARTYPE, F, ALIGNMENT >, op_trans > &proxy, const vector< SCALARTYPE, VEC_ALIGNMENT > &vec, TAG const &tag)
 Convenience functions for result = solve(trans(mat), vec, some_tag()); Creates a temporary result vector and forwards the request to inplace_solve().
template<typename SCALARTYPE , typename F , unsigned int ALIGNMENT>
void lu_factorize (matrix< SCALARTYPE, F, ALIGNMENT > &mat)
 LU factorization of a dense matrix.
template<typename SCALARTYPE , typename F1 , typename F2 , unsigned int ALIGNMENT_A, unsigned int ALIGNMENT_B>
void lu_substitute (matrix< SCALARTYPE, F1, ALIGNMENT_A > const &A, matrix< SCALARTYPE, F2, ALIGNMENT_B > &B)
 LU substitution for the system LU = rhs.
template<typename SCALARTYPE , typename F , unsigned int ALIGNMENT, unsigned int VEC_ALIGNMENT>
void lu_substitute (matrix< SCALARTYPE, F, ALIGNMENT > const &mat, vector< SCALARTYPE, VEC_ALIGNMENT > &vec)
 LU substitution for the system LU = rhs.
template<typename MatrixType , typename VectorType , typename PreconditionerType >
VectorType solve (const MatrixType &matrix, VectorType const &rhs, gmres_tag const &tag, PreconditionerType const &precond)
 Implementation of the GMRES solver.
template<typename MatrixType , typename VectorType >
VectorType solve (const MatrixType &matrix, VectorType const &rhs, gmres_tag const &tag)
 Convenience overload of the solve() function using GMRES. Per default, no preconditioner is used.
template<typename T >
void ilut_inc_row_iterator_to_row_index (T &row_iter, unsigned int k)
 Increments a row iterator (iteration along increasing row indices) up to a certain row index k.
template<typename ScalarType >
void ilut_inc_row_iterator_to_row_index (viennacl::tools::sparse_matrix_adapter< ScalarType > &row_iter, unsigned int k)
 Increments a row iterator (iteration along increasing row indices) up to a certain row index k.
template<typename ScalarType >
void ilut_inc_row_iterator_to_row_index (viennacl::tools::const_sparse_matrix_adapter< ScalarType > &row_iter, unsigned int k)
 Increments a row iterator (iteration along increasing row indices) up to a certain row index k.
template<typename MatrixType , typename LUType >
void precondition (MatrixType const &input, LUType &output, ilut_tag const &tag)
 Implementation of a ILU-preconditioner with threshold.
template<typename MatrixType , typename VectorType >
void ilu_inplace_solve (MatrixType const &mat, VectorType &vec, viennacl::linalg::unit_lower_tag)
 Generic inplace solution of a unit lower triangular system.
template<typename MatrixType , typename VectorType >
void ilu_inplace_solve (MatrixType const &mat, VectorType &vec, viennacl::linalg::upper_tag)
 Generic inplace solution of a upper triangular system.
template<typename MatrixType , typename VectorType >
void ilu_lu_substitute (MatrixType const &mat, VectorType &vec)
 Generic LU substitution.
template<typename VectorT1 , typename VectorT2 >
VectorT1::value_type inner_prod (VectorT1 const &v1, VectorT2 const &v2, typename viennacl::tools::enable_if< viennacl::is_stl< typename viennacl::traits::tag_of< VectorT1 >::type >::value >::type *dummy=0)
template<typename ScalarType , unsigned int alignment1, unsigned int alignment2>
viennacl::scalar_expression
< const viennacl::vector
< ScalarType, alignment1 >
, const viennacl::vector
< ScalarType, alignment2 >
, viennacl::op_inner_prod > 
inner_prod (viennacl::vector< ScalarType, alignment1 > const &vector1, viennacl::vector< ScalarType, alignment2 > const &vector2, typename viennacl::tools::enable_if< viennacl::is_viennacl< typename viennacl::traits::tag_of< viennacl::vector< ScalarType, alignment1 > >::type >::value >::type *dummy=0)
template<class TYPE , typename F , unsigned int ALIGNMENT>
void add (const viennacl::matrix< TYPE, F, ALIGNMENT > &mat1, const viennacl::matrix< TYPE, F, ALIGNMENT > &mat2, viennacl::matrix< TYPE, F, ALIGNMENT > &result)
 Adds two dense matrices and writes the result to a third matrix.
template<class TYPE , typename F , unsigned int ALIGNMENT>
void inplace_add (viennacl::matrix< TYPE, F, ALIGNMENT > &result, const viennacl::matrix< TYPE, F, ALIGNMENT > &mat2)
 Adds a dense matrix to another.
template<class TYPE , typename F , unsigned int ALIGNMENT>
void sub (const viennacl::matrix< TYPE, F, ALIGNMENT > &mat1, const viennacl::matrix< TYPE, F, ALIGNMENT > &mat2, viennacl::matrix< TYPE, F, ALIGNMENT > &result)
 Adds two dense matrices and writes the result to a third matrix.
template<class TYPE , typename F , unsigned int ALIGNMENT>
void inplace_sub (viennacl::matrix< TYPE, F, ALIGNMENT > &result, const viennacl::matrix< TYPE, F, ALIGNMENT > &mat2)
 Adds a dense matrix to another.
template<class SCALARTYPE , typename F , unsigned int ALIGNMENT>
void inplace_mult (viennacl::matrix< SCALARTYPE, F, ALIGNMENT > &result, SCALARTYPE val)
 Multiplies a dense matrix by a scalar.
template<class SCALARTYPE , typename F , unsigned int ALIGNMENT>
void inplace_mult (viennacl::matrix< SCALARTYPE, F, ALIGNMENT > &result, viennacl::scalar< SCALARTYPE > const &val)
 Multiplies a dense matrix by a scalar.
template<class SCALARTYPE , typename F , unsigned int ALIGNMENT>
void inplace_divide (viennacl::matrix< SCALARTYPE, F, ALIGNMENT > &result, viennacl::scalar< SCALARTYPE > const &val)
 Multiplies a dense matrix by a scalar.
template<class TYPE , typename F , unsigned int ALIGNMENT, unsigned int VECTOR_ALIGNMENT>
void prod_impl (const viennacl::matrix< TYPE, F, ALIGNMENT > &mat, const viennacl::vector< TYPE, VECTOR_ALIGNMENT > &vec, viennacl::vector< TYPE, VECTOR_ALIGNMENT > &result)
 Carries out matrix-vector multiplication.
template<class SCALARTYPE , typename F , unsigned int ALIGNMENT, unsigned int VECTOR_ALIGNMENT>
viennacl::vector_expression
< const
viennacl::matrix_expression
< const matrix< SCALARTYPE, F,
ALIGNMENT >, const matrix
< SCALARTYPE, F, ALIGNMENT >
, op_trans >, const
viennacl::vector< SCALARTYPE,
VECTOR_ALIGNMENT >, op_prod > 
prod_impl (const viennacl::matrix_expression< const matrix< SCALARTYPE, F, ALIGNMENT >, const matrix< SCALARTYPE, F, ALIGNMENT >, op_trans > &proxy, const viennacl::vector< SCALARTYPE, VECTOR_ALIGNMENT > &vec)
 Returns a proxy class that represents matrix-vector multiplication with a transposed matrix.
template<class SCALARTYPE , typename F , unsigned int ALIGNMENT, unsigned int VECTOR_ALIGNMENT>
void prod_impl (const viennacl::matrix_expression< const matrix< SCALARTYPE, F, ALIGNMENT >, const matrix< SCALARTYPE, F, ALIGNMENT >, op_trans > &mat, const viennacl::vector< SCALARTYPE, VECTOR_ALIGNMENT > &vec, viennacl::vector< SCALARTYPE, VECTOR_ALIGNMENT > &result)
 Unwraps the transposed matrix proxy and forwards to trans_prod_impl().
template<class SCALARTYPE , typename F , unsigned int ALIGNMENT, unsigned int VECTOR_ALIGNMENT>
void trans_prod_impl (const matrix< SCALARTYPE, F, ALIGNMENT > &mat, const viennacl::vector< SCALARTYPE, VECTOR_ALIGNMENT > &vec, viennacl::vector< SCALARTYPE, VECTOR_ALIGNMENT > &result)
 Carries out matrix-vector multiplication with a transposed matrix.
template<class TYPE , typename F1 , typename F2 , typename F3 , unsigned int ALIGNMENT>
void prod_impl (const viennacl::matrix< TYPE, F1, ALIGNMENT > &A, const viennacl::matrix< TYPE, F2, ALIGNMENT > &B, viennacl::matrix< TYPE, F3, ALIGNMENT > &C, int block_size=15)
 Carries out matrix-matrix multiplication.
template<class TYPE , typename F1 , typename F2 , typename F3 , unsigned int ALIGNMENT>
void prod_impl (const viennacl::matrix_expression< const matrix< TYPE, F1, ALIGNMENT >, const matrix< TYPE, F1, ALIGNMENT >, op_trans > &A, const viennacl::matrix< TYPE, F2, ALIGNMENT > &B, viennacl::matrix< TYPE, F3, ALIGNMENT > &C)
 Carries out matrix-matrix multiplication.
template<class TYPE , typename F1 , typename F2 , typename F3 , unsigned int ALIGNMENT>
void prod_impl (const viennacl::matrix< TYPE, F1, ALIGNMENT > &A, const viennacl::matrix_expression< const matrix< TYPE, F2, ALIGNMENT >, const matrix< TYPE, F2, ALIGNMENT >, op_trans > &B, viennacl::matrix< TYPE, F3, ALIGNMENT > &C)
 Carries out matrix-matrix multiplication.
template<class TYPE , typename F1 , typename F2 , typename F3 , unsigned int ALIGNMENT>
void prod_impl (const viennacl::matrix_expression< const matrix< TYPE, F1, ALIGNMENT >, const matrix< TYPE, F1, ALIGNMENT >, op_trans > &A, const viennacl::matrix_expression< const matrix< TYPE, F2, ALIGNMENT >, const matrix< TYPE, F2, ALIGNMENT >, op_trans > &B, viennacl::matrix< TYPE, F3, ALIGNMENT > &C)
 Carries out matrix-matrix multiplication.
template<class SCALARTYPE , unsigned int VA1, unsigned int VA2>
viennacl::matrix_expression
< const viennacl::vector
< SCALARTYPE, VA1 >, const
viennacl::vector< SCALARTYPE,
VA2 >, op_prod > 
outer_prod (const viennacl::vector< SCALARTYPE, VA1 > &vec1, const viennacl::vector< SCALARTYPE, VA2 > &vec2)
 Returns a proxy class for the operation mat += vec1 * vec2^T, i.e. a rank 1 update.
template<class SCALARTYPE , typename F , unsigned int ALIGNMENT>
void rank_1_update (viennacl::matrix< SCALARTYPE, F, ALIGNMENT > &mat1, const viennacl::vector< SCALARTYPE, ALIGNMENT > &vec1, const viennacl::vector< SCALARTYPE, ALIGNMENT > &vec2)
 The implementation of the operation mat += vec1 * vec2^T, i.e. a rank 1 update.
template<class SCALARTYPE , typename F , unsigned int ALIGNMENT>
void scaled_rank_1_update (viennacl::matrix< SCALARTYPE, F, ALIGNMENT > &mat1, SCALARTYPE val, const viennacl::vector< SCALARTYPE, ALIGNMENT > &vec1, const viennacl::vector< SCALARTYPE, ALIGNMENT > &vec2)
 The implementation of the operation mat += alpha * vec1 * vec2^T, i.e. a scaled rank 1 update.
template<typename VectorT >
VectorT::value_type norm_1 (VectorT const &v1, typename viennacl::tools::enable_if< viennacl::is_stl< typename viennacl::traits::tag_of< VectorT >::type >::value >::type *dummy=0)
template<typename ScalarType , unsigned int alignment>
viennacl::scalar_expression
< const viennacl::vector
< ScalarType, alignment >
, const viennacl::vector
< ScalarType, alignment >
, viennacl::op_norm_1 > 
norm_1 (viennacl::vector< ScalarType, alignment > const &vector, typename viennacl::tools::enable_if< viennacl::is_viennacl< typename viennacl::traits::tag_of< viennacl::vector< ScalarType, alignment > >::type >::value >::type *dummy=0)
template<typename VectorT >
VectorT::value_type norm_2 (VectorT const &v1, typename viennacl::tools::enable_if< viennacl::is_stl< typename viennacl::traits::tag_of< VectorT >::type >::value >::type *dummy=0)
template<typename ScalarType , unsigned int alignment>
viennacl::scalar_expression
< const viennacl::vector
< ScalarType, alignment >
, const viennacl::vector
< ScalarType, alignment >
, viennacl::op_norm_2 > 
norm_2 (viennacl::vector< ScalarType, alignment > const &v, typename viennacl::tools::enable_if< viennacl::is_viennacl< typename viennacl::traits::tag_of< viennacl::vector< ScalarType, alignment > >::type >::value >::type *dummy=0)
template<typename VectorT >
VectorT::value_type norm_inf (VectorT const &v1, typename viennacl::tools::enable_if< viennacl::is_stl< typename viennacl::traits::tag_of< VectorT >::type >::value >::type *dummy=0)
template<typename ScalarType , unsigned int alignment>
viennacl::scalar_expression
< const viennacl::vector
< ScalarType, alignment >
, const viennacl::vector
< ScalarType, alignment >
, viennacl::op_norm_inf > 
norm_inf (viennacl::vector< ScalarType, alignment > const &v1, typename viennacl::tools::enable_if< viennacl::is_viennacl< typename viennacl::traits::tag_of< viennacl::vector< ScalarType, alignment > >::type >::value >::type *dummy=0)
template<typename T , typename A1 , typename A2 , typename VectorT >
VectorT prod_impl (std::vector< std::vector< T, A1 >, A2 > const &matrix, VectorT const &vector)
template<typename KEY , typename DATA , typename COMPARE , typename AMAP , typename AVEC , typename VectorT >
VectorT prod_impl (std::vector< std::map< KEY, DATA, COMPARE, AMAP >, AVEC > const &matrix, VectorT const &vector)
template<typename MatrixT , typename VectorT >
VectorT prod (MatrixT const &matrix, VectorT const &vector, typename viennacl::tools::enable_if< viennacl::is_stl< typename viennacl::traits::tag_of< MatrixT >::type >::value >::type *dummy=0)
template<typename MatrixT , typename NumericT , unsigned int ALIGNMENT>
viennacl::vector_expression
< const MatrixT, const
viennacl::vector< NumericT,
ALIGNMENT >, viennacl::op_prod > 
prod (MatrixT const &matrix, viennacl::vector< NumericT, ALIGNMENT > const &vector, typename viennacl::tools::enable_if< viennacl::is_viennacl< typename viennacl::traits::tag_of< MatrixT >::type >::value >::type *dummy=0)
template<typename MatrixT , typename NumericT , typename F , unsigned int ALIGNMENT>
viennacl::matrix_expression
< const MatrixT, const
viennacl::matrix< NumericT, F,
ALIGNMENT >, viennacl::op_prod > 
prod (MatrixT const &matrix_A, viennacl::matrix< NumericT, F, ALIGNMENT > const &matrix_B, typename viennacl::tools::enable_if< viennacl::is_viennacl< typename viennacl::traits::tag_of< MatrixT >::type >::value >::type *dummy=0)
template<typename MatrixT , typename NumericT , typename F , unsigned int ALIGNMENT>
viennacl::matrix_expression
< const MatrixT, const
viennacl::matrix_expression
< const viennacl::matrix
< NumericT, F, ALIGNMENT >
, const viennacl::matrix
< NumericT, F, ALIGNMENT >
, viennacl::op_trans >
, viennacl::op_prod > 
prod (MatrixT const &matrix_A, const viennacl::matrix_expression< const viennacl::matrix< NumericT, F, ALIGNMENT >, const viennacl::matrix< NumericT, F, ALIGNMENT >, viennacl::op_trans > &matrix_B, typename viennacl::tools::enable_if< viennacl::is_viennacl< typename viennacl::traits::tag_of< MatrixT >::type >::value >::type *dummy=0)
template<class SCALARTYPE , unsigned int ALIGNMENT>
void add (const viennacl::vector< SCALARTYPE, ALIGNMENT > &vec1, const viennacl::vector< SCALARTYPE, ALIGNMENT > &vec2, viennacl::vector< SCALARTYPE, ALIGNMENT > &result)
 Addition of two vectors. Try to use the overloaded operators for vector instead, unless you want to fine-tune the number of GPU threads involved.
template<class SCALARTYPE , unsigned int ALIGNMENT>
void inplace_add (viennacl::vector< SCALARTYPE, ALIGNMENT > &vec1, const viennacl::vector< SCALARTYPE, ALIGNMENT > &vec2)
 Inplace addition of two vectors. Try to use the overloaded operators for vector instead, unless you want to fine-tune the number of GPU threads involved.
template<class SCALARTYPE , unsigned int ALIGNMENT>
void sub (const viennacl::vector< SCALARTYPE, ALIGNMENT > &vec1, const viennacl::vector< SCALARTYPE, ALIGNMENT > &vec2, viennacl::vector< SCALARTYPE, ALIGNMENT > &result)
 Subtraction of two vectors. Try to use the overloaded operators for vector instead, unless you want to fine-tune the number of GPU threads involved.
template<class SCALARTYPE , unsigned int ALIGNMENT>
void inplace_sub (viennacl::vector< SCALARTYPE, ALIGNMENT > &vec1, const viennacl::vector< SCALARTYPE, ALIGNMENT > &vec2)
 Inplace addition of two vectors. Try to use the overloaded operators for vector instead, unless you want to fine-tune the number of GPU threads involved.
template<class SCALARTYPE , unsigned int ALIGNMENT>
void mult (const viennacl::vector< SCALARTYPE, ALIGNMENT > &vec, scalar< SCALARTYPE > const &alpha, viennacl::vector< SCALARTYPE, ALIGNMENT > &result)
 Scales a vector. Try to use the overloaded operators for vector instead, unless you want to fine-tune the number of GPU threads involved.
template<class SCALARTYPE , unsigned int ALIGNMENT>
void mult (const viennacl::vector< SCALARTYPE, ALIGNMENT > &vec, SCALARTYPE alpha, viennacl::vector< SCALARTYPE, ALIGNMENT > &result)
 Scales a vector. Try to use the overloaded operators for vector instead, unless you want to fine-tune the number of GPU threads involved.
template<class SCALARTYPE , unsigned int ALIGNMENT>
void inplace_mult (viennacl::vector< SCALARTYPE, ALIGNMENT > &vec, scalar< SCALARTYPE > const &alpha)
 Scales a vector inplace. Try to use the overloaded operators for vector instead, unless you want to fine-tune the number of GPU threads involved.
template<class SCALARTYPE , unsigned int ALIGNMENT>
void inplace_mult (viennacl::vector< SCALARTYPE, ALIGNMENT > &vec, SCALARTYPE alpha)
 Scales a vector inplace. Try to use the overloaded operators for vector instead, unless you want to fine-tune the number of GPU threads involved.
template<class SCALARTYPE , unsigned int ALIGNMENT>
void divide (const viennacl::vector< SCALARTYPE, ALIGNMENT > &vec, scalar< SCALARTYPE > const &alpha, viennacl::vector< SCALARTYPE, ALIGNMENT > &result)
 Scales a vector. Try to use the overloaded operators for vector instead, unless you want to fine-tune the number of GPU threads involved.
template<class SCALARTYPE , unsigned int ALIGNMENT>
void inplace_divide (viennacl::vector< SCALARTYPE, ALIGNMENT > &vec, scalar< SCALARTYPE > const &alpha)
 Scales a vector inplace. Try to use the overloaded operators for vector instead, unless you want to fine-tune the number of GPU threads involved.
template<class SCALARTYPE , unsigned int ALIGNMENT>
void mul_add (const viennacl::vector< SCALARTYPE, ALIGNMENT > &vec1, scalar< SCALARTYPE > const &alpha, const viennacl::vector< SCALARTYPE, ALIGNMENT > &vec2, viennacl::vector< SCALARTYPE, ALIGNMENT > &result)
 Multiply-add operation. Try to use the overloaded operators for vector instead, unless you want to fine-tune the number of GPU threads involved.
template<class SCALARTYPE , unsigned int ALIGNMENT>
void mul_add (const viennacl::vector< SCALARTYPE, ALIGNMENT > &vec1, SCALARTYPE alpha, const viennacl::vector< SCALARTYPE, ALIGNMENT > &vec2, viennacl::vector< SCALARTYPE, ALIGNMENT > &result)
 Multiply-add operation. Try to use the overloaded operators for vector instead, unless you want to fine-tune the number of GPU threads involved.
template<class SCALARTYPE , unsigned int ALIGNMENT>
void inplace_mul_add (viennacl::vector< SCALARTYPE, ALIGNMENT > &vec1, const viennacl::vector< SCALARTYPE, ALIGNMENT > &vec2, scalar< SCALARTYPE > const &alpha)
 Inplace Multiply-add operation. Try to use the overloaded operators for vector instead, unless you want to fine-tune the number of GPU threads involved.
template<class SCALARTYPE , unsigned int ALIGNMENT>
void inplace_mul_add (viennacl::vector< SCALARTYPE, ALIGNMENT > &vec1, const viennacl::vector< SCALARTYPE, ALIGNMENT > &vec2, SCALARTYPE alpha)
 Inplace Multiply-add operation. Try to use the overloaded operators for vector instead, unless you want to fine-tune the number of GPU threads involved.
template<class SCALARTYPE , unsigned int ALIGNMENT>
void mul_sub (const viennacl::vector< SCALARTYPE, ALIGNMENT > &vec1, scalar< SCALARTYPE > const &alpha, const viennacl::vector< SCALARTYPE, ALIGNMENT > &vec2, viennacl::vector< SCALARTYPE, ALIGNMENT > &result)
 Multiply-subtract operation. Try to use the overloaded operators for vector instead, unless you want to fine-tune the number of GPU threads involved.
template<class SCALARTYPE , unsigned int ALIGNMENT>
void inplace_mul_sub (viennacl::vector< SCALARTYPE, ALIGNMENT > &vec1, const viennacl::vector< SCALARTYPE, ALIGNMENT > &vec2, scalar< SCALARTYPE > const &alpha)
 Inplace Multiply-subtract operation. Try to use the overloaded operators for vector instead, unless you want to fine-tune the number of GPU threads involved.
template<class SCALARTYPE , unsigned int ALIGNMENT>
void inplace_div_add (viennacl::vector< SCALARTYPE, ALIGNMENT > &vec1, const viennacl::vector< SCALARTYPE, ALIGNMENT > &vec2, scalar< SCALARTYPE > const &alpha)
 Inplace divide-add operation. Try to use the overloaded operators for vector instead, unless you want to fine-tune the number of GPU threads involved.
template<class SCALARTYPE , unsigned int ALIGNMENT>
void inplace_div_sub (viennacl::vector< SCALARTYPE, ALIGNMENT > &vec1, const viennacl::vector< SCALARTYPE, ALIGNMENT > &vec2, scalar< SCALARTYPE > const &alpha)
 Inplace divide-subtract operation. Try to use the overloaded operators for vector instead, unless you want to fine-tune the number of GPU threads involved.
template<class SCALARTYPE , unsigned int ALIGNMENT>
cl_uint index_norm_inf (const viennacl::vector< SCALARTYPE, ALIGNMENT > &vcl_vec)
 Computes the index of the first entry that is equal to the supremum-norm in modulus.
template<class SCALARTYPE , unsigned int ALIGNMENT>
void plane_rotation (viennacl::vector< SCALARTYPE, ALIGNMENT > &vec1, const viennacl::vector< SCALARTYPE, ALIGNMENT > &vec2, SCALARTYPE alpha, SCALARTYPE beta)
 Computes a plane rotation of two vectors.

Function Documentation

void viennacl::linalg::add ( const viennacl::matrix< TYPE, F, ALIGNMENT > &  mat1,
const viennacl::matrix< TYPE, F, ALIGNMENT > &  mat2,
viennacl::matrix< TYPE, F, ALIGNMENT > &  result 
)

Adds two dense matrices and writes the result to a third matrix.

This is the implementation of the convenience expression result = mat1 + mat2;

Parameters:
mat1 The left hand side operand
mat2 The right hand side operand
result The resulting matrix
void viennacl::linalg::add ( const viennacl::vector< SCALARTYPE, ALIGNMENT > &  vec1,
const viennacl::vector< SCALARTYPE, ALIGNMENT > &  vec2,
viennacl::vector< SCALARTYPE, ALIGNMENT > &  result 
)

Addition of two vectors. Try to use the overloaded operators for vector instead, unless you want to fine-tune the number of GPU threads involved.

Parameters:
vec1 The first addend.
vec2 The second addend.
result The result vector.
void viennacl::linalg::divide ( const viennacl::vector< SCALARTYPE, ALIGNMENT > &  vec,
scalar< SCALARTYPE > const &  alpha,
viennacl::vector< SCALARTYPE, ALIGNMENT > &  result 
)

Scales a vector. Try to use the overloaded operators for vector instead, unless you want to fine-tune the number of GPU threads involved.

Computes result = vec / alpha, where alpha is a gpu scalar

Parameters:
vec The vector to be scaled.
alpha The (inverse) scaling factor.
result The result vector.
void viennacl::linalg::ilu_inplace_solve ( MatrixType const &  mat,
VectorType &  vec,
viennacl::linalg::upper_tag   
)

Generic inplace solution of a upper triangular system.

Parameters:
mat The system matrix
vec The right hand side vector
void viennacl::linalg::ilu_inplace_solve ( MatrixType const &  mat,
VectorType &  vec,
viennacl::linalg::unit_lower_tag   
)

Generic inplace solution of a unit lower triangular system.

Parameters:
mat The system matrix
vec The right hand side vector
void viennacl::linalg::ilu_lu_substitute ( MatrixType const &  mat,
VectorType &  vec 
)

Generic LU substitution.

Parameters:
mat The system matrix
vec The right hand side vector
void viennacl::linalg::ilut_inc_row_iterator_to_row_index ( viennacl::tools::const_sparse_matrix_adapter< ScalarType > &  row_iter,
unsigned int  k 
)

Increments a row iterator (iteration along increasing row indices) up to a certain row index k.

Specialization for the const sparse matrix adapter shipped with ViennaCL

Parameters:
row_iter The row iterator
k The final row index
void viennacl::linalg::ilut_inc_row_iterator_to_row_index ( T &  row_iter,
unsigned int  k 
)

Increments a row iterator (iteration along increasing row indices) up to a certain row index k.

Generic implementation using the iterator concept from boost::numeric::ublas. Could not find a better way for sparse matrices...

Parameters:
row_iter The row iterator
k The final row index
void viennacl::linalg::ilut_inc_row_iterator_to_row_index ( viennacl::tools::sparse_matrix_adapter< ScalarType > &  row_iter,
unsigned int  k 
)

Increments a row iterator (iteration along increasing row indices) up to a certain row index k.

Specialization for the sparse matrix adapter shipped with ViennaCL

Parameters:
row_iter The row iterator
k The final row index
cl_uint viennacl::linalg::index_norm_inf ( const viennacl::vector< SCALARTYPE, ALIGNMENT > &  vcl_vec  ) 

Computes the index of the first entry that is equal to the supremum-norm in modulus.

Parameters:
vcl_vec The vector
Returns:
The result. Note that the result must be a CPU scalar (unsigned int), since gpu scalars are floating point types.
VectorT1::value_type viennacl::linalg::inner_prod ( VectorT1 const &  v1,
VectorT2 const &  v2,
typename viennacl::tools::enable_if< viennacl::is_stl< typename viennacl::traits::tag_of< VectorT1 >::type >::value >::type *  dummy = 0 
)
viennacl::scalar_expression< const viennacl::vector<ScalarType, alignment1>, const viennacl::vector<ScalarType, alignment2>, viennacl::op_inner_prod > viennacl::linalg::inner_prod ( viennacl::vector< ScalarType, alignment1 > const &  vector1,
viennacl::vector< ScalarType, alignment2 > const &  vector2,
typename viennacl::tools::enable_if< viennacl::is_viennacl< typename viennacl::traits::tag_of< viennacl::vector< ScalarType, alignment1 > >::type >::value >::type *  dummy = 0 
)
void inner_prod_impl ( const viennacl::vector< SCALARTYPE, ALIGNMENT > &  vec1,
const viennacl::vector< SCALARTYPE, ALIGNMENT > &  vec2,
scalar< SCALARTYPE > &  result 
)

Computes the inner product of two vectors - implementation. Library users should call inner_prod(vec1, vec2).

Parameters:
vec1 The first vector
vec2 The second vector
result The result scalar (on the gpu)
viennacl::scalar_expression< const viennacl::vector< SCALARTYPE, ALIGNMENT1 >, const viennacl::vector< SCALARTYPE, ALIGNMENT2 >, viennacl::op_inner_prod > inner_prod_impl ( const viennacl::vector< SCALARTYPE, ALIGNMENT1 > &  vec1,
const viennacl::vector< SCALARTYPE, ALIGNMENT2 > &  vec2 
)

Computes the inner product of two vectors.

Parameters:
vec1 The first vector
vec2 The second vector
Returns:
The result
void viennacl::linalg::inplace_add ( viennacl::vector< SCALARTYPE, ALIGNMENT > &  vec1,
const viennacl::vector< SCALARTYPE, ALIGNMENT > &  vec2 
)

Inplace addition of two vectors. Try to use the overloaded operators for vector instead, unless you want to fine-tune the number of GPU threads involved.

Computes vec1 += vec2.

Parameters:
vec1 The result.
vec2 The addend
void viennacl::linalg::inplace_add ( viennacl::matrix< TYPE, F, ALIGNMENT > &  result,
const viennacl::matrix< TYPE, F, ALIGNMENT > &  mat2 
)

Adds a dense matrix to another.

This is the implementation of the convenience expression result += mat1;

Parameters:
mat1 The left hand side operand
mat2 The right hand side operand
result The resulting matrix
void viennacl::linalg::inplace_div_add ( viennacl::vector< SCALARTYPE, ALIGNMENT > &  vec1,
const viennacl::vector< SCALARTYPE, ALIGNMENT > &  vec2,
scalar< SCALARTYPE > const &  alpha 
)

Inplace divide-add operation. Try to use the overloaded operators for vector instead, unless you want to fine-tune the number of GPU threads involved.

Computes vec1 += vec2 / alpha, where alpha is a gpu scalar

Parameters:
vec1 The first vector
vec2 The vector update
alpha The scaling factor for the second vector.
void viennacl::linalg::inplace_div_sub ( viennacl::vector< SCALARTYPE, ALIGNMENT > &  vec1,
const viennacl::vector< SCALARTYPE, ALIGNMENT > &  vec2,
scalar< SCALARTYPE > const &  alpha 
)

Inplace divide-subtract operation. Try to use the overloaded operators for vector instead, unless you want to fine-tune the number of GPU threads involved.

Computes vec1 -= vec2 / alpha, where alpha is a gpu scalar

Parameters:
vec1 The first vector
vec2 The vector update
alpha The scaling factor for the second vector.
void viennacl::linalg::inplace_divide ( viennacl::vector< SCALARTYPE, ALIGNMENT > &  vec,
scalar< SCALARTYPE > const &  alpha 
)

Scales a vector inplace. Try to use the overloaded operators for vector instead, unless you want to fine-tune the number of GPU threads involved.

Computes result *= alpha, where alpha is a gpu scalar

Parameters:
vec The vector to be scaled.
alpha The (inverse) scaling factor.
void viennacl::linalg::inplace_divide ( viennacl::matrix< SCALARTYPE, F, ALIGNMENT > &  result,
viennacl::scalar< SCALARTYPE > const &  val 
)

Multiplies a dense matrix by a scalar.

This is the implementation of the convenience expression result += mat1;

Parameters:
mat1 The left hand side operand
mat2 The right hand side operand
result The resulting matrix
void viennacl::linalg::inplace_mul_add ( viennacl::vector< SCALARTYPE, ALIGNMENT > &  vec1,
const viennacl::vector< SCALARTYPE, ALIGNMENT > &  vec2,
scalar< SCALARTYPE > const &  alpha 
)

Inplace Multiply-add operation. Try to use the overloaded operators for vector instead, unless you want to fine-tune the number of GPU threads involved.

Computes vec1 += alpha * vec2, where alpha is a gpu scalar

Parameters:
vec1 The first added
alpha The scaling factor for the first addend.
vec2 The second added.
void viennacl::linalg::inplace_mul_add ( viennacl::vector< SCALARTYPE, ALIGNMENT > &  vec1,
const viennacl::vector< SCALARTYPE, ALIGNMENT > &  vec2,
SCALARTYPE  alpha 
)

Inplace Multiply-add operation. Try to use the overloaded operators for vector instead, unless you want to fine-tune the number of GPU threads involved.

Computes vec1 += alpha * vec2, where alpha is a cpu scalar

Parameters:
vec1 The first added
vec2 The second added.
alpha The scaling factor for the first addend.
void viennacl::linalg::inplace_mul_sub ( viennacl::vector< SCALARTYPE, ALIGNMENT > &  vec1,
const viennacl::vector< SCALARTYPE, ALIGNMENT > &  vec2,
scalar< SCALARTYPE > const &  alpha 
)

Inplace Multiply-subtract operation. Try to use the overloaded operators for vector instead, unless you want to fine-tune the number of GPU threads involved.

Computes vec1 -= alpha * vec2, where alpha is a gpu scalar

Parameters:
vec1 The result vector which is updated
vec2 The second operand.
alpha The scaling factor for the vector update.
void viennacl::linalg::inplace_mult ( viennacl::vector< SCALARTYPE, ALIGNMENT > &  vec,
scalar< SCALARTYPE > const &  alpha 
)

Scales a vector inplace. Try to use the overloaded operators for vector instead, unless you want to fine-tune the number of GPU threads involved.

Computes result *= alpha, where alpha is a gpu scalar

Parameters:
vec The vector to be scaled.
alpha The scaling factor.
void viennacl::linalg::inplace_mult ( viennacl::vector< SCALARTYPE, ALIGNMENT > &  vec,
SCALARTYPE  alpha 
)

Scales a vector inplace. Try to use the overloaded operators for vector instead, unless you want to fine-tune the number of GPU threads involved.

Computes result *= alpha, where alpha is a cpu scalar

Parameters:
vec The vector to be scaled.
alpha The scaling factor.
void viennacl::linalg::inplace_mult ( viennacl::matrix< SCALARTYPE, F, ALIGNMENT > &  result,
viennacl::scalar< SCALARTYPE > const &  val 
)

Multiplies a dense matrix by a scalar.

This is the implementation of the convenience expression result += mat1;

Parameters:
mat1 The left hand side operand
mat2 The right hand side operand
result The resulting matrix
void viennacl::linalg::inplace_mult ( viennacl::matrix< SCALARTYPE, F, ALIGNMENT > &  result,
SCALARTYPE  val 
)

Multiplies a dense matrix by a scalar.

This is the implementation of the convenience expression result += mat1;

Parameters:
mat1 The left hand side operand
mat2 The right hand side operand
result The resulting matrix
void viennacl::linalg::inplace_solve ( compressed_matrix< SCALARTYPE, MAT_ALIGNMENT > const &  L,
vector< SCALARTYPE, VEC_ALIGNMENT > &  vec,
viennacl::linalg::unit_lower_tag   
)

Inplace solution of a lower triangular compressed_matrix with unit diagonal. Typically used for LU substitutions.

Parameters:
L The matrix
vec The vector
void viennacl::linalg::inplace_solve ( compressed_matrix< SCALARTYPE, MAT_ALIGNMENT > const &  U,
vector< SCALARTYPE, VEC_ALIGNMENT > &  vec,
viennacl::linalg::upper_tag   
)

Inplace solution of a upper triangular compressed_matrix. Typically used for LU substitutions.

Parameters:
U The upper triangular matrix
vec The vector
void viennacl::linalg::inplace_solve ( const matrix< SCALARTYPE, F1, A1 > &  mat,
matrix< SCALARTYPE, F2, A2 > &  B,
SOLVERTAG   
)

Direct inplace solver for dense upper triangular systems.

Parameters:
mat The system matrix
B The matrix of row vectors, where the solution is directly written to
void viennacl::linalg::inplace_solve ( const matrix< SCALARTYPE, F1, A1 > &  mat,
const matrix_expression< const matrix< SCALARTYPE, F2, A2 >, const matrix< SCALARTYPE, F2, A2 >, op_trans > &  B,
SOLVERTAG   
)

Direct inplace solver for dense upper triangular systems.

Parameters:
mat The system matrix
B The (transposed) matrix of row vectors, where the solution is directly written to
void viennacl::linalg::inplace_solve ( const matrix_expression< const matrix< SCALARTYPE, F1, A1 >, const matrix< SCALARTYPE, F1, A1 >, op_trans > &  proxy,
matrix< SCALARTYPE, F2, A2 > &  B,
SOLVERTAG   
)

Direct inplace solver for dense upper triangular systems that stem from transposed lower triangular systems.

Parameters:
proxy The system matrix proxy
B The matrix holding the load vectors, where the solution is directly written to
void viennacl::linalg::inplace_solve ( const matrix_expression< const matrix< SCALARTYPE, F1, A1 >, const matrix< SCALARTYPE, F1, A1 >, op_trans > &  proxy,
const matrix_expression< const matrix< SCALARTYPE, F2, A2 >, const matrix< SCALARTYPE, F2, A2 >, op_trans > &  B,
SOLVERTAG   
)

Direct inplace solver for dense upper triangular systems that stem from transposed lower triangular systems.

Parameters:
proxy The system matrix proxy
B The matrix holding the load vectors, where the solution is directly written to
void viennacl::linalg::inplace_solve ( const matrix< SCALARTYPE, F, ALIGNMENT > &  mat,
vector< SCALARTYPE, VEC_ALIGNMENT > &  vec,
SOLVERTAG   
)
void viennacl::linalg::inplace_solve ( const matrix_expression< const matrix< SCALARTYPE, F, ALIGNMENT >, const matrix< SCALARTYPE, F, ALIGNMENT >, op_trans > &  proxy,
vector< SCALARTYPE, VEC_ALIGNMENT > &  vec,
SOLVERTAG   
)

Direct inplace solver for dense upper triangular systems that stem from transposed lower triangular systems.

Parameters:
proxy The system matrix proxy
vec The load vector, where the solution is directly written to
void viennacl::linalg::inplace_sub ( viennacl::vector< SCALARTYPE, ALIGNMENT > &  vec1,
const viennacl::vector< SCALARTYPE, ALIGNMENT > &  vec2 
)

Inplace addition of two vectors. Try to use the overloaded operators for vector instead, unless you want to fine-tune the number of GPU threads involved.

Computes vec1 -= vec2.

Parameters:
vec1 The result.
vec2 The subtracted vector
void viennacl::linalg::inplace_sub ( viennacl::matrix< TYPE, F, ALIGNMENT > &  result,
const viennacl::matrix< TYPE, F, ALIGNMENT > &  mat2 
)

Adds a dense matrix to another.

This is the implementation of the convenience expression result += mat1;

Parameters:
mat1 The left hand side operand
mat2 The right hand side operand
result The resulting matrix
void viennacl::linalg::lu_factorize ( matrix< SCALARTYPE, F, ALIGNMENT > &  mat  ) 

LU factorization of a dense matrix.

Parameters:
mat The system matrix, where the LU matrices are directly written to. The implicit unit diagonal of L is not written.
void viennacl::linalg::lu_substitute ( matrix< SCALARTYPE, F1, ALIGNMENT_A > const &  A,
matrix< SCALARTYPE, F2, ALIGNMENT_B > &  B 
)

LU substitution for the system LU = rhs.

Parameters:
A The system matrix, where the LU matrices are directly written to. The implicit unit diagonal of L is not written.
B The matrix of load vectors, where the solution is directly written to
void viennacl::linalg::lu_substitute ( matrix< SCALARTYPE, F, ALIGNMENT > const &  mat,
vector< SCALARTYPE, VEC_ALIGNMENT > &  vec 
)

LU substitution for the system LU = rhs.

Parameters:
mat The system matrix, where the LU matrices are directly written to. The implicit unit diagonal of L is not written.
vec The load vector, where the solution is directly written to
void viennacl::linalg::mul_add ( const viennacl::vector< SCALARTYPE, ALIGNMENT > &  vec1,
scalar< SCALARTYPE > const &  alpha,
const viennacl::vector< SCALARTYPE, ALIGNMENT > &  vec2,
viennacl::vector< SCALARTYPE, ALIGNMENT > &  result 
)

Multiply-add operation. Try to use the overloaded operators for vector instead, unless you want to fine-tune the number of GPU threads involved.

Computes result = alpha * vec1 + vec2, where alpha is a gpu scalar

Parameters:
vec1 The first added
alpha The scaling factor for the first addend.
vec2 The second added.
result The result vector.
void viennacl::linalg::mul_add ( const viennacl::vector< SCALARTYPE, ALIGNMENT > &  vec1,
SCALARTYPE  alpha,
const viennacl::vector< SCALARTYPE, ALIGNMENT > &  vec2,
viennacl::vector< SCALARTYPE, ALIGNMENT > &  result 
)

Multiply-add operation. Try to use the overloaded operators for vector instead, unless you want to fine-tune the number of GPU threads involved.

Computes result = alpha * vec1 + vec2, where alpha is a cpu scalar

Parameters:
vec1 The first added
alpha The scaling factor for the first addend.
vec2 The second added.
result The result vector.
void viennacl::linalg::mul_sub ( const viennacl::vector< SCALARTYPE, ALIGNMENT > &  vec1,
scalar< SCALARTYPE > const &  alpha,
const viennacl::vector< SCALARTYPE, ALIGNMENT > &  vec2,
viennacl::vector< SCALARTYPE, ALIGNMENT > &  result 
)

Multiply-subtract operation. Try to use the overloaded operators for vector instead, unless you want to fine-tune the number of GPU threads involved.

Computes result = alpha * vec1 - vec2, where alpha is a gpu scalar

Parameters:
vec1 The first vector operand
alpha The scaling factor for the first vector.
vec2 The second operand.
result The result vector.
void viennacl::linalg::mult ( const viennacl::vector< SCALARTYPE, ALIGNMENT > &  vec,
scalar< SCALARTYPE > const &  alpha,
viennacl::vector< SCALARTYPE, ALIGNMENT > &  result 
)

Scales a vector. Try to use the overloaded operators for vector instead, unless you want to fine-tune the number of GPU threads involved.

Computes result = vec * alpha, where alpha is a gpu scalar

Parameters:
vec The vector to be scaled.
alpha The scaling factor.
result The result vector.
void viennacl::linalg::mult ( const viennacl::vector< SCALARTYPE, ALIGNMENT > &  vec,
SCALARTYPE  alpha,
viennacl::vector< SCALARTYPE, ALIGNMENT > &  result 
)

Scales a vector. Try to use the overloaded operators for vector instead, unless you want to fine-tune the number of GPU threads involved.

Computes result = vec * alpha, where alpha is a cpu scalar

Parameters:
vec The vector to be scaled.
alpha The scaling factor.
result The result vector.
VectorT::value_type viennacl::linalg::norm_1 ( VectorT const &  v1,
typename viennacl::tools::enable_if< viennacl::is_stl< typename viennacl::traits::tag_of< VectorT >::type >::value >::type *  dummy = 0 
)
viennacl::scalar_expression< const viennacl::vector<ScalarType, alignment>, const viennacl::vector<ScalarType, alignment>, viennacl::op_norm_1 > viennacl::linalg::norm_1 ( viennacl::vector< ScalarType, alignment > const &  vector,
typename viennacl::tools::enable_if< viennacl::is_viennacl< typename viennacl::traits::tag_of< viennacl::vector< ScalarType, alignment > >::type >::value >::type *  dummy = 0 
)
void norm_1_impl ( const viennacl::vector< SCALARTYPE, ALIGNMENT > &  vcl_vec,
scalar< SCALARTYPE > &  result 
)

Computes the l^1-norm of a vector.

Parameters:
vcl_vec The vector
result The result scalar
VectorT::value_type viennacl::linalg::norm_2 ( VectorT const &  v1,
typename viennacl::tools::enable_if< viennacl::is_stl< typename viennacl::traits::tag_of< VectorT >::type >::value >::type *  dummy = 0 
)
viennacl::scalar_expression< const viennacl::vector<ScalarType, alignment>, const viennacl::vector<ScalarType, alignment>, viennacl::op_norm_2 > viennacl::linalg::norm_2 ( viennacl::vector< ScalarType, alignment > const &  v,
typename viennacl::tools::enable_if< viennacl::is_viennacl< typename viennacl::traits::tag_of< viennacl::vector< ScalarType, alignment > >::type >::value >::type *  dummy = 0 
)
void norm_2_impl ( const viennacl::vector< SCALARTYPE, ALIGNMENT > &  vcl_vec,
scalar< SCALARTYPE > &  result 
)

Computes the l^2-norm of a vector - implementation.

Parameters:
vcl_vec The vector
result The result scalar
viennacl::scalar_expression< const viennacl::vector<ScalarType, alignment>, const viennacl::vector<ScalarType, alignment>, viennacl::op_norm_inf > viennacl::linalg::norm_inf ( viennacl::vector< ScalarType, alignment > const &  v1,
typename viennacl::tools::enable_if< viennacl::is_viennacl< typename viennacl::traits::tag_of< viennacl::vector< ScalarType, alignment > >::type >::value >::type *  dummy = 0 
)
VectorT::value_type viennacl::linalg::norm_inf ( VectorT const &  v1,
typename viennacl::tools::enable_if< viennacl::is_stl< typename viennacl::traits::tag_of< VectorT >::type >::value >::type *  dummy = 0 
)
void norm_inf_impl ( const viennacl::vector< SCALARTYPE, ALIGNMENT > &  vcl_vec,
scalar< SCALARTYPE > &  result 
)

Computes the supremum-norm of a vector.

Parameters:
vcl_vec The vector
result The result scalar
viennacl::matrix_expression< const viennacl::vector<SCALARTYPE, VA1>, const viennacl::vector<SCALARTYPE, VA2>, op_prod> viennacl::linalg::outer_prod ( const viennacl::vector< SCALARTYPE, VA1 > &  vec1,
const viennacl::vector< SCALARTYPE, VA2 > &  vec2 
)

Returns a proxy class for the operation mat += vec1 * vec2^T, i.e. a rank 1 update.

Parameters:
vec1 The first vector
vec2 The second vector
void viennacl::linalg::plane_rotation ( viennacl::vector< SCALARTYPE, ALIGNMENT > &  vec1,
const viennacl::vector< SCALARTYPE, ALIGNMENT > &  vec2,
SCALARTYPE  alpha,
SCALARTYPE  beta 
)

Computes a plane rotation of two vectors.

Computes (x,y) <- (alpha * x + beta * y, -beta * x + alpha * y)

Parameters:
vec1 The first vector
vec2 The second vector
alpha The first transformation coefficient
beta The second transformation coefficient
void viennacl::linalg::precondition ( MatrixType const &  input,
LUType &  output,
ilut_tag const &  tag 
)

Implementation of a ILU-preconditioner with threshold.

refer to Algorithm 10.6 by Saad's book (1996 edition)

Parameters:
input The input matrix. Type requirements: const_iterator1 for iteration along rows, const_iterator2 for iteration along columns
output The output matrix. Type requirements: const_iterator1 for iteration along rows, const_iterator2 for iteration along columns and write access via operator()
tag An ilut_tag in order to dispatch among several other preconditioners.
VectorT viennacl::linalg::prod ( MatrixT const &  matrix,
VectorT const &  vector,
typename viennacl::tools::enable_if< viennacl::is_stl< typename viennacl::traits::tag_of< MatrixT >::type >::value >::type *  dummy = 0 
)
viennacl::matrix_expression< const MatrixT, const viennacl::matrix<NumericT, F, ALIGNMENT>, viennacl::op_prod > viennacl::linalg::prod ( MatrixT const &  matrix_A,
viennacl::matrix< NumericT, F, ALIGNMENT > const &  matrix_B,
typename viennacl::tools::enable_if< viennacl::is_viennacl< typename viennacl::traits::tag_of< MatrixT >::type >::value >::type *  dummy = 0 
)
viennacl::vector_expression< const MatrixT, const viennacl::vector<NumericT, ALIGNMENT>, viennacl::op_prod > viennacl::linalg::prod ( MatrixT const &  matrix,
viennacl::vector< NumericT, ALIGNMENT > const &  vector,
typename viennacl::tools::enable_if< viennacl::is_viennacl< typename viennacl::traits::tag_of< MatrixT >::type >::value >::type *  dummy = 0 
)
viennacl::matrix_expression< const MatrixT, const viennacl::matrix_expression< const viennacl::matrix<NumericT, F, ALIGNMENT>, const viennacl::matrix<NumericT, F, ALIGNMENT>, viennacl::op_trans >, viennacl::op_prod > viennacl::linalg::prod ( MatrixT const &  matrix_A,
const viennacl::matrix_expression< const viennacl::matrix< NumericT, F, ALIGNMENT >, const viennacl::matrix< NumericT, F, ALIGNMENT >, viennacl::op_trans > &  matrix_B,
typename viennacl::tools::enable_if< viennacl::is_viennacl< typename viennacl::traits::tag_of< MatrixT >::type >::value >::type *  dummy = 0 
)
void viennacl::linalg::prod_impl ( const viennacl::matrix_expression< const matrix< TYPE, F1, ALIGNMENT >, const matrix< TYPE, F1, ALIGNMENT >, op_trans > &  A,
const viennacl::matrix_expression< const matrix< TYPE, F2, ALIGNMENT >, const matrix< TYPE, F2, ALIGNMENT >, op_trans > &  B,
viennacl::matrix< TYPE, F3, ALIGNMENT > &  C 
)

Carries out matrix-matrix multiplication.

Implementation of C = prod(trans(A), trans(B));

void viennacl::linalg::prod_impl ( const viennacl::matrix< TYPE, F, ALIGNMENT > &  mat,
const viennacl::vector< TYPE, VECTOR_ALIGNMENT > &  vec,
viennacl::vector< TYPE, VECTOR_ALIGNMENT > &  result 
)

Carries out matrix-vector multiplication.

Implementation of the convenience expression result = prod(mat, vec);

Parameters:
mat The matrix
vec The vector
result The result vector
VectorT viennacl::linalg::prod_impl ( std::vector< std::vector< T, A1 >, A2 > const &  matrix,
VectorT const &  vector 
)
void viennacl::linalg::prod_impl ( const viennacl::matrix< TYPE, F1, ALIGNMENT > &  A,
const viennacl::matrix_expression< const matrix< TYPE, F2, ALIGNMENT >, const matrix< TYPE, F2, ALIGNMENT >, op_trans > &  B,
viennacl::matrix< TYPE, F3, ALIGNMENT > &  C 
)

Carries out matrix-matrix multiplication.

Implementation of C = prod(A, trans(B));

void viennacl::linalg::prod_impl ( const viennacl::coordinate_matrix< TYPE, ALIGNMENT > &  mat,
const viennacl::vector< TYPE, VECTOR_ALIGNMENT > &  vec,
viennacl::vector< TYPE, VECTOR_ALIGNMENT > &  result 
)

Carries out matrix-vector multiplication with a coordinate_matrix.

Implementation of the convenience expression result = prod(mat, vec);

Parameters:
mat The matrix
vec The vector
result The result vector
vector_expression< const compressed_matrix< SCALARTYPE, ALIGNMENT >, const vector< SCALARTYPE, VECTOR_ALIGNMENT >, op_prod > prod_impl ( const compressed_matrix< SCALARTYPE, ALIGNMENT > &  mat,
const vector< SCALARTYPE, VECTOR_ALIGNMENT > &  vec 
)

Returns a proxy class that represents matrix-vector multiplication with a compressed_matrix.

This is used for the convenience expression result = prod(mat, vec);

Parameters:
mat The matrix
vec The vector
void viennacl::linalg::prod_impl ( const viennacl::compressed_matrix< TYPE, ALIGNMENT > &  mat,
const viennacl::vector< TYPE, VECTOR_ALIGNMENT > &  vec,
viennacl::vector< TYPE, VECTOR_ALIGNMENT > &  result,
size_t  NUM_THREADS = 0 
)

Carries out matrix-vector multiplication with a compressed_matrix.

Implementation of the convenience expression result = prod(mat, vec);

Parameters:
mat The matrix
vec The vector
result The result vector
NUM_THREADS Number of threads per work group. Can be used for fine-tuning.
VectorT viennacl::linalg::prod_impl ( std::vector< std::map< KEY, DATA, COMPARE, AMAP >, AVEC > const &  matrix,
VectorT const &  vector 
)
void viennacl::linalg::prod_impl ( const viennacl::matrix< TYPE, F1, ALIGNMENT > &  A,
const viennacl::matrix< TYPE, F2, ALIGNMENT > &  B,
viennacl::matrix< TYPE, F3, ALIGNMENT > &  C,
int  block_size = 15 
)

Carries out matrix-matrix multiplication.

Implementation of C = prod(A, B);

vector_expression< const coordinate_matrix< SCALARTYPE, ALIGNMENT >, const vector< SCALARTYPE, VECTOR_ALIGNMENT >, op_prod > prod_impl ( const coordinate_matrix< SCALARTYPE, ALIGNMENT > &  mat,
const vector< SCALARTYPE, VECTOR_ALIGNMENT > &  vec 
)

Returns a proxy class that represents matrix-vector multiplication with a compressed_matrix.

This is used for the convenience expression result = prod(mat, vec);

Parameters:
mat The matrix
vec The vector
viennacl::vector_expression< const viennacl::matrix< SCALARTYPE, F, ALIGNMENT >, const viennacl::vector< SCALARTYPE, VECTOR_ALIGNMENT >, op_prod > prod_impl ( const viennacl::matrix< SCALARTYPE, F, ALIGNMENT > &  mat,
const viennacl::vector< SCALARTYPE, VECTOR_ALIGNMENT > &  vec 
)

Returns a proxy class that represents matrix-vector multiplication.

This is used for the convenience expression result = prod(mat, vec);

Parameters:
mat The matrix
vec The vector
viennacl::vector_expression<const viennacl::coordinate_matrix<SCALARTYPE, ALIGNMENT>, const viennacl::vector<SCALARTYPE, VECTOR_ALIGNMENT>, viennacl::op_prod > viennacl::linalg::prod_impl ( const viennacl::coordinate_matrix< SCALARTYPE, ALIGNMENT > &  mat,
const viennacl::vector< SCALARTYPE, VECTOR_ALIGNMENT > &  vec,
size_t  NUM_THREADS 
)

Returns a proxy class that represents matrix-vector multiplication with a coordinate_matrix.

This is used for the convenience expression result = prod(mat, vec);

Parameters:
mat The matrix
vec The vector
NUM_THREADS Number of threads per work group. Can be used for fine-tuning.
viennacl::vector_expression<const viennacl::matrix_expression< const matrix<SCALARTYPE, F, ALIGNMENT>, const matrix<SCALARTYPE, F, ALIGNMENT>, op_trans>, const viennacl::vector<SCALARTYPE, VECTOR_ALIGNMENT>, op_prod > viennacl::linalg::prod_impl ( const viennacl::matrix_expression< const matrix< SCALARTYPE, F, ALIGNMENT >, const matrix< SCALARTYPE, F, ALIGNMENT >, op_trans > &  proxy,
const viennacl::vector< SCALARTYPE, VECTOR_ALIGNMENT > &  vec 
)

Returns a proxy class that represents matrix-vector multiplication with a transposed matrix.

This is used for the convenience expression result = trans(mat) * vec;

Parameters:
proxy The transposed matrix proxy
vec The vector
void viennacl::linalg::prod_impl ( const viennacl::matrix_expression< const matrix< SCALARTYPE, F, ALIGNMENT >, const matrix< SCALARTYPE, F, ALIGNMENT >, op_trans > &  mat,
const viennacl::vector< SCALARTYPE, VECTOR_ALIGNMENT > &  vec,
viennacl::vector< SCALARTYPE, VECTOR_ALIGNMENT > &  result 
)

Unwraps the transposed matrix proxy and forwards to trans_prod_impl().

void viennacl::linalg::prod_impl ( const viennacl::matrix_expression< const matrix< TYPE, F1, ALIGNMENT >, const matrix< TYPE, F1, ALIGNMENT >, op_trans > &  A,
const viennacl::matrix< TYPE, F2, ALIGNMENT > &  B,
viennacl::matrix< TYPE, F3, ALIGNMENT > &  C 
)

Carries out matrix-matrix multiplication.

Implementation of C = prod(trans(A), B);

void viennacl::linalg::rank_1_update ( viennacl::matrix< SCALARTYPE, F, ALIGNMENT > &  mat1,
const viennacl::vector< SCALARTYPE, ALIGNMENT > &  vec1,
const viennacl::vector< SCALARTYPE, ALIGNMENT > &  vec2 
)

The implementation of the operation mat += vec1 * vec2^T, i.e. a rank 1 update.

Implementation of the convenience expression result += outer_prod(vec1, vec2);

Parameters:
mat1 The matrix to be updated
vec1 The first vector
vec2 The second vector
void viennacl::linalg::scaled_rank_1_update ( viennacl::matrix< SCALARTYPE, F, ALIGNMENT > &  mat1,
SCALARTYPE  val,
const viennacl::vector< SCALARTYPE, ALIGNMENT > &  vec1,
const viennacl::vector< SCALARTYPE, ALIGNMENT > &  vec2 
)

The implementation of the operation mat += alpha * vec1 * vec2^T, i.e. a scaled rank 1 update.

Implementation of the convenience expression result += alpha * outer_prod(vec1, vec2);

Parameters:
mat1 The matrix to be updated
val The scaling factor
vec1 The first vector
vec2 The second vector
vector<SCALARTYPE, VEC_ALIGNMENT> viennacl::linalg::solve ( compressed_matrix< SCALARTYPE, MAT_ALIGNMENT > const &  L,
const vector< SCALARTYPE, VEC_ALIGNMENT > &  vec,
const viennacl::linalg::unit_lower_tag tag 
)

Convenience functions for result = solve(trans(mat), vec, unit_lower_tag()); Creates a temporary result vector and forwards the request to inplace_solve().

Parameters:
L The lower triangular sparse matrix
vec The load vector, where the solution is directly written to
tag Dispatch tag
VectorType viennacl::linalg::solve ( const MatrixType &  matrix,
VectorType const &  rhs,
gmres_tag const &  tag,
PreconditionerType const &  precond 
)

Implementation of the GMRES solver.

Following the algorithm proposed by Walker in "A Simpler GMRES"

Parameters:
matrix The system matrix
rhs The load vector
tag Solver configuration tag
precond A preconditioner. Precondition operation is done via member function apply()
Returns:
The result vector
VectorType viennacl::linalg::solve ( const MatrixType &  matrix,
VectorType const &  rhs,
cg_tag const &  tag 
)

Implementation of the conjugate gradient solver without preconditioner.

Following the algorithm in the book by Y. Saad "Iterative Methods for sparse linear systems"

Parameters:
matrix The system matrix
rhs The load vector
tag Solver configuration tag
Returns:
The result vector
matrix<SCALARTYPE, F2, ALIGNMENT_B> viennacl::linalg::solve ( const matrix< SCALARTYPE, F1, ALIGNMENT_A > &  A,
const matrix< SCALARTYPE, F2, ALIGNMENT_B > &  B,
TAG const &  tag 
)

Convenience functions for C = solve(A, B, some_tag()); Creates a temporary result matrix and forwards the request to inplace_solve().

Parameters:
A The system matrix
B The matrix of load vectors
tag Dispatch tag
VectorType viennacl::linalg::solve ( const MatrixType &  matrix,
VectorType const &  rhs,
cg_tag const &  tag,
viennacl::linalg::no_precond   
)
matrix<SCALARTYPE, F2, ALIGNMENT_B> viennacl::linalg::solve ( const matrix_expression< const matrix< SCALARTYPE, F1, ALIGNMENT_A >, const matrix< SCALARTYPE, F1, ALIGNMENT_A >, op_trans > &  proxy,
const matrix< SCALARTYPE, F2, ALIGNMENT_B > &  B,
TAG const &  tag 
)

Convenience functions for result = solve(trans(mat), B, some_tag()); Creates a temporary result matrix and forwards the request to inplace_solve().

Parameters:
proxy The transposed system matrix proxy
B The matrix of load vectors
tag Dispatch tag
vector<SCALARTYPE, VEC_ALIGNMENT> viennacl::linalg::solve ( const matrix_expression< const matrix< SCALARTYPE, F, ALIGNMENT >, const matrix< SCALARTYPE, F, ALIGNMENT >, op_trans > &  proxy,
const vector< SCALARTYPE, VEC_ALIGNMENT > &  vec,
TAG const &  tag 
)

Convenience functions for result = solve(trans(mat), vec, some_tag()); Creates a temporary result vector and forwards the request to inplace_solve().

Parameters:
proxy The transposed system matrix proxy
vec The load vector, where the solution is directly written to
tag Dispatch tag
VectorType viennacl::linalg::solve ( const MatrixType &  matrix,
VectorType const &  rhs,
bicgstab_tag const &  tag,
viennacl::linalg::no_precond   
)
VectorType viennacl::linalg::solve ( const MatrixType &  matrix,
VectorType const &  rhs,
bicgstab_tag const &  tag 
)

Implementation of the stabilized Bi-conjugate gradient solver.

Following the description in "Iterative Methods for Sparse Linear Systems" by Y. Saad

Parameters:
matrix The system matrix
rhs The load vector
tag Solver configuration tag
Returns:
The result vector
matrix<SCALARTYPE, F2, ALIGNMENT_B> viennacl::linalg::solve ( const matrix< SCALARTYPE, F1, ALIGNMENT_A > &  A,
const matrix_expression< const matrix< SCALARTYPE, F2, ALIGNMENT_B >, const matrix< SCALARTYPE, F2, ALIGNMENT_B >, op_trans > &  proxy,
TAG const &  tag 
)

Convenience functions for C = solve(A, B^T, some_tag()); Creates a temporary result matrix and forwards the request to inplace_solve().

Parameters:
A The system matrix
proxy The transposed load vector
tag Dispatch tag
VectorType viennacl::linalg::solve ( const MatrixType &  matrix,
VectorType const &  rhs,
bicgstab_tag const &  tag,
PreconditionerType const &  precond 
)

Implementation of the preconditioned stabilized Bi-conjugate gradient solver.

Following the description of the unpreconditioned case in "Iterative Methods for Sparse Linear Systems" by Y. Saad

Parameters:
matrix The system matrix
rhs The load vector
tag Solver configuration tag
precond A preconditioner. Precondition operation is done via member function apply()
Returns:
The result vector
VectorType viennacl::linalg::solve ( const MatrixType &  matrix,
VectorType const &  rhs,
cg_tag const &  tag,
PreconditionerType const &  precond 
)

Implementation of the preconditioned conjugate gradient solver.

Following Algorithm 9.1 in "Iterative Methods for Sparse Linear Systems" by Y. Saad

Parameters:
matrix The system matrix
rhs The load vector
tag Solver configuration tag
precond A preconditioner. Precondition operation is done via member function apply()
Returns:
The result vector
vector<SCALARTYPE, VEC_ALIGNMENT> viennacl::linalg::solve ( const matrix< SCALARTYPE, F, ALIGNMENT > &  mat,
const vector< SCALARTYPE, VEC_ALIGNMENT > &  vec,
TAG const &  tag 
)

Convenience functions for result = solve(mat, vec, some_tag()); Creates a temporary result vector and forwards the request to inplace_solve().

Parameters:
mat The system matrix
vec The load vector
tag Dispatch tag
vector<SCALARTYPE, VEC_ALIGNMENT> viennacl::linalg::solve ( compressed_matrix< SCALARTYPE, MAT_ALIGNMENT > const &  L,
const vector< SCALARTYPE, VEC_ALIGNMENT > &  vec,
viennacl::linalg::upper_tag const &  tag 
)

Convenience functions for result = solve(trans(mat), vec, unit_lower_tag()); Creates a temporary result vector and forwards the request to inplace_solve().

Parameters:
L The lower triangular sparse matrix
vec The load vector, where the solution is directly written to
tag Dispatch tag
matrix<SCALARTYPE, F2, ALIGNMENT_B> viennacl::linalg::solve ( const matrix_expression< const matrix< SCALARTYPE, F1, ALIGNMENT_A >, const matrix< SCALARTYPE, F1, ALIGNMENT_A >, op_trans > &  proxy_A,
const matrix_expression< const matrix< SCALARTYPE, F2, ALIGNMENT_B >, const matrix< SCALARTYPE, F2, ALIGNMENT_B >, op_trans > &  proxy_B,
TAG const &  tag 
)

Convenience functions for result = solve(trans(mat), vec, some_tag()); Creates a temporary result vector and forwards the request to inplace_solve().

Parameters:
proxy_A The transposed system matrix proxy
proxy_B The transposed matrix of load vectors, where the solution is directly written to
tag Dispatch tag
VectorType viennacl::linalg::solve ( const MatrixType &  matrix,
VectorType const &  rhs,
gmres_tag const &  tag 
)

Convenience overload of the solve() function using GMRES. Per default, no preconditioner is used.

void viennacl::linalg::sub ( const viennacl::matrix< TYPE, F, ALIGNMENT > &  mat1,
const viennacl::matrix< TYPE, F, ALIGNMENT > &  mat2,
viennacl::matrix< TYPE, F, ALIGNMENT > &  result 
)

Adds two dense matrices and writes the result to a third matrix.

This is the implementation of the convenience expression result = mat1 + mat2;

Parameters:
mat1 The left hand side operand
mat2 The right hand side operand
result The resulting matrix
void viennacl::linalg::sub ( const viennacl::vector< SCALARTYPE, ALIGNMENT > &  vec1,
const viennacl::vector< SCALARTYPE, ALIGNMENT > &  vec2,
viennacl::vector< SCALARTYPE, ALIGNMENT > &  result 
)

Subtraction of two vectors. Try to use the overloaded operators for vector instead, unless you want to fine-tune the number of GPU threads involved.

result = vec1 - vec2

Parameters:
vec1 The first operand.
vec2 The second operand.
result The result vector.
void viennacl::linalg::trans_prod_impl ( const matrix< SCALARTYPE, F, ALIGNMENT > &  mat,
const viennacl::vector< SCALARTYPE, VECTOR_ALIGNMENT > &  vec,
viennacl::vector< SCALARTYPE, VECTOR_ALIGNMENT > &  result 
)

Carries out matrix-vector multiplication with a transposed matrix.

Implementation of the convenience expression result = trans(mat) * vec;

Parameters:
mat The matrix
vec The vector
result The result vector