ViennaCL - The Vienna Computing Library  1.5.1
Namespaces | Data Structures | Typedefs | Functions | Variables
viennacl::linalg Namespace Reference

Provides all linear algebra operations which are not covered by operator overloads. More...

Namespaces

namespace  cuda
 

Holds all CUDA compute kernels used by ViennaCL.


namespace  detail
 

Namespace holding implementation details for linear algebra routines. Usually not of interest for a library user.


namespace  host_based
 

Holds all compute kernels with conventional host-based execution (buffers in CPU RAM).


namespace  kernels
 

Namespace containing the OpenCL kernels. Deprecated, will be moved to viennacl::linalg::opencl in future releases.


namespace  opencl
 

Holds all routines providing OpenCL linear algebra operations.


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  amg_precond
 AMG preconditioner class, can be supplied to solve()-routines. More...
class  amg_precond< compressed_matrix< ScalarType, MAT_ALIGNMENT > >
 AMG preconditioner class, can be supplied to solve()-routines. 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  block_ilu_precond
 A block ILU preconditioner class, can be supplied to solve()-routines. More...
class  block_ilu_precond< compressed_matrix< ScalarType, MAT_ALIGNMENT >, ILUTag >
 ILUT preconditioner class, can be supplied to solve()-routines. More...
class  ilu0_tag
 A tag for incomplete LU factorization with static pattern (ILU0) More...
class  ilu0_precond
 ILU0 preconditioner class, can be supplied to solve()-routines. More...
class  ilu0_precond< compressed_matrix< ScalarType, MAT_ALIGNMENT > >
 ILU0 preconditioner class, can be supplied to solve()-routines. 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  gmres_tag
 A tag for the solver GMRES. Used for supplying solver parameters and for dispatching the solve() function. More...
class  ichol0_tag
 A tag for incomplete Cholesky factorization with static pattern (ILU0) More...
class  ichol0_precond
 Incomplete Cholesky preconditioner class with static pattern (ICHOL0), can be supplied to solve()-routines. More...
class  ichol0_precond< compressed_matrix< ScalarType, MAT_ALIGNMENT > >
 ILU0 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. Generic version for non-ViennaCL matrices. More...
class  jacobi_precond< MatrixType, true >
 Jacobi preconditioner class, can be supplied to solve()-routines. More...
class  lanczos_tag
 A tag for the lanczos algorithm. More...
class  mixed_precision_cg_tag
 A tag for the conjugate gradient Used for supplying solver parameters and for dispatching the solve() function. More...
class  nmf_config
 Configuration class for the nonnegative-matrix-factorization algorithm. Specify tolerances, maximum iteration counts, etc., here. More...
class  power_iter_tag
 A tag for the power iteration algorithm. More...
class  row_scaling_tag
 A tag for a row scaling preconditioner which merely normalizes the equation system such that each row of the system matrix has unit norm. More...
class  row_scaling
 Jacobi-type preconditioner class, can be supplied to solve()-routines. This is a diagonal preconditioner with the diagonal entries being (configurable) row norms of the matrix. More...
class  row_scaling< MatrixType, true >
 Jacobi preconditioner class, can be supplied to solve()-routines. More...
class  spai_precond
 Implementation of the SParse Approximate Inverse Algorithm for a generic, uBLAS-compatible matrix type. More...
class  spai_precond< viennacl::compressed_matrix< ScalarType, MAT_ALIGNMENT > >
 Implementation of the SParse Approximate Inverse Algorithm for a ViennaCL compressed_matrix. More...
class  fspai_precond
 Implementation of the Factored SParse Approximate Inverse Algorithm for a generic, uBLAS-compatible matrix type. More...
class  fspai_precond< viennacl::compressed_matrix< ScalarType, MAT_ALIGNMENT > >
 Implementation of the Factored SParse Approximate Inverse Algorithm for a ViennaCL compressed_matrix. More...

Typedefs

typedef detail::amg::amg_tag amg_tag
typedef
viennacl::linalg::detail::spai::spai_tag 
spai_tag
typedef
viennacl::linalg::detail::spai::fspai_tag 
fspai_tag

Functions

template<class SCALARTYPE , unsigned int ALIGNMENT>
void convolve_i (viennacl::vector< SCALARTYPE, ALIGNMENT > &input1, viennacl::vector< SCALARTYPE, ALIGNMENT > &input2, viennacl::vector< SCALARTYPE, ALIGNMENT > &output)
template<typename T >
viennacl::vector_expression
< const vector_base< T >
, const vector_base< T >
, op_element_binary< op_prod > > 
element_prod (vector_base< T > const &v1, vector_base< T > const &v2)
template<typename T >
viennacl::vector_expression
< const vector_base< T >
, const vector_base< T >
, op_element_binary< op_div > > 
element_div (vector_base< T > const &v1, vector_base< T > const &v2)
template<typename T >
void inner_prod_impl (vector_base< T > const &vec1, vector_base< T > const &vec2, scalar< T > &result)
 Computes the inner product of two vectors - dispatcher interface.
template<typename LHS , typename RHS , typename OP , typename T >
void inner_prod_impl (viennacl::vector_expression< LHS, RHS, OP > const &vec1, vector_base< T > const &vec2, scalar< T > &result)
template<typename T , typename LHS , typename RHS , typename OP >
void inner_prod_impl (vector_base< T > const &vec1, viennacl::vector_expression< LHS, RHS, OP > const &vec2, scalar< T > &result)
template<typename LHS1 , typename RHS1 , typename OP1 , typename LHS2 , typename RHS2 , typename OP2 , typename T >
void inner_prod_impl (viennacl::vector_expression< LHS1, RHS1, OP1 > const &vec1, viennacl::vector_expression< LHS2, RHS2, OP2 > const &vec2, scalar< T > &result)
template<typename T >
void inner_prod_cpu (vector_base< T > const &vec1, vector_base< T > const &vec2, T &result)
 Computes the inner product of two vectors with the final reduction step on the CPU - dispatcher interface.
template<typename LHS , typename RHS , typename OP , typename T >
void inner_prod_cpu (viennacl::vector_expression< LHS, RHS, OP > const &vec1, vector_base< T > const &vec2, T &result)
template<typename T , typename LHS , typename RHS , typename OP >
void inner_prod_cpu (vector_base< T > const &vec1, viennacl::vector_expression< LHS, RHS, OP > const &vec2, T &result)
template<typename LHS1 , typename RHS1 , typename OP1 , typename LHS2 , typename RHS2 , typename OP2 , typename S3 >
void inner_prod_cpu (viennacl::vector_expression< LHS1, RHS1, OP1 > const &vec1, viennacl::vector_expression< LHS2, RHS2, OP2 > const &vec2, S3 &result)
template<typename T >
void norm_1_impl (vector_base< T > const &vec, scalar< T > &result)
 Computes the l^1-norm of a vector - dispatcher interface.
template<typename LHS , typename RHS , typename OP , typename T >
void norm_1_impl (viennacl::vector_expression< LHS, RHS, OP > const &vec, scalar< T > &result)
template<typename T >
void norm_1_cpu (vector_base< T > const &vec, T &result)
 Computes the l^1-norm of a vector with final reduction on the CPU.
template<typename LHS , typename RHS , typename OP , typename S2 >
void norm_1_cpu (viennacl::vector_expression< LHS, RHS, OP > const &vec, S2 &result)
 Computes the l^1-norm of a vector with final reduction on the CPU - interface for a vector expression. Creates a temporary.
template<typename T >
void norm_2_impl (vector_base< T > const &vec, scalar< T > &result)
 Computes the l^2-norm of a vector - dispatcher interface.
template<typename LHS , typename RHS , typename OP , typename T >
void norm_2_impl (viennacl::vector_expression< LHS, RHS, OP > const &vec, scalar< T > &result)
 Computes the l^2-norm of a vector - interface for a vector expression. Creates a temporary.
template<typename T >
void norm_2_cpu (vector_base< T > const &vec, T &result)
 Computes the l^2-norm of a vector with final reduction on the CPU - dispatcher interface.
template<typename LHS , typename RHS , typename OP , typename S2 >
void norm_2_cpu (viennacl::vector_expression< LHS, RHS, OP > const &vec, S2 &result)
 Computes the l^2-norm of a vector with final reduction on the CPU - interface for a vector expression. Creates a temporary.
template<typename T >
void norm_inf_impl (vector_base< T > const &vec, scalar< T > &result)
 Computes the supremum-norm of a vector.
template<typename LHS , typename RHS , typename OP , typename T >
void norm_inf_impl (viennacl::vector_expression< LHS, RHS, OP > const &vec, scalar< T > &result)
 Computes the supremum norm of a vector - interface for a vector expression. Creates a temporary.
template<typename T >
void norm_inf_cpu (vector_base< T > const &vec, T &result)
 Computes the supremum-norm of a vector with final reduction on the CPU.
template<typename LHS , typename RHS , typename OP , typename S2 >
void norm_inf_cpu (viennacl::vector_expression< LHS, RHS, OP > const &vec, S2 &result)
 Computes the supremum norm of a vector with final reduction on the CPU - interface for a vector expression. Creates a temporary.
template<typename T , typename F >
void norm_frobenius_impl (matrix_base< T, F > const &A, scalar< T > &result)
 Computes the Frobenius norm of a matrix - dispatcher interface.
template<typename T , typename F >
void norm_frobenius_cpu (matrix_base< T, F > const &A, T &result)
 Computes the Frobenius norm of a vector with final reduction on the CPU.
template<typename T >
vcl_size_t index_norm_inf (vector_base< T > const &vec)
 Computes the index of the first entry that is equal to the supremum-norm in modulus.
template<typename LHS , typename RHS , typename OP >
vcl_size_t index_norm_inf (viennacl::vector_expression< LHS, RHS, OP > const &vec)
 Computes the supremum norm of a vector with final reduction on the CPU - interface for a vector expression. Creates a temporary.
template<typename NumericT , typename F >
void prod_impl (const matrix_base< NumericT, F > &mat, const vector_base< NumericT > &vec, vector_base< NumericT > &result)
 Carries out matrix-vector multiplication.
template<typename NumericT , typename F >
void prod_impl (const matrix_expression< const matrix_base< NumericT, F >, const matrix_base< NumericT, F >, op_trans > &mat_trans, const vector_base< NumericT > &vec, vector_base< NumericT > &result)
 Carries out matrix-vector multiplication with a transposed matrix.
template<typename SparseMatrixType , class SCALARTYPE , unsigned int ALIGNMENT>
viennacl::enable_if
< viennacl::is_any_sparse_matrix
< SparseMatrixType >::value,
vector_expression< const
SparseMatrixType, const vector
< SCALARTYPE, ALIGNMENT >
, op_prod > >::type 
prod_impl (const SparseMatrixType &mat, const vector< SCALARTYPE, ALIGNMENT > &vec)
template<typename InternalType1 , typename InternalType2 >
void amg_setup (InternalType1 &A, InternalType1 &P, InternalType2 &Pointvector, amg_tag &tag)
 Setup AMG preconditioner.
template<typename MatrixType , typename InternalType1 , typename InternalType2 >
void amg_init (MatrixType const &mat, InternalType1 &A, InternalType1 &P, InternalType2 &Pointvector, amg_tag &tag)
 Initialize AMG preconditioner.
template<typename InternalType1 , typename InternalType2 >
void amg_transform_cpu (InternalType1 &A, InternalType1 &P, InternalType1 &R, InternalType2 &A_setup, InternalType2 &P_setup, amg_tag &tag)
 Save operators after setup phase for CPU computation.
template<typename InternalType1 , typename InternalType2 >
void amg_transform_gpu (InternalType1 &A, InternalType1 &P, InternalType1 &R, InternalType2 &A_setup, InternalType2 &P_setup, amg_tag &tag, viennacl::context ctx)
 Save operators after setup phase for GPU computation.
template<typename InternalVectorType , typename SparseMatrixType >
void amg_setup_apply (InternalVectorType &result, InternalVectorType &rhs, InternalVectorType &residual, SparseMatrixType const &A, amg_tag const &tag)
 Setup data structures for precondition phase.
template<typename InternalVectorType , typename SparseMatrixType >
void amg_setup_apply (InternalVectorType &result, InternalVectorType &rhs, InternalVectorType &residual, SparseMatrixType const &A, amg_tag const &tag, viennacl::context ctx)
 Setup data structures for precondition phase for later use on the GPU.
template<typename ScalarType , typename SparseMatrixType >
void amg_lu (boost::numeric::ublas::compressed_matrix< ScalarType > &op, boost::numeric::ublas::permutation_matrix<> &Permutation, SparseMatrixType const &A)
 Pre-compute LU factorization for direct solve (ublas library).
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 VectorT >
std::vector< typename
viennacl::result_of::cpu_value_type
< typename VectorT::value_type >
::type > 
bisect (VectorT const &alphas, VectorT const &betas)
 Implementation of the bisect-algorithm for the calculation of the eigenvalues of a tridiagonal matrix. Experimental - interface might change.
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 SCALARTYPE , unsigned int ALIGNMENT>
void prod_impl (const viennacl::circulant_matrix< SCALARTYPE, ALIGNMENT > &mat, const viennacl::vector_base< SCALARTYPE > &vec, viennacl::vector_base< SCALARTYPE > &result)
 Carries out matrix-vector multiplication with a circulant_matrix.
template<typename ScalarType >
void precondition (viennacl::compressed_matrix< ScalarType > &A, ilu0_tag const &)
 Implementation of a ILU-preconditioner with static pattern. Optimized version for CSR matrices.
template<typename ScalarType , typename SizeType , typename SparseVector >
ScalarType setup_w (viennacl::compressed_matrix< ScalarType > const &A, SizeType row, SparseVector &w)
 Dispatcher overload for extracting the row of nonzeros of a compressed matrix.
template<typename ScalarType , typename SizeType , typename SparseVector >
ScalarType setup_w (std::vector< std::map< SizeType, ScalarType > > const &A, SizeType row, SparseVector &w)
 Dispatcher overload for extracting the row of nonzeros of a STL-grown sparse matrix.
template<typename SparseMatrixType , typename ScalarType , typename SizeType >
void precondition (SparseMatrixType const &A, std::vector< std::map< SizeType, ScalarType > > &output, ilut_tag const &tag)
 Implementation of a ILU-preconditioner with threshold. Optimized implementation for compressed_matrix.
template<typename NumericT , typename F1 , typename F2 , typename SOLVERTAG >
void inplace_solve (const matrix_base< NumericT, F1 > &A, matrix_base< NumericT, F2 > &B, SOLVERTAG)
 Direct inplace solver for dense triangular systems. Matlab notation: A \ B.
template<typename NumericT , typename F1 , typename F2 , typename SOLVERTAG >
void inplace_solve (const matrix_base< NumericT, F1 > &A, matrix_expression< const matrix_base< NumericT, F2 >, const matrix_base< NumericT, F2 >, op_trans > proxy_B, SOLVERTAG)
 Direct inplace solver for dense triangular systems with transposed right hand side.
template<typename NumericT , typename F1 , typename F2 , typename SOLVERTAG >
void inplace_solve (const matrix_expression< const matrix_base< NumericT, F1 >, const matrix_base< NumericT, F1 >, op_trans > &proxy_A, matrix_base< NumericT, F2 > &B, SOLVERTAG)
 Direct inplace solver for dense triangular systems that stem from transposed triangular systems.
template<typename NumericT , typename F1 , typename F2 , typename SOLVERTAG >
void inplace_solve (const matrix_expression< const matrix_base< NumericT, F1 >, const matrix_base< NumericT, F1 >, op_trans > &proxy_A, matrix_expression< const matrix_base< NumericT, F2 >, const matrix_base< NumericT, F2 >, op_trans > proxy_B, SOLVERTAG)
 Direct inplace solver for dense transposed triangular systems with transposed right hand side. Matlab notation: A' \ B'.
template<typename NumericT , typename F , typename SOLVERTAG >
void inplace_solve (const matrix_base< NumericT, F > &mat, vector_base< NumericT > &vec, SOLVERTAG)
template<typename NumericT , typename F , typename SOLVERTAG >
void inplace_solve (const matrix_expression< const matrix_base< NumericT, F >, const matrix_base< NumericT, F >, op_trans > &proxy, vector_base< NumericT > &vec, SOLVERTAG)
 Direct inplace solver for dense upper triangular systems that stem from transposed lower triangular systems.
template<typename NumericT , typename F1 , typename F2 , typename SOLVERTAG >
matrix< NumericT, F2 > solve (const matrix_base< NumericT, F1 > &A, const matrix_base< NumericT, F2 > &B, SOLVERTAG tag)
 Convenience functions for C = solve(A, B, some_tag()); Creates a temporary result matrix and forwards the request to inplace_solve()
template<typename NumericT , typename F1 , typename F2 , typename SOLVERTAG >
matrix< NumericT, F2 > solve (const matrix_base< NumericT, F1 > &A, const matrix_expression< const matrix_base< NumericT, F2 >, const matrix_base< NumericT, F2 >, op_trans > &proxy, SOLVERTAG 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 NumericT , typename F1 , typename SOLVERTAG >
vector< NumericT > solve (const matrix_base< NumericT, F1 > &mat, const vector_base< NumericT > &vec, SOLVERTAG 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 NumericT , typename F1 , typename F2 , typename SOLVERTAG >
matrix< NumericT, F2 > solve (const matrix_expression< const matrix_base< NumericT, F1 >, const matrix_base< NumericT, F1 >, op_trans > &proxy, const matrix_base< NumericT, F2 > &B, SOLVERTAG 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 NumericT , typename F1 , typename F2 , typename SOLVERTAG >
matrix< NumericT, F2 > solve (const matrix_expression< const matrix_base< NumericT, F1 >, const matrix_base< NumericT, F1 >, op_trans > &proxy_A, const matrix_expression< const matrix_base< NumericT, F2 >, const matrix_base< NumericT, F2 >, op_trans > &proxy_B, SOLVERTAG 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 NumericT , typename F1 , typename SOLVERTAG >
vector< NumericT > solve (const matrix_expression< const matrix_base< NumericT, F1 >, const matrix_base< NumericT, F1 >, op_trans > &proxy, const vector_base< NumericT > &vec, SOLVERTAG 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 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<class SCALARTYPE , unsigned int ALIGNMENT>
void prod_impl (const viennacl::hankel_matrix< SCALARTYPE, ALIGNMENT > &mat, const viennacl::vector_base< SCALARTYPE > &vec, viennacl::vector_base< SCALARTYPE > &result)
 Carries out matrix-vector multiplication with a hankel_matrix.
template<typename ScalarType >
void precondition (viennacl::compressed_matrix< ScalarType > &A, ichol0_tag const &)
 Implementation of a ILU-preconditioner with static pattern. Optimized version for CSR matrices.
template<typename VectorT1 , typename VectorT2 >
viennacl::enable_if
< viennacl::is_stl< typename
viennacl::traits::tag_of
< VectorT1 >::type >::value,
typename VectorT1::value_type >
::type 
inner_prod (VectorT1 const &v1, VectorT2 const &v2)
template<typename NumericT >
viennacl::scalar_expression
< const vector_base< NumericT >
, const vector_base< NumericT >
, viennacl::op_inner_prod
inner_prod (vector_base< NumericT > const &vector1, vector_base< NumericT > const &vector2)
template<typename LHS , typename RHS , typename OP , typename NumericT >
viennacl::scalar_expression
< const
viennacl::vector_expression
< LHS, RHS, OP >, const
vector_base< NumericT >
, viennacl::op_inner_prod
inner_prod (viennacl::vector_expression< LHS, RHS, OP > const &vector1, vector_base< NumericT > const &vector2)
template<typename NumericT , typename LHS , typename RHS , typename OP >
viennacl::scalar_expression
< const vector_base< NumericT >
, const
viennacl::vector_expression
< LHS, RHS, OP >
, viennacl::op_inner_prod
inner_prod (vector_base< NumericT > const &vector1, viennacl::vector_expression< LHS, RHS, OP > const &vector2)
template<typename LHS1 , typename RHS1 , typename OP1 , typename LHS2 , typename RHS2 , typename OP2 >
viennacl::scalar_expression
< const
viennacl::vector_expression
< LHS1, RHS1, OP1 >, const
viennacl::vector_expression
< LHS2, RHS2, OP2 >
, viennacl::op_inner_prod
inner_prod (viennacl::vector_expression< LHS1, RHS1, OP1 > const &vector1, viennacl::vector_expression< LHS2, RHS2, OP2 > const &vector2)
template<typename NumericT >
viennacl::vector_expression
< const vector_base< NumericT >
, const vector_tuple< NumericT >
, viennacl::op_inner_prod
inner_prod (vector_base< NumericT > const &x, vector_tuple< NumericT > const &y_tuple)
template<typename MatrixT >
std::vector< typename
viennacl::result_of::cpu_value_type
< typename MatrixT::value_type >
::type > 
eig (MatrixT const &matrix, lanczos_tag const &tag)
 Implementation of the calculation of eigenvalues using lanczos.
template<typename SCALARTYPE >
void lu_factorize (matrix< SCALARTYPE, viennacl::row_major > &A)
 LU factorization of a row-major dense matrix.
template<typename SCALARTYPE >
void lu_factorize (matrix< SCALARTYPE, viennacl::column_major > &A)
 LU factorization of a column-major 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 &A, vector< SCALARTYPE, VEC_ALIGNMENT > &vec)
 LU substitution for the system LU = rhs.
template<typename NumericT , typename F , typename ScalarType1 >
void am (matrix_base< NumericT, F > &mat1, matrix_base< NumericT, F > const &mat2, ScalarType1 const &alpha, vcl_size_t len_alpha, bool reciprocal_alpha, bool flip_sign_alpha)
template<typename NumericT , typename F , typename ScalarType1 , typename ScalarType2 >
void ambm (matrix_base< NumericT, F > &mat1, matrix_base< NumericT, F > const &mat2, ScalarType1 const &alpha, vcl_size_t len_alpha, bool reciprocal_alpha, bool flip_sign_alpha, matrix_base< NumericT, F > const &mat3, ScalarType2 const &beta, vcl_size_t len_beta, bool reciprocal_beta, bool flip_sign_beta)
template<typename NumericT , typename F , typename ScalarType1 , typename ScalarType2 >
void ambm_m (matrix_base< NumericT, F > &mat1, matrix_base< NumericT, F > const &mat2, ScalarType1 const &alpha, vcl_size_t len_alpha, bool reciprocal_alpha, bool flip_sign_alpha, matrix_base< NumericT, F > const &mat3, ScalarType2 const &beta, vcl_size_t len_beta, bool reciprocal_beta, bool flip_sign_beta)
template<typename NumericT , typename F >
void matrix_assign (matrix_base< NumericT, F > &mat, NumericT s, bool clear=false)
template<typename NumericT , typename F >
void matrix_diagonal_assign (matrix_base< NumericT, F > &mat, NumericT s)
template<typename NumericT , typename F >
void matrix_diag_from_vector (const vector_base< NumericT > &v, int k, matrix_base< NumericT, F > &A)
 Dispatcher interface for A = diag(v, k)
template<typename NumericT , typename F >
void matrix_diag_to_vector (const matrix_base< NumericT, F > &A, int k, vector_base< NumericT > &v)
 Dispatcher interface for v = diag(A, k)
template<typename NumericT , typename F >
void matrix_row (const matrix_base< NumericT, F > &A, unsigned int i, vector_base< NumericT > &v)
template<typename NumericT , typename F >
void matrix_column (const matrix_base< NumericT, F > &A, unsigned int j, vector_base< NumericT > &v)
template<typename NumericT , typename F1 , typename F2 , typename F3 , typename ScalarType >
void prod_impl (const matrix_base< NumericT, F1 > &A, const matrix_base< NumericT, F2 > &B, matrix_base< NumericT, F3 > &C, ScalarType alpha, ScalarType beta)
 Carries out matrix-matrix multiplication.
template<typename NumericT , typename F1 , typename F2 , typename F3 , typename ScalarType >
void prod_impl (const viennacl::matrix_expression< const matrix_base< NumericT, F1 >, const matrix_base< NumericT, F1 >, op_trans > &A, const matrix_base< NumericT, F2 > &B, matrix_base< NumericT, F3 > &C, ScalarType alpha, ScalarType beta)
 Carries out matrix-matrix multiplication.
template<typename NumericT , typename F1 , typename F2 , typename F3 , typename ScalarType >
void prod_impl (const matrix_base< NumericT, F1 > &A, const viennacl::matrix_expression< const matrix_base< NumericT, F2 >, const matrix_base< NumericT, F2 >, op_trans > &B, matrix_base< NumericT, F3 > &C, ScalarType alpha, ScalarType beta)
 Carries out matrix-matrix multiplication.
template<typename NumericT , typename F1 , typename F2 , typename F3 , typename ScalarType >
void prod_impl (const viennacl::matrix_expression< const matrix_base< NumericT, F1 >, const matrix_base< NumericT, F1 >, op_trans > &A, const viennacl::matrix_expression< const matrix_base< NumericT, F2 >, const matrix_base< NumericT, F2 >, op_trans > &B, matrix_base< NumericT, F3 > &C, ScalarType alpha, ScalarType beta)
 Carries out matrix-matrix multiplication.
template<typename T , typename F , typename OP >
void element_op (matrix_base< T, F > &A, matrix_expression< const matrix_base< T, F >, const matrix_base< T, F >, OP > const &proxy)
 Implementation of the element-wise operation A = B .* C and A = B ./ C for matrices (using MATLAB syntax). Don't use this function directly, use element_prod() and element_div().
template<typename NumericT >
viennacl::matrix_expression
< const vector_base< NumericT >
, const vector_base< NumericT >
, op_prod
outer_prod (const vector_base< NumericT > &vec1, const vector_base< NumericT > &vec2)
 Returns a proxy class for the operation mat += vec1 * vec2^T, i.e. a rank 1 update.
template<typename NumericT , typename F , typename S1 >
void scaled_rank_1_update (matrix_base< NumericT, F > &mat1, S1 const &alpha, vcl_size_t len_alpha, bool reciprocal_alpha, bool flip_sign_alpha, const vector_base< NumericT > &vec1, const vector_base< NumericT > &vec2)
 The implementation of the operation mat += alpha * vec1 * vec2^T, i.e. a scaled rank 1 update.
template<typename MatrixType , typename VectorType >
VectorType solve (const MatrixType &matrix, VectorType const &rhs, mixed_precision_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, mixed_precision_cg_tag const &tag, viennacl::linalg::no_precond)
template<typename ScalarType >
void nmf (viennacl::matrix< ScalarType > const &V, viennacl::matrix< ScalarType > &W, viennacl::matrix< ScalarType > &H, nmf_config const &conf)
 The nonnegative matrix factorization (approximation) algorithm as suggested by Lee and Seung. Factorizes a matrix V with nonnegative entries into matrices W and H such that ||V - W*H|| is minimized.
template<typename T , typename A >
norm_1 (std::vector< T, A > const &v1)
template<typename ScalarType >
viennacl::scalar_expression
< const viennacl::vector_base
< ScalarType >, const
viennacl::vector_base
< ScalarType >
, viennacl::op_norm_1
norm_1 (viennacl::vector_base< ScalarType > const &vector)
template<typename LHS , typename RHS , typename OP >
viennacl::scalar_expression
< const
viennacl::vector_expression
< const LHS, const RHS, OP >
, const
viennacl::vector_expression
< const LHS, const RHS, OP >
, viennacl::op_norm_1
norm_1 (viennacl::vector_expression< const LHS, const RHS, OP > const &vector)
template<typename T , typename A >
norm_2 (std::vector< T, A > const &v1)
template<typename ScalarType >
viennacl::scalar_expression
< const viennacl::vector_base
< ScalarType >, const
viennacl::vector_base
< ScalarType >
, viennacl::op_norm_2
norm_2 (viennacl::vector_base< ScalarType > const &v)
template<typename LHS , typename RHS , typename OP >
viennacl::scalar_expression
< const
viennacl::vector_expression
< const LHS, const RHS, OP >
, const
viennacl::vector_expression
< const LHS, const RHS, OP >
, viennacl::op_norm_2
norm_2 (viennacl::vector_expression< const LHS, const RHS, OP > const &vector)
template<typename NumericT , typename F >
scalar_expression< const
matrix_base< NumericT, F >
, const matrix_base< NumericT,
F >, op_norm_frobenius
norm_frobenius (const matrix< NumericT, F > &A)
template<typename T , typename A >
norm_inf (std::vector< T, A > const &v1)
template<typename ScalarType >
viennacl::scalar_expression
< const viennacl::vector_base
< ScalarType >, const
viennacl::vector_base
< ScalarType >
, viennacl::op_norm_inf
norm_inf (viennacl::vector_base< ScalarType > const &v1)
template<typename LHS , typename RHS , typename OP >
viennacl::scalar_expression
< const
viennacl::vector_expression
< const LHS, const RHS, OP >
, const
viennacl::vector_expression
< const LHS, const RHS, OP >
, viennacl::op_norm_inf
norm_inf (viennacl::vector_expression< const LHS, const RHS, OP > const &vector)
template<typename MatrixT >
viennacl::result_of::cpu_value_type
< typename MatrixT::value_type >
::type 
eig (MatrixT const &matrix, power_iter_tag const &tag)
 Implementation of the calculation of eigenvalues using poweriteration.
template<typename T , typename A1 , typename A2 , typename VectorT >
VectorT prod (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 (std::vector< std::map< KEY, DATA, COMPARE, AMAP >, AVEC > const &matrix, VectorT const &vector)
template<typename NumericT , typename F1 , typename F2 >
viennacl::matrix_expression
< const viennacl::matrix_base
< NumericT, F1 >, const
viennacl::matrix_base
< NumericT, F2 >
, viennacl::op_mat_mat_prod
prod (viennacl::matrix_base< NumericT, F1 > const &A, viennacl::matrix_base< NumericT, F2 > const &B)
template<typename NumericT , typename F1 , typename F2 >
viennacl::matrix_expression
< const viennacl::matrix_base
< NumericT, F1 >, const
viennacl::matrix_expression
< const viennacl::matrix_base
< NumericT, F2 >, const
viennacl::matrix_base
< NumericT, F2 >, op_trans >
, viennacl::op_mat_mat_prod
prod (viennacl::matrix_base< NumericT, F1 > const &A, viennacl::matrix_expression< const viennacl::matrix_base< NumericT, F2 >, const viennacl::matrix_base< NumericT, F2 >, op_trans > const &B)
template<typename NumericT , typename F1 , typename F2 >
viennacl::matrix_expression
< const
viennacl::matrix_expression
< const viennacl::matrix_base
< NumericT, F1 >, const
viennacl::matrix_base
< NumericT, F1 >, op_trans >
, const viennacl::matrix_base
< NumericT, F2 >
, viennacl::op_mat_mat_prod
prod (viennacl::matrix_expression< const viennacl::matrix_base< NumericT, F1 >, const viennacl::matrix_base< NumericT, F1 >, op_trans > const &A, viennacl::matrix_base< NumericT, F2 > const &B)
template<typename NumericT , typename F1 , typename F2 >
viennacl::matrix_expression
< const
viennacl::matrix_expression
< const viennacl::matrix_base
< NumericT, F1 >, const
viennacl::matrix_base
< NumericT, F1 >, op_trans >
, const
viennacl::matrix_expression
< const viennacl::matrix_base
< NumericT, F2 >, const
viennacl::matrix_base
< NumericT, F2 >, op_trans >
, viennacl::op_mat_mat_prod
prod (viennacl::matrix_expression< const viennacl::matrix_base< NumericT, F1 >, const viennacl::matrix_base< NumericT, F1 >, op_trans > const &A, viennacl::matrix_expression< const viennacl::matrix_base< NumericT, F2 >, const viennacl::matrix_base< NumericT, F2 >, op_trans > const &B)
template<typename NumericT , typename F >
viennacl::vector_expression
< const viennacl::matrix_base
< NumericT, F >, const
viennacl::vector_base
< NumericT >
, viennacl::op_prod
prod (viennacl::matrix_base< NumericT, F > const &matrix, viennacl::vector_base< NumericT > const &vector)
template<typename NumericT , typename F >
viennacl::vector_expression
< const
viennacl::matrix_expression
< const viennacl::matrix_base
< NumericT, F >, const
viennacl::matrix_base
< NumericT, F >, op_trans >
, const viennacl::vector_base
< NumericT >
, viennacl::op_prod
prod (viennacl::matrix_expression< const viennacl::matrix_base< NumericT, F >, const viennacl::matrix_base< NumericT, F >, op_trans > const &matrix, viennacl::vector_base< NumericT > const &vector)
template<typename SparseMatrixType , class SCALARTYPE >
viennacl::enable_if
< viennacl::is_any_sparse_matrix
< SparseMatrixType >::value,
vector_expression< const
SparseMatrixType, const
vector_base< SCALARTYPE >
, op_prod > >::type 
prod (const SparseMatrixType &mat, const vector_base< SCALARTYPE > &vec)
template<typename SparseMatrixType , typename SCALARTYPE , typename F1 >
viennacl::enable_if
< viennacl::is_any_sparse_matrix
< SparseMatrixType >::value,
viennacl::matrix_expression
< const SparseMatrixType,
const matrix_base< SCALARTYPE,
F1 >, op_prod > >::type 
prod (const SparseMatrixType &sp_mat, const viennacl::matrix_base< SCALARTYPE, F1 > &d_mat)
template<typename SparseMatrixType , typename SCALARTYPE , typename F1 >
viennacl::enable_if
< viennacl::is_any_sparse_matrix
< SparseMatrixType >::value,
viennacl::matrix_expression
< const SparseMatrixType,
const
viennacl::matrix_expression
< const viennacl::matrix_base
< SCALARTYPE, F1 >, const
viennacl::matrix_base
< SCALARTYPE, F1 >, op_trans >
, viennacl::op_prod > >::type 
prod (const SparseMatrixType &A, viennacl::matrix_expression< const viennacl::matrix_base< SCALARTYPE, F1 >, const viennacl::matrix_base< SCALARTYPE, F1 >, op_trans > const &B)
template<typename StructuredMatrixType , class SCALARTYPE >
viennacl::enable_if
< viennacl::is_any_dense_structured_matrix
< StructuredMatrixType >
::value, vector_expression
< const StructuredMatrixType,
const vector_base< SCALARTYPE >
, op_prod > >::type 
prod (const StructuredMatrixType &mat, const vector_base< SCALARTYPE > &vec)
template<typename SCALARTYPE , typename F , unsigned int ALIGNMENT>
void qr_method_nsm (viennacl::matrix< SCALARTYPE, F, ALIGNMENT > &A, viennacl::matrix< SCALARTYPE, F, ALIGNMENT > &Q, boost::numeric::ublas::vector< SCALARTYPE > &D, boost::numeric::ublas::vector< SCALARTYPE > &E)
template<typename SCALARTYPE , typename F , unsigned int ALIGNMENT>
void qr_method_sym (viennacl::matrix< SCALARTYPE, F, ALIGNMENT > &A, viennacl::matrix< SCALARTYPE, F, ALIGNMENT > &Q, boost::numeric::ublas::vector< SCALARTYPE > &D)
template<typename MatrixType , typename VectorType >
void recoverQ (MatrixType const &A, VectorType const &betas, MatrixType &Q, MatrixType &R)
template<typename MatrixType , typename VectorType1 , typename VectorType2 >
void inplace_qr_apply_trans_Q (MatrixType const &A, VectorType1 const &betas, VectorType2 &b)
 Computes Q^T b, where Q is an implicit orthogonal matrix defined via its Householder reflectors stored in A.
template<typename T , typename F , unsigned int ALIGNMENT, typename VectorType1 , unsigned int A2>
void inplace_qr_apply_trans_Q (viennacl::matrix< T, F, ALIGNMENT > const &A, VectorType1 const &betas, viennacl::vector< T, A2 > &b)
template<typename T , typename F , unsigned int ALIGNMENT>
std::vector< T > inplace_qr (viennacl::matrix< T, F, ALIGNMENT > &A, vcl_size_t block_size=16)
 Overload of inplace-QR factorization of a ViennaCL matrix A.
template<typename MatrixType >
std::vector< typename
MatrixType::value_type > 
inplace_qr (MatrixType &A, vcl_size_t block_size=16)
 Overload of inplace-QR factorization for a general Boost.uBLAS compatible matrix A.
template<typename S1 , typename S2 , typename ScalarType1 >
viennacl::enable_if
< viennacl::is_scalar< S1 >
::value &&viennacl::is_scalar
< S2 >::value
&&viennacl::is_any_scalar
< ScalarType1 >::value >::type 
as (S1 &s1, S2 const &s2, ScalarType1 const &alpha, vcl_size_t len_alpha, bool reciprocal_alpha, bool flip_sign_alpha)
 Interface for the generic operation s1 = s2 @ alpha, where s1 and s2 are GPU scalars, @ denotes multiplication or division, and alpha is either a GPU or a CPU scalar.
template<typename S1 , typename S2 , typename ScalarType1 , typename S3 , typename ScalarType2 >
viennacl::enable_if
< viennacl::is_scalar< S1 >
::value &&viennacl::is_scalar
< S2 >::value
&&viennacl::is_scalar< S3 >
::value
&&viennacl::is_any_scalar
< ScalarType1 >::value
&&viennacl::is_any_scalar
< ScalarType2 >::value >::type 
asbs (S1 &s1, S2 const &s2, ScalarType1 const &alpha, vcl_size_t len_alpha, bool reciprocal_alpha, bool flip_sign_alpha, S3 const &s3, ScalarType2 const &beta, vcl_size_t len_beta, bool reciprocal_beta, bool flip_sign_beta)
 Interface for the generic operation s1 = s2 @ alpha + s3 @ beta, where s1, s2 and s3 are GPU scalars, @ denotes multiplication or division, and alpha, beta are either a GPU or a CPU scalar.
template<typename S1 , typename S2 , typename ScalarType1 , typename S3 , typename ScalarType2 >
viennacl::enable_if
< viennacl::is_scalar< S1 >
::value &&viennacl::is_scalar
< S2 >::value
&&viennacl::is_scalar< S3 >
::value
&&viennacl::is_any_scalar
< ScalarType1 >::value
&&viennacl::is_any_scalar
< ScalarType2 >::value >::type 
asbs_s (S1 &s1, S2 const &s2, ScalarType1 const &alpha, vcl_size_t len_alpha, bool reciprocal_alpha, bool flip_sign_alpha, S3 const &s3, ScalarType2 const &beta, vcl_size_t len_beta, bool reciprocal_beta, bool flip_sign_beta)
 Interface for the generic operation s1 += s2 @ alpha + s3 @ beta, where s1, s2 and s3 are GPU scalars, @ denotes multiplication or division, and alpha, beta are either a GPU or a CPU scalar.
template<typename S1 , typename S2 >
viennacl::enable_if
< viennacl::is_scalar< S1 >
::value &&viennacl::is_scalar
< S2 >::value >::type 
swap (S1 &s1, S2 &s2)
 Swaps the contents of two scalars.
template<typename SparseMatrixType , class ScalarType >
viennacl::enable_if
< viennacl::is_any_sparse_matrix
< SparseMatrixType >::value >
::type 
prod_impl (const SparseMatrixType &mat, const viennacl::vector_base< ScalarType > &vec, viennacl::vector_base< ScalarType > &result)
 Carries out matrix-vector multiplication involving a sparse matrix type.
template<typename SparseMatrixType , class ScalarType , typename F1 , typename F2 >
viennacl::enable_if
< viennacl::is_any_sparse_matrix
< SparseMatrixType >::value >
::type 
prod_impl (const SparseMatrixType &sp_mat, const viennacl::matrix_base< ScalarType, F1 > &d_mat, viennacl::matrix_base< ScalarType, F2 > &result)
 Carries out matrix-matrix multiplication first matrix being sparse.
template<typename SparseMatrixType , class ScalarType , typename F1 , typename F2 >
viennacl::enable_if
< viennacl::is_any_sparse_matrix
< SparseMatrixType >::value >
::type 
prod_impl (const SparseMatrixType &sp_mat, const viennacl::matrix_expression< const viennacl::matrix_base< ScalarType, F1 >, const viennacl::matrix_base< ScalarType, F1 >, viennacl::op_trans > &d_mat, viennacl::matrix_base< ScalarType, F2 > &result)
 Carries out matrix-matrix multiplication first matrix being sparse, and the second transposed.
template<typename SparseMatrixType , class ScalarType , typename SOLVERTAG >
viennacl::enable_if
< viennacl::is_any_sparse_matrix
< SparseMatrixType >::value >
::type 
inplace_solve (const SparseMatrixType &mat, viennacl::vector_base< ScalarType > &vec, SOLVERTAG tag)
 Carries out triangular inplace solves.
template<typename SparseMatrixType , class ScalarType , typename SOLVERTAG >
viennacl::enable_if
< viennacl::is_any_sparse_matrix
< SparseMatrixType >::value >
::type 
inplace_solve (const matrix_expression< const SparseMatrixType, const SparseMatrixType, op_trans > &mat, viennacl::vector_base< ScalarType > &vec, SOLVERTAG tag)
 Carries out transposed triangular inplace solves.
template<typename SCALARTYPE , unsigned int ALIGNMENT>
void svd (viennacl::matrix< SCALARTYPE, row_major, ALIGNMENT > &A, viennacl::matrix< SCALARTYPE, row_major, ALIGNMENT > &QL, viennacl::matrix< SCALARTYPE, row_major, ALIGNMENT > &QR)
 Computes the singular value decomposition of a matrix A. Experimental in 1.3.x.
template<class SCALARTYPE , unsigned int ALIGNMENT>
void prod_impl (const viennacl::toeplitz_matrix< SCALARTYPE, ALIGNMENT > &mat, const viennacl::vector_base< SCALARTYPE > &vec, viennacl::vector_base< SCALARTYPE > &result)
 Carries out matrix-vector multiplication with a toeplitz_matrix.
template<typename ScalarType >
void inplace_tred2 (boost::numeric::ublas::matrix< ScalarType > const &A, vcl_size_t block_size=1)
 Inplace reduction of a hermitian (or real symmetric) to tridiagonal form using householder similarity transforms (preserving eigenvalues)
template<class SCALARTYPE , unsigned int ALIGNMENT>
void prod_impl (const viennacl::vandermonde_matrix< SCALARTYPE, ALIGNMENT > &mat, const viennacl::vector_base< SCALARTYPE > &vec, viennacl::vector_base< SCALARTYPE > &result)
 Carries out matrix-vector multiplication with a vandermonde_matrix.
template<typename T , typename ScalarType1 >
void av (vector_base< T > &vec1, vector_base< T > const &vec2, ScalarType1 const &alpha, vcl_size_t len_alpha, bool reciprocal_alpha, bool flip_sign_alpha)
template<typename T , typename ScalarType1 , typename ScalarType2 >
void avbv (vector_base< T > &vec1, vector_base< T > const &vec2, ScalarType1 const &alpha, vcl_size_t len_alpha, bool reciprocal_alpha, bool flip_sign_alpha, vector_base< T > const &vec3, ScalarType2 const &beta, vcl_size_t len_beta, bool reciprocal_beta, bool flip_sign_beta)
template<typename T , typename ScalarType1 , typename ScalarType2 >
void avbv_v (vector_base< T > &vec1, vector_base< T > const &vec2, ScalarType1 const &alpha, vcl_size_t len_alpha, bool reciprocal_alpha, bool flip_sign_alpha, vector_base< T > const &vec3, ScalarType2 const &beta, vcl_size_t len_beta, bool reciprocal_beta, bool flip_sign_beta)
template<typename T >
void vector_assign (vector_base< T > &vec1, const T &alpha, bool up_to_internal_size=false)
 Assign a constant value to a vector (-range/-slice)
template<typename T >
void vector_swap (vector_base< T > &vec1, vector_base< T > &vec2)
 Swaps the contents of two vectors, data is copied.
template<typename T , typename OP >
void element_op (vector_base< T > &vec1, vector_expression< const vector_base< T >, const vector_base< T >, OP > const &proxy)
 Implementation of the element-wise operation v1 = v2 .* v3 and v1 = v2 ./ v3 (using MATLAB syntax)
template<typename T >
void inner_prod_impl (vector_base< T > const &x, vector_tuple< T > const &y_tuple, vector_base< T > &result)
 Computes the inner products <x, y1>, <x, y2>, ..., <x, y_N> and writes the result to a (sub-)vector.
template<typename LHS , typename RHS , typename OP , typename S2 >
void norm_1_impl (viennacl::vector_expression< LHS, RHS, OP > const &vec, S2 &result)
 Computes the l^1-norm of a vector - interface for a vector expression. Creates a temporary.
template<typename T >
void plane_rotation (vector_base< T > &vec1, vector_base< T > &vec2, T alpha, T beta)
 Computes a plane rotation of two vectors.

Variables

const char * double_float_conversion_program = "};\n"
const std::string SVD_BIDIAG_PACK_KERNEL = "bidiag_pack"
const std::string SVD_HOUSEHOLDER_UPDATE_A_LEFT_KERNEL = "house_update_A_left"
const std::string SVD_HOUSEHOLDER_UPDATE_A_RIGHT_KERNEL = "house_update_A_right"
const std::string SVD_HOUSEHOLDER_UPDATE_QL_KERNEL = "house_update_QL"
const std::string SVD_HOUSEHOLDER_UPDATE_QR_KERNEL = "house_update_QR"
const std::string SVD_COPY_COL_KERNEL = "copy_col"
const std::string SVD_COPY_ROW_KERNEL = "copy_row"
const std::string SVD_MATRIX_TRANSPOSE_KERNEL = "transpose_inplace"
const std::string SVD_INVERSE_SIGNS_KERNEL = "inverse_signs"
const std::string SVD_GIVENS_PREV_KERNEL = "givens_prev"
const std::string SVD_GIVENS_NEXT_KERNEL = "givens_next"
const std::string SVD_FINAL_ITER_UPDATE_KERNEL = "final_iter_update"
const std::string SVD_UPDATE_QR_COLUMN_KERNEL = "update_qr_column"

Detailed Description

Provides all linear algebra operations which are not covered by operator overloads.


Typedef Documentation


Function Documentation

void viennacl::linalg::am ( matrix_base< NumericT, F > &  mat1,
matrix_base< NumericT, F > const &  mat2,
ScalarType1 const &  alpha,
vcl_size_t  len_alpha,
bool  reciprocal_alpha,
bool  flip_sign_alpha 
)
void viennacl::linalg::ambm ( matrix_base< NumericT, F > &  mat1,
matrix_base< NumericT, F > const &  mat2,
ScalarType1 const &  alpha,
vcl_size_t  len_alpha,
bool  reciprocal_alpha,
bool  flip_sign_alpha,
matrix_base< NumericT, F > const &  mat3,
ScalarType2 const &  beta,
vcl_size_t  len_beta,
bool  reciprocal_beta,
bool  flip_sign_beta 
)
void viennacl::linalg::ambm_m ( matrix_base< NumericT, F > &  mat1,
matrix_base< NumericT, F > const &  mat2,
ScalarType1 const &  alpha,
vcl_size_t  len_alpha,
bool  reciprocal_alpha,
bool  flip_sign_alpha,
matrix_base< NumericT, F > const &  mat3,
ScalarType2 const &  beta,
vcl_size_t  len_beta,
bool  reciprocal_beta,
bool  flip_sign_beta 
)
void viennacl::linalg::amg_init ( MatrixType const &  mat,
InternalType1 &  A,
InternalType1 &  P,
InternalType2 &  Pointvector,
amg_tag &  tag 
)

Initialize AMG preconditioner.

Parameters:
matSystem matrix
AOperator matrices on all levels
PProlongation/Interpolation operators on all levels
PointvectorVector of points on all levels
tagAMG preconditioner tag
void viennacl::linalg::amg_lu ( boost::numeric::ublas::compressed_matrix< ScalarType > &  op,
boost::numeric::ublas::permutation_matrix<> &  Permutation,
SparseMatrixType const &  A 
)

Pre-compute LU factorization for direct solve (ublas library).

Speeds up precondition phase as this is computed only once overall instead of once per iteration.

Parameters:
opOperator matrix for direct solve
PermutationPermutation matrix which saves the factorization result
AOperator matrix on coarsest level
void viennacl::linalg::amg_setup ( InternalType1 &  A,
InternalType1 &  P,
InternalType2 &  Pointvector,
amg_tag &  tag 
)

Setup AMG preconditioner.

Parameters:
AOperator matrices on all levels
PProlongation/Interpolation operators on all levels
PointvectorVector of points on all levels
tagAMG preconditioner tag
void viennacl::linalg::amg_setup_apply ( InternalVectorType &  result,
InternalVectorType &  rhs,
InternalVectorType &  residual,
SparseMatrixType const &  A,
amg_tag const &  tag 
)

Setup data structures for precondition phase.

Parameters:
resultResult vector on all levels
rhsRHS vector on all levels
residualResidual vector on all levels
AOperators matrices on all levels from setup phase
tagAMG preconditioner tag
void viennacl::linalg::amg_setup_apply ( InternalVectorType &  result,
InternalVectorType &  rhs,
InternalVectorType &  residual,
SparseMatrixType const &  A,
amg_tag const &  tag,
viennacl::context  ctx 
)

Setup data structures for precondition phase for later use on the GPU.

Parameters:
resultResult vector on all levels
rhsRHS vector on all levels
residualResidual vector on all levels
AOperators matrices on all levels from setup phase
tagAMG preconditioner tag
ctxOptional context in which the auxiliary objects are created (one out of multiple OpenCL contexts, CUDA, host)
void viennacl::linalg::amg_transform_cpu ( InternalType1 &  A,
InternalType1 &  P,
InternalType1 &  R,
InternalType2 &  A_setup,
InternalType2 &  P_setup,
amg_tag &  tag 
)

Save operators after setup phase for CPU computation.

Parameters:
AOperator matrices on all levels on the CPU
PProlongation/Interpolation operators on all levels on the CPU
RRestriction operators on all levels on the CPU
A_setupOperators matrices on all levels from setup phase
P_setupProlongation/Interpolation operators on all levels from setup phase
tagAMG preconditioner tag
void viennacl::linalg::amg_transform_gpu ( InternalType1 &  A,
InternalType1 &  P,
InternalType1 &  R,
InternalType2 &  A_setup,
InternalType2 &  P_setup,
amg_tag &  tag,
viennacl::context  ctx 
)

Save operators after setup phase for GPU computation.

Parameters:
AOperator matrices on all levels on the GPU
PProlongation/Interpolation operators on all levels on the GPU
RRestriction operators on all levels on the GPU
A_setupOperators matrices on all levels from setup phase
P_setupProlongation/Interpolation operators on all levels from setup phase
tagAMG preconditioner tag
ctxOptional context in which the auxiliary objects are created (one out of multiple OpenCL contexts, CUDA, host)
viennacl::enable_if< viennacl::is_scalar<S1>::value && viennacl::is_scalar<S2>::value && viennacl::is_any_scalar<ScalarType1>::value >::type viennacl::linalg::as ( S1 &  s1,
S2 const &  s2,
ScalarType1 const &  alpha,
vcl_size_t  len_alpha,
bool  reciprocal_alpha,
bool  flip_sign_alpha 
)

Interface for the generic operation s1 = s2 @ alpha, where s1 and s2 are GPU scalars, @ denotes multiplication or division, and alpha is either a GPU or a CPU scalar.

Parameters:
s1The first (GPU) scalar
s2The second (GPU) scalar
alphaThe scalar alpha in the operation
len_alphaIf alpha is obtained from summing over a small GPU vector (e.g. the final summation after a multi-group reduction), then supply the length of the array here
reciprocal_alphaIf true, then s2 / alpha instead of s2 * alpha is computed
flip_sign_alphaIf true, then (-alpha) is used instead of alpha
viennacl::enable_if< viennacl::is_scalar<S1>::value && viennacl::is_scalar<S2>::value && viennacl::is_scalar<S3>::value && viennacl::is_any_scalar<ScalarType1>::value && viennacl::is_any_scalar<ScalarType2>::value >::type viennacl::linalg::asbs ( S1 &  s1,
S2 const &  s2,
ScalarType1 const &  alpha,
vcl_size_t  len_alpha,
bool  reciprocal_alpha,
bool  flip_sign_alpha,
S3 const &  s3,
ScalarType2 const &  beta,
vcl_size_t  len_beta,
bool  reciprocal_beta,
bool  flip_sign_beta 
)

Interface for the generic operation s1 = s2 @ alpha + s3 @ beta, where s1, s2 and s3 are GPU scalars, @ denotes multiplication or division, and alpha, beta are either a GPU or a CPU scalar.

Parameters:
s1The first (GPU) scalar
s2The second (GPU) scalar
alphaThe scalar alpha in the operation
len_alphaIf alpha is a small GPU vector, which needs to be summed in order to obtain the final scalar, then supply the length of the array here
reciprocal_alphaIf true, then s2 / alpha instead of s2 * alpha is computed
flip_sign_alphaIf true, then (-alpha) is used instead of alpha
s3The third (GPU) scalar
betaThe scalar beta in the operation
len_betaIf beta is obtained from summing over a small GPU vector (e.g. the final summation after a multi-group reduction), then supply the length of the array here
reciprocal_betaIf true, then s2 / beta instead of s2 * beta is computed
flip_sign_betaIf true, then (-beta) is used instead of beta
viennacl::enable_if< viennacl::is_scalar<S1>::value && viennacl::is_scalar<S2>::value && viennacl::is_scalar<S3>::value && viennacl::is_any_scalar<ScalarType1>::value && viennacl::is_any_scalar<ScalarType2>::value >::type viennacl::linalg::asbs_s ( S1 &  s1,
S2 const &  s2,
ScalarType1 const &  alpha,
vcl_size_t  len_alpha,
bool  reciprocal_alpha,
bool  flip_sign_alpha,
S3 const &  s3,
ScalarType2 const &  beta,
vcl_size_t  len_beta,
bool  reciprocal_beta,
bool  flip_sign_beta 
)

Interface for the generic operation s1 += s2 @ alpha + s3 @ beta, where s1, s2 and s3 are GPU scalars, @ denotes multiplication or division, and alpha, beta are either a GPU or a CPU scalar.

Parameters:
s1The first (GPU) scalar
s2The second (GPU) scalar
alphaThe scalar alpha in the operation
len_alphaIf alpha is a small GPU vector, which needs to be summed in order to obtain the final scalar, then supply the length of the array here
reciprocal_alphaIf true, then s2 / alpha instead of s2 * alpha is computed
flip_sign_alphaIf true, then (-alpha) is used instead of alpha
s3The third (GPU) scalar
betaThe scalar beta in the operation
len_betaIf beta is obtained from summing over a small GPU vector (e.g. the final summation after a multi-group reduction), then supply the length of the array here
reciprocal_betaIf true, then s2 / beta instead of s2 * beta is computed
flip_sign_betaIf true, then (-beta) is used instead of beta
void viennacl::linalg::av ( vector_base< T > &  vec1,
vector_base< T > const &  vec2,
ScalarType1 const &  alpha,
vcl_size_t  len_alpha,
bool  reciprocal_alpha,
bool  flip_sign_alpha 
)
void viennacl::linalg::avbv ( vector_base< T > &  vec1,
vector_base< T > const &  vec2,
ScalarType1 const &  alpha,
vcl_size_t  len_alpha,
bool  reciprocal_alpha,
bool  flip_sign_alpha,
vector_base< T > const &  vec3,
ScalarType2 const &  beta,
vcl_size_t  len_beta,
bool  reciprocal_beta,
bool  flip_sign_beta 
)
void viennacl::linalg::avbv_v ( vector_base< T > &  vec1,
vector_base< T > const &  vec2,
ScalarType1 const &  alpha,
vcl_size_t  len_alpha,
bool  reciprocal_alpha,
bool  flip_sign_alpha,
vector_base< T > const &  vec3,
ScalarType2 const &  beta,
vcl_size_t  len_beta,
bool  reciprocal_beta,
bool  flip_sign_beta 
)
std::vector< typename viennacl::result_of::cpu_value_type<typename VectorT::value_type>::type > viennacl::linalg::bisect ( VectorT const &  alphas,
VectorT const &  betas 
)

Implementation of the bisect-algorithm for the calculation of the eigenvalues of a tridiagonal matrix. Experimental - interface might change.

Parameters:
alphasElements of the main diagonal
betasElements of the secondary diagonal
Returns:
Returns the eigenvalues of the tridiagonal matrix defined by alpha and beta
void viennacl::linalg::convolve_i ( viennacl::vector< SCALARTYPE, ALIGNMENT > &  input1,
viennacl::vector< SCALARTYPE, ALIGNMENT > &  input2,
viennacl::vector< SCALARTYPE, ALIGNMENT > &  output 
)
viennacl::result_of::cpu_value_type<typename MatrixT::value_type>::type viennacl::linalg::eig ( MatrixT const &  matrix,
power_iter_tag const &  tag 
)

Implementation of the calculation of eigenvalues using poweriteration.

Parameters:
matrixThe system matrix
tagTag with termination factor
Returns:
Returns the largest eigenvalue computed by the power iteration method
std::vector< typename viennacl::result_of::cpu_value_type<typename MatrixT::value_type>::type > viennacl::linalg::eig ( MatrixT const &  matrix,
lanczos_tag const &  tag 
)

Implementation of the calculation of eigenvalues using lanczos.

Parameters:
matrixThe system matrix
tagTag with several options for the lanczos algorithm
Returns:
Returns the n largest eigenvalues (n defined in the lanczos_tag)
viennacl::vector_expression<const vector_base<T>, const vector_base<T>, op_element_binary<op_div> > viennacl::linalg::element_div ( vector_base< T > const &  v1,
vector_base< T > const &  v2 
)
void viennacl::linalg::element_op ( vector_base< T > &  vec1,
vector_expression< const vector_base< T >, const vector_base< T >, OP > const &  proxy 
)

Implementation of the element-wise operation v1 = v2 .* v3 and v1 = v2 ./ v3 (using MATLAB syntax)

Parameters:
vec1The result vector (or -range, or -slice)
proxyThe proxy object holding v2, v3 and the operation
void viennacl::linalg::element_op ( matrix_base< T, F > &  A,
matrix_expression< const matrix_base< T, F >, const matrix_base< T, F >, OP > const &  proxy 
)

Implementation of the element-wise operation A = B .* C and A = B ./ C for matrices (using MATLAB syntax). Don't use this function directly, use element_prod() and element_div().

Parameters:
AThe result matrix (or -range, or -slice)
proxyThe proxy object holding B, C, and the operation
viennacl::vector_expression<const vector_base<T>, const vector_base<T>, op_element_binary<op_prod> > viennacl::linalg::element_prod ( vector_base< T > const &  v1,
vector_base< T > const &  v2 
)
vcl_size_t index_norm_inf ( vector_base< T > const &  vec)

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

Parameters:
vecThe vector
Returns:
The result. Note that the result must be a CPU scalar
vcl_size_t index_norm_inf ( viennacl::vector_expression< LHS, RHS, OP > const &  vec)

Computes the supremum norm of a vector with final reduction on the CPU - interface for a vector expression. Creates a temporary.

Parameters:
vecThe vector expression
viennacl::enable_if< viennacl::is_stl< typename viennacl::traits::tag_of< VectorT1 >::type >::value, typename VectorT1::value_type>::type viennacl::linalg::inner_prod ( VectorT1 const &  v1,
VectorT2 const &  v2 
)
viennacl::scalar_expression< const vector_base<NumericT>, const vector_base<NumericT>, viennacl::op_inner_prod > viennacl::linalg::inner_prod ( vector_base< NumericT > const &  vector1,
vector_base< NumericT > const &  vector2 
)
viennacl::scalar_expression< const viennacl::vector_expression<LHS, RHS, OP>, const vector_base<NumericT>, viennacl::op_inner_prod > viennacl::linalg::inner_prod ( viennacl::vector_expression< LHS, RHS, OP > const &  vector1,
vector_base< NumericT > const &  vector2 
)
viennacl::scalar_expression< const vector_base<NumericT>, const viennacl::vector_expression<LHS, RHS, OP>, viennacl::op_inner_prod > viennacl::linalg::inner_prod ( vector_base< NumericT > const &  vector1,
viennacl::vector_expression< LHS, RHS, OP > const &  vector2 
)
viennacl::scalar_expression< const viennacl::vector_expression<LHS1, RHS1, OP1>, const viennacl::vector_expression<LHS2, RHS2, OP2>, viennacl::op_inner_prod > viennacl::linalg::inner_prod ( viennacl::vector_expression< LHS1, RHS1, OP1 > const &  vector1,
viennacl::vector_expression< LHS2, RHS2, OP2 > const &  vector2 
)
viennacl::vector_expression< const vector_base<NumericT>, const vector_tuple<NumericT>, viennacl::op_inner_prod > viennacl::linalg::inner_prod ( vector_base< NumericT > const &  x,
vector_tuple< NumericT > const &  y_tuple 
)
void inner_prod_cpu ( vector_base< T > const &  vec1,
vector_base< T > const &  vec2,
T &  result 
)

Computes the inner product of two vectors with the final reduction step on the CPU - dispatcher interface.

Parameters:
vec1The first vector
vec2The second vector
resultThe result scalar (on the gpu)
void inner_prod_cpu ( viennacl::vector_expression< LHS, RHS, OP > const &  vec1,
vector_base< T > const &  vec2,
T &  result 
)
void inner_prod_cpu ( vector_base< T > const &  vec1,
viennacl::vector_expression< LHS, RHS, OP > const &  vec2,
T &  result 
)
void inner_prod_cpu ( viennacl::vector_expression< LHS1, RHS1, OP1 > const &  vec1,
viennacl::vector_expression< LHS2, RHS2, OP2 > const &  vec2,
S3 &  result 
)
void viennacl::linalg::inner_prod_impl ( vector_base< T > const &  x,
vector_tuple< T > const &  y_tuple,
vector_base< T > &  result 
)

Computes the inner products <x, y1>, <x, y2>, ..., <x, y_N> and writes the result to a (sub-)vector.

Parameters:
xThe common vector
y_tupleA collection of vector, all of the same size.
resultThe result scalar (on the gpu). Needs to match the number of elements in y_tuple
void inner_prod_impl ( vector_base< T > const &  vec1,
vector_base< T > const &  vec2,
scalar< T > &  result 
)

Computes the inner product of two vectors - dispatcher interface.

Parameters:
vec1The first vector
vec2The second vector
resultThe result scalar (on the gpu)
void inner_prod_impl ( viennacl::vector_expression< LHS, RHS, OP > const &  vec1,
vector_base< T > const &  vec2,
scalar< T > &  result 
)
void inner_prod_impl ( vector_base< T > const &  vec1,
viennacl::vector_expression< LHS, RHS, OP > const &  vec2,
scalar< T > &  result 
)
void inner_prod_impl ( viennacl::vector_expression< LHS1, RHS1, OP1 > const &  vec1,
viennacl::vector_expression< LHS2, RHS2, OP2 > const &  vec2,
scalar< T > &  result 
)
std::vector<T> viennacl::linalg::inplace_qr ( viennacl::matrix< T, F, ALIGNMENT > &  A,
vcl_size_t  block_size = 16 
)

Overload of inplace-QR factorization of a ViennaCL matrix A.

Parameters:
AA dense ViennaCL matrix to be factored
block_sizeThe block size to be used.
std::vector<typename MatrixType::value_type> viennacl::linalg::inplace_qr ( MatrixType &  A,
vcl_size_t  block_size = 16 
)

Overload of inplace-QR factorization for a general Boost.uBLAS compatible matrix A.

Parameters:
AA dense compatible to Boost.uBLAS
block_sizeThe block size to be used.
void viennacl::linalg::inplace_qr_apply_trans_Q ( MatrixType const &  A,
VectorType1 const &  betas,
VectorType2 &  b 
)

Computes Q^T b, where Q is an implicit orthogonal matrix defined via its Householder reflectors stored in A.

Parameters:
AA matrix holding the Householder reflectors in the lower triangular part. Typically obtained from calling inplace_qr() on the original matrix
betasThe scalars beta_i for each Householder reflector (I - beta_i v_i v_i^T)
bThe vector b to which the result Q^T b is directly written to
void viennacl::linalg::inplace_qr_apply_trans_Q ( viennacl::matrix< T, F, ALIGNMENT > const &  A,
VectorType1 const &  betas,
viennacl::vector< T, A2 > &  b 
)
void viennacl::linalg::inplace_solve ( const matrix_base< NumericT, F1 > &  A,
matrix_base< NumericT, F2 > &  B,
SOLVERTAG   
)

Direct inplace solver for dense triangular systems. Matlab notation: A \ B.

Parameters:
AThe system matrix
BThe matrix of row vectors, where the solution is directly written to
void viennacl::linalg::inplace_solve ( const matrix_base< NumericT, F1 > &  A,
matrix_expression< const matrix_base< NumericT, F2 >, const matrix_base< NumericT, F2 >, op_trans >  proxy_B,
SOLVERTAG   
)

Direct inplace solver for dense triangular systems with transposed right hand side.

Parameters:
AThe system matrix
proxy_BThe transposed matrix of row vectors, where the solution is directly written to
void viennacl::linalg::inplace_solve ( const matrix_expression< const matrix_base< NumericT, F1 >, const matrix_base< NumericT, F1 >, op_trans > &  proxy_A,
matrix_base< NumericT, F2 > &  B,
SOLVERTAG   
)

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

Parameters:
proxy_AThe system matrix proxy
BThe matrix holding the load vectors, where the solution is directly written to
void viennacl::linalg::inplace_solve ( const matrix_expression< const matrix_base< NumericT, F1 >, const matrix_base< NumericT, F1 >, op_trans > &  proxy_A,
matrix_expression< const matrix_base< NumericT, F2 >, const matrix_base< NumericT, F2 >, op_trans >  proxy_B,
SOLVERTAG   
)

Direct inplace solver for dense transposed triangular systems with transposed right hand side. Matlab notation: A' \ B'.

Parameters:
proxy_AThe system matrix proxy
proxy_BThe matrix holding the load vectors, where the solution is directly written to
void viennacl::linalg::inplace_solve ( const matrix_base< NumericT, F > &  mat,
vector_base< NumericT > &  vec,
SOLVERTAG   
)
viennacl::enable_if< viennacl::is_any_sparse_matrix<SparseMatrixType>::value>::type viennacl::linalg::inplace_solve ( const SparseMatrixType &  mat,
viennacl::vector_base< ScalarType > &  vec,
SOLVERTAG  tag 
)

Carries out triangular inplace solves.

Parameters:
matThe matrix
vecThe vector
tagThe solver tag (lower_tag, unit_lower_tag, unit_upper_tag, or upper_tag)
void viennacl::linalg::inplace_solve ( const matrix_expression< const matrix_base< NumericT, F >, const matrix_base< NumericT, F >, op_trans > &  proxy,
vector_base< NumericT > &  vec,
SOLVERTAG   
)

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

Parameters:
proxyThe system matrix proxy
vecThe load vector, where the solution is directly written to
viennacl::enable_if< viennacl::is_any_sparse_matrix<SparseMatrixType>::value>::type viennacl::linalg::inplace_solve ( const matrix_expression< const SparseMatrixType, const SparseMatrixType, op_trans > &  mat,
viennacl::vector_base< ScalarType > &  vec,
SOLVERTAG  tag 
)

Carries out transposed triangular inplace solves.

Parameters:
matThe matrix
vecThe vector
tagThe solver tag (lower_tag, unit_lower_tag, unit_upper_tag, or upper_tag)
void viennacl::linalg::inplace_tred2 ( boost::numeric::ublas::matrix< ScalarType > const &  A,
vcl_size_t  block_size = 1 
)

Inplace reduction of a hermitian (or real symmetric) to tridiagonal form using householder similarity transforms (preserving eigenvalues)

Parameters:
AA dense matrix to be tridiagonalized
block_sizeThe block size to be used
void viennacl::linalg::lu_factorize ( matrix< SCALARTYPE, viennacl::row_major > &  A)

LU factorization of a row-major dense matrix.

Parameters:
AThe system matrix, where the LU matrices are directly written to. The implicit unit diagonal of L is not written.
void viennacl::linalg::lu_factorize ( matrix< SCALARTYPE, viennacl::column_major > &  A)

LU factorization of a column-major dense matrix.

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

LU substitution for the system LU = rhs.

Parameters:
AThe system matrix, where the LU matrices are directly written to. The implicit unit diagonal of L is not written.
vecThe load vector, where the solution is directly written to
void viennacl::linalg::matrix_assign ( matrix_base< NumericT, F > &  mat,
NumericT  s,
bool  clear = false 
)
void viennacl::linalg::matrix_column ( const matrix_base< NumericT, F > &  A,
unsigned int  j,
vector_base< NumericT > &  v 
)
void viennacl::linalg::matrix_diag_from_vector ( const vector_base< NumericT > &  v,
int  k,
matrix_base< NumericT, F > &  A 
)

Dispatcher interface for A = diag(v, k)

void viennacl::linalg::matrix_diag_to_vector ( const matrix_base< NumericT, F > &  A,
int  k,
vector_base< NumericT > &  v 
)

Dispatcher interface for v = diag(A, k)

void viennacl::linalg::matrix_diagonal_assign ( matrix_base< NumericT, F > &  mat,
NumericT  s 
)
void viennacl::linalg::matrix_row ( const matrix_base< NumericT, F > &  A,
unsigned int  i,
vector_base< NumericT > &  v 
)
void viennacl::linalg::nmf ( viennacl::matrix< ScalarType > const &  V,
viennacl::matrix< ScalarType > &  W,
viennacl::matrix< ScalarType > &  H,
nmf_config const &  conf 
)

The nonnegative matrix factorization (approximation) algorithm as suggested by Lee and Seung. Factorizes a matrix V with nonnegative entries into matrices W and H such that ||V - W*H|| is minimized.

Parameters:
VInput matrix
WFirst factor
HSecond factor
confA configuration object holding tolerances and the like
T viennacl::linalg::norm_1 ( std::vector< T, A > const &  v1)
viennacl::scalar_expression<const viennacl::vector_expression<const LHS, const RHS, OP>, const viennacl::vector_expression<const LHS, const RHS, OP>, viennacl::op_norm_1> viennacl::linalg::norm_1 ( viennacl::vector_expression< const LHS, const RHS, OP > const &  vector)
void norm_1_cpu ( vector_base< T > const &  vec,
T &  result 
)

Computes the l^1-norm of a vector with final reduction on the CPU.

Parameters:
vecThe vector
resultThe result scalar
void norm_1_cpu ( viennacl::vector_expression< LHS, RHS, OP > const &  vec,
S2 &  result 
)

Computes the l^1-norm of a vector with final reduction on the CPU - interface for a vector expression. Creates a temporary.

Parameters:
vecThe vector expression
resultThe result scalar
void viennacl::linalg::norm_1_impl ( viennacl::vector_expression< LHS, RHS, OP > const &  vec,
S2 &  result 
)

Computes the l^1-norm of a vector - interface for a vector expression. Creates a temporary.

Parameters:
vecThe vector expression
resultThe result scalar
void norm_1_impl ( vector_base< T > const &  vec,
scalar< T > &  result 
)

Computes the l^1-norm of a vector - dispatcher interface.

Parameters:
vecThe vector
resultThe result scalar
void viennacl::linalg::norm_1_impl ( viennacl::vector_expression< LHS, RHS, OP > const &  vec,
scalar< T > &  result 
)
T viennacl::linalg::norm_2 ( std::vector< T, A > const &  v1)
viennacl::scalar_expression<const viennacl::vector_expression<const LHS, const RHS, OP>, const viennacl::vector_expression<const LHS, const RHS, OP>, viennacl::op_norm_2> viennacl::linalg::norm_2 ( viennacl::vector_expression< const LHS, const RHS, OP > const &  vector)
void norm_2_cpu ( vector_base< T > const &  vec,
T &  result 
)

Computes the l^2-norm of a vector with final reduction on the CPU - dispatcher interface.

Parameters:
vecThe vector
resultThe result scalar
void norm_2_cpu ( viennacl::vector_expression< LHS, RHS, OP > const &  vec,
S2 &  result 
)

Computes the l^2-norm of a vector with final reduction on the CPU - interface for a vector expression. Creates a temporary.

Parameters:
vecThe vector expression
resultThe result scalar
void norm_2_impl ( vector_base< T > const &  vec,
scalar< T > &  result 
)

Computes the l^2-norm of a vector - dispatcher interface.

Parameters:
vecThe vector
resultThe result scalar
void norm_2_impl ( viennacl::vector_expression< LHS, RHS, OP > const &  vec,
scalar< T > &  result 
)

Computes the l^2-norm of a vector - interface for a vector expression. Creates a temporary.

Parameters:
vecThe vector expression
resultThe result scalar
scalar_expression< const matrix_base<NumericT, F>, const matrix_base<NumericT, F>, op_norm_frobenius> viennacl::linalg::norm_frobenius ( const matrix< NumericT, F > &  A)
void norm_frobenius_cpu ( matrix_base< T, F > const &  A,
T &  result 
)

Computes the Frobenius norm of a vector with final reduction on the CPU.

Parameters:
AThe matrix
resultThe result scalar
void norm_frobenius_impl ( matrix_base< T, F > const &  A,
scalar< T > &  result 
)

Computes the Frobenius norm of a matrix - dispatcher interface.

Parameters:
AThe matrix
resultThe result scalar
T viennacl::linalg::norm_inf ( std::vector< T, A > const &  v1)
viennacl::scalar_expression<const viennacl::vector_expression<const LHS, const RHS, OP>, const viennacl::vector_expression<const LHS, const RHS, OP>, viennacl::op_norm_inf> viennacl::linalg::norm_inf ( viennacl::vector_expression< const LHS, const RHS, OP > const &  vector)
void norm_inf_cpu ( vector_base< T > const &  vec,
T &  result 
)

Computes the supremum-norm of a vector with final reduction on the CPU.

Parameters:
vecThe vector
resultThe result scalar
void norm_inf_cpu ( viennacl::vector_expression< LHS, RHS, OP > const &  vec,
S2 &  result 
)

Computes the supremum norm of a vector with final reduction on the CPU - interface for a vector expression. Creates a temporary.

Parameters:
vecThe vector expression
resultThe result scalar
void norm_inf_impl ( vector_base< T > const &  vec,
scalar< T > &  result 
)

Computes the supremum-norm of a vector.

Parameters:
vecThe vector
resultThe result scalar
void norm_inf_impl ( viennacl::vector_expression< LHS, RHS, OP > const &  vec,
scalar< T > &  result 
)

Computes the supremum norm of a vector - interface for a vector expression. Creates a temporary.

Parameters:
vecThe vector expression
resultThe result scalar
viennacl::matrix_expression<const vector_base<NumericT>, const vector_base<NumericT>, op_prod> viennacl::linalg::outer_prod ( const vector_base< NumericT > &  vec1,
const vector_base< NumericT > &  vec2 
)

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

Parameters:
vec1The first vector
vec2The second vector
void viennacl::linalg::plane_rotation ( vector_base< T > &  vec1,
vector_base< T > &  vec2,
alpha,
beta 
)

Computes a plane rotation of two vectors.

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

Parameters:
vec1The first vector
vec2The second vector
alphaThe first transformation coefficient (CPU scalar)
betaThe second transformation coefficient (CPU scalar)
void viennacl::linalg::precondition ( viennacl::compressed_matrix< ScalarType > &  A,
ichol0_tag const &   
)

Implementation of a ILU-preconditioner with static pattern. Optimized version for CSR matrices.

Refer to Chih-Jen Lin and Jorge J. Moré, Incomplete Cholesky Factorizations with Limited Memory, SIAM J. Sci. Comput., 21(1), 24–45 for one of many descriptions of incomplete Cholesky Factorizations

Parameters:
AThe input matrix in CSR format // param tag An ichol0_tag in order to dispatch among several other preconditioners.
void viennacl::linalg::precondition ( viennacl::compressed_matrix< ScalarType > &  A,
ilu0_tag const &   
)

Implementation of a ILU-preconditioner with static pattern. Optimized version for CSR matrices.

refer to the Algorithm in Saad's book (1996 edition)

Parameters:
AThe sparse matrix matrix. The result is directly written to A.
void viennacl::linalg::precondition ( SparseMatrixType const &  A,
std::vector< std::map< SizeType, ScalarType > > &  output,
ilut_tag const &  tag 
)

Implementation of a ILU-preconditioner with threshold. Optimized implementation for compressed_matrix.

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

Parameters:
AThe input matrix. Either a compressed_matrix or of type std::vector< std::map<T, U> >
outputThe output matrix. Type requirements: const_iterator1 for iteration along rows, const_iterator2 for iteration along columns and write access via operator()
tagAn ilut_tag in order to dispatch among several other preconditioners.
VectorT viennacl::linalg::prod ( std::vector< std::vector< T, A1 >, A2 > const &  matrix,
VectorT const &  vector 
)
VectorT viennacl::linalg::prod ( std::vector< std::map< KEY, DATA, COMPARE, AMAP >, AVEC > const &  matrix,
VectorT const &  vector 
)
viennacl::matrix_expression< const viennacl::matrix_base<NumericT, F1>, const viennacl::matrix_base<NumericT, F2>, viennacl::op_mat_mat_prod > viennacl::linalg::prod ( viennacl::matrix_base< NumericT, F1 > const &  A,
viennacl::matrix_base< NumericT, F2 > const &  B 
)
viennacl::matrix_expression< const viennacl::matrix_base<NumericT, F1>, const viennacl::matrix_expression<const viennacl::matrix_base<NumericT, F2>, const viennacl::matrix_base<NumericT, F2>, op_trans>, viennacl::op_mat_mat_prod > viennacl::linalg::prod ( viennacl::matrix_base< NumericT, F1 > const &  A,
viennacl::matrix_expression< const viennacl::matrix_base< NumericT, F2 >, const viennacl::matrix_base< NumericT, F2 >, op_trans > const &  B 
)
viennacl::matrix_expression< const viennacl::matrix_expression<const viennacl::matrix_base<NumericT, F1>, const viennacl::matrix_base<NumericT, F1>, op_trans>, const viennacl::matrix_base<NumericT, F2>, viennacl::op_mat_mat_prod > viennacl::linalg::prod ( viennacl::matrix_expression< const viennacl::matrix_base< NumericT, F1 >, const viennacl::matrix_base< NumericT, F1 >, op_trans > const &  A,
viennacl::matrix_base< NumericT, F2 > const &  B 
)
viennacl::matrix_expression< const viennacl::matrix_expression<const viennacl::matrix_base<NumericT, F1>, const viennacl::matrix_base<NumericT, F1>, op_trans>, const viennacl::matrix_expression<const viennacl::matrix_base<NumericT, F2>, const viennacl::matrix_base<NumericT, F2>, op_trans>, viennacl::op_mat_mat_prod > viennacl::linalg::prod ( viennacl::matrix_expression< const viennacl::matrix_base< NumericT, F1 >, const viennacl::matrix_base< NumericT, F1 >, op_trans > const &  A,
viennacl::matrix_expression< const viennacl::matrix_base< NumericT, F2 >, const viennacl::matrix_base< NumericT, F2 >, op_trans > const &  B 
)
viennacl::vector_expression< const viennacl::matrix_base<NumericT, F>, const viennacl::vector_base<NumericT>, viennacl::op_prod > viennacl::linalg::prod ( viennacl::matrix_base< NumericT, F > const &  matrix,
viennacl::vector_base< NumericT > const &  vector 
)
viennacl::vector_expression< const viennacl::matrix_expression<const viennacl::matrix_base<NumericT, F>, const viennacl::matrix_base<NumericT, F>, op_trans>, const viennacl::vector_base<NumericT>, viennacl::op_prod > viennacl::linalg::prod ( viennacl::matrix_expression< const viennacl::matrix_base< NumericT, F >, const viennacl::matrix_base< NumericT, F >, op_trans > const &  matrix,
viennacl::vector_base< NumericT > const &  vector 
)
viennacl::enable_if< viennacl::is_any_sparse_matrix<SparseMatrixType>::value, vector_expression<const SparseMatrixType, const vector_base<SCALARTYPE>, op_prod > >::type viennacl::linalg::prod ( const SparseMatrixType &  mat,
const vector_base< SCALARTYPE > &  vec 
)
viennacl::enable_if< viennacl::is_any_sparse_matrix<SparseMatrixType>::value, viennacl::matrix_expression<const SparseMatrixType, const matrix_base < SCALARTYPE, F1 >, op_prod > >::type viennacl::linalg::prod ( const SparseMatrixType &  sp_mat,
const viennacl::matrix_base< SCALARTYPE, F1 > &  d_mat 
)
viennacl::enable_if< viennacl::is_any_sparse_matrix<SparseMatrixType>::value, viennacl::matrix_expression< const SparseMatrixType, const viennacl::matrix_expression<const viennacl::matrix_base<SCALARTYPE, F1>, const viennacl::matrix_base<SCALARTYPE, F1>, op_trans>, viennacl::op_prod > >::type viennacl::linalg::prod ( const SparseMatrixType &  A,
viennacl::matrix_expression< const viennacl::matrix_base< SCALARTYPE, F1 >, const viennacl::matrix_base< SCALARTYPE, F1 >, op_trans > const &  B 
)
viennacl::enable_if< viennacl::is_any_dense_structured_matrix<StructuredMatrixType>::value, vector_expression<const StructuredMatrixType, const vector_base<SCALARTYPE>, op_prod > >::type viennacl::linalg::prod ( const StructuredMatrixType &  mat,
const vector_base< SCALARTYPE > &  vec 
)
void viennacl::linalg::prod_impl ( const viennacl::vandermonde_matrix< SCALARTYPE, ALIGNMENT > &  mat,
const viennacl::vector_base< SCALARTYPE > &  vec,
viennacl::vector_base< SCALARTYPE > &  result 
)

Carries out matrix-vector multiplication with a vandermonde_matrix.

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

Parameters:
matThe matrix
vecThe vector
resultThe result vector
void viennacl::linalg::prod_impl ( const viennacl::hankel_matrix< SCALARTYPE, ALIGNMENT > &  mat,
const viennacl::vector_base< SCALARTYPE > &  vec,
viennacl::vector_base< SCALARTYPE > &  result 
)

Carries out matrix-vector multiplication with a hankel_matrix.

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

Parameters:
matThe matrix
vecThe vector
resultThe result vector
void viennacl::linalg::prod_impl ( const viennacl::toeplitz_matrix< SCALARTYPE, ALIGNMENT > &  mat,
const viennacl::vector_base< SCALARTYPE > &  vec,
viennacl::vector_base< SCALARTYPE > &  result 
)

Carries out matrix-vector multiplication with a toeplitz_matrix.

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

Parameters:
matThe matrix
vecThe vector
resultThe result vector
void viennacl::linalg::prod_impl ( const viennacl::circulant_matrix< SCALARTYPE, ALIGNMENT > &  mat,
const viennacl::vector_base< SCALARTYPE > &  vec,
viennacl::vector_base< SCALARTYPE > &  result 
)

Carries out matrix-vector multiplication with a circulant_matrix.

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

Parameters:
matThe matrix
vecThe vector
resultThe result vector
viennacl::enable_if< viennacl::is_any_sparse_matrix<SparseMatrixType>::value>::type viennacl::linalg::prod_impl ( const SparseMatrixType &  mat,
const viennacl::vector_base< ScalarType > &  vec,
viennacl::vector_base< ScalarType > &  result 
)

Carries out matrix-vector multiplication involving a sparse matrix type.

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

Parameters:
matThe matrix
vecThe vector
resultThe result vector
viennacl::enable_if< viennacl::is_any_sparse_matrix<SparseMatrixType>::value>::type viennacl::linalg::prod_impl ( const SparseMatrixType &  sp_mat,
const viennacl::matrix_base< ScalarType, F1 > &  d_mat,
viennacl::matrix_base< ScalarType, F2 > &  result 
)

Carries out matrix-matrix multiplication first matrix being sparse.

Implementation of the convenience expression result = prod(sp_mat, d_mat);

Parameters:
sp_matThe sparse matrix
d_matThe dense matrix
resultThe result matrix (dense)
viennacl::enable_if< viennacl::is_any_sparse_matrix<SparseMatrixType>::value>::type viennacl::linalg::prod_impl ( const SparseMatrixType &  sp_mat,
const viennacl::matrix_expression< const viennacl::matrix_base< ScalarType, F1 >, const viennacl::matrix_base< ScalarType, F1 >, viennacl::op_trans > &  d_mat,
viennacl::matrix_base< ScalarType, F2 > &  result 
)

Carries out matrix-matrix multiplication first matrix being sparse, and the second transposed.

Implementation of the convenience expression result = prod(sp_mat, d_mat);

Parameters:
sp_matThe sparse matrix
d_matThe dense matrix (transposed)
resultThe result matrix (dense)
void viennacl::linalg::prod_impl ( const matrix_base< NumericT, F1 > &  A,
const matrix_base< NumericT, F2 > &  B,
matrix_base< NumericT, F3 > &  C,
ScalarType  alpha,
ScalarType  beta 
)

Carries out matrix-matrix multiplication.

Implementation of C = prod(A, B);

void viennacl::linalg::prod_impl ( const viennacl::matrix_expression< const matrix_base< NumericT, F1 >, const matrix_base< NumericT, F1 >, op_trans > &  A,
const matrix_base< NumericT, F2 > &  B,
matrix_base< NumericT, F3 > &  C,
ScalarType  alpha,
ScalarType  beta 
)

Carries out matrix-matrix multiplication.

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

void viennacl::linalg::prod_impl ( const matrix_base< NumericT, F1 > &  A,
const viennacl::matrix_expression< const matrix_base< NumericT, F2 >, const matrix_base< NumericT, F2 >, op_trans > &  B,
matrix_base< NumericT, F3 > &  C,
ScalarType  alpha,
ScalarType  beta 
)

Carries out matrix-matrix multiplication.

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

void viennacl::linalg::prod_impl ( const viennacl::matrix_expression< const matrix_base< NumericT, F1 >, const matrix_base< NumericT, F1 >, op_trans > &  A,
const viennacl::matrix_expression< const matrix_base< NumericT, F2 >, const matrix_base< NumericT, F2 >, op_trans > &  B,
matrix_base< NumericT, F3 > &  C,
ScalarType  alpha,
ScalarType  beta 
)

Carries out matrix-matrix multiplication.

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

void prod_impl ( const matrix_base< NumericT, F > &  mat,
const vector_base< NumericT > &  vec,
vector_base< NumericT > &  result 
)

Carries out matrix-vector multiplication.

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

Parameters:
matThe matrix
vecThe vector
resultThe result vector
void prod_impl ( const matrix_expression< const matrix_base< NumericT, F >, const matrix_base< NumericT, F >, op_trans > &  mat_trans,
const vector_base< NumericT > &  vec,
vector_base< NumericT > &  result 
)

Carries out matrix-vector multiplication with a transposed matrix.

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

Parameters:
mat_transThe transposed matrix proxy
vecThe vector
resultThe result vector
viennacl::enable_if< viennacl::is_any_sparse_matrix<SparseMatrixType>::value, vector_expression<const SparseMatrixType, const vector<SCALARTYPE, ALIGNMENT>, op_prod > >::type viennacl::linalg::prod_impl ( const SparseMatrixType &  mat,
const vector< SCALARTYPE, ALIGNMENT > &  vec 
)
void viennacl::linalg::qr_method_nsm ( viennacl::matrix< SCALARTYPE, F, ALIGNMENT > &  A,
viennacl::matrix< SCALARTYPE, F, ALIGNMENT > &  Q,
boost::numeric::ublas::vector< SCALARTYPE > &  D,
boost::numeric::ublas::vector< SCALARTYPE > &  E 
)
void viennacl::linalg::qr_method_sym ( viennacl::matrix< SCALARTYPE, F, ALIGNMENT > &  A,
viennacl::matrix< SCALARTYPE, F, ALIGNMENT > &  Q,
boost::numeric::ublas::vector< SCALARTYPE > &  D 
)
void viennacl::linalg::recoverQ ( MatrixType const &  A,
VectorType const &  betas,
MatrixType &  Q,
MatrixType &  R 
)
void viennacl::linalg::scaled_rank_1_update ( matrix_base< NumericT, F > &  mat1,
S1 const &  alpha,
vcl_size_t  len_alpha,
bool  reciprocal_alpha,
bool  flip_sign_alpha,
const vector_base< NumericT > &  vec1,
const vector_base< NumericT > &  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:
mat1The matrix to be updated
alphaThe scaling factor (either a viennacl::scalar<>, float, or double)
len_alphaLength of the buffer for an eventual final reduction step (currently always '1')
reciprocal_alphaUse 1/alpha instead of alpha
flip_sign_alphaUse -alpha instead of alpha
vec1The first vector
vec2The second vector
ScalarType viennacl::linalg::setup_w ( viennacl::compressed_matrix< ScalarType > const &  A,
SizeType  row,
SparseVector &  w 
)

Dispatcher overload for extracting the row of nonzeros of a compressed matrix.

ScalarType viennacl::linalg::setup_w ( std::vector< std::map< SizeType, ScalarType > > const &  A,
SizeType  row,
SparseVector &  w 
)

Dispatcher overload for extracting the row of nonzeros of a STL-grown sparse matrix.

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:
matrixThe system matrix
rhsThe load vector
tagSolver configuration tag
Returns:
The result vector
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:
matrixThe system matrix
rhsThe load vector
tagSolver configuration tag
Returns:
The result vector
VectorType viennacl::linalg::solve ( const MatrixType &  matrix,
VectorType const &  rhs,
mixed_precision_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:
matrixThe system matrix
rhsThe load vector
tagSolver configuration tag
Returns:
The result vector
VectorType viennacl::linalg::solve ( const MatrixType &  matrix,
VectorType const &  rhs,
cg_tag const &  tag,
viennacl::linalg::no_precond   
)
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:
matrixThe system matrix
rhsThe load vector
tagSolver configuration tag
precondA preconditioner. Precondition operation is done via member function apply()
Returns:
The result vector
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:
matrixThe system matrix
rhsThe load vector
tagSolver configuration tag
precondA preconditioner. Precondition operation is done via member function apply()
Returns:
The result vector
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,
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:
matrixThe system matrix
rhsThe load vector
tagSolver configuration tag
precondA preconditioner. Precondition operation is done via member function apply()
Returns:
The result vector
VectorType viennacl::linalg::solve ( const MatrixType &  matrix,
VectorType const &  rhs,
mixed_precision_cg_tag const &  tag,
viennacl::linalg::no_precond   
)
matrix<NumericT, F2> viennacl::linalg::solve ( const matrix_base< NumericT, F1 > &  A,
const matrix_base< NumericT, F2 > &  B,
SOLVERTAG  tag 
)

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

Parameters:
AThe system matrix
BThe matrix of load vectors
tagDispatch tag
matrix<NumericT, F2> viennacl::linalg::solve ( const matrix_base< NumericT, F1 > &  A,
const matrix_expression< const matrix_base< NumericT, F2 >, const matrix_base< NumericT, F2 >, op_trans > &  proxy,
SOLVERTAG  tag 
)

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

Parameters:
AThe system matrix
proxyThe transposed load vector
tagDispatch tag
vector<NumericT> viennacl::linalg::solve ( const matrix_base< NumericT, F1 > &  mat,
const vector_base< NumericT > &  vec,
SOLVERTAG const &  tag 
)

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

Parameters:
matThe system matrix
vecThe load vector
tagDispatch tag
matrix<NumericT, F2> viennacl::linalg::solve ( const matrix_expression< const matrix_base< NumericT, F1 >, const matrix_base< NumericT, F1 >, op_trans > &  proxy,
const matrix_base< NumericT, F2 > &  B,
SOLVERTAG  tag 
)

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

Parameters:
proxyThe transposed system matrix proxy
BThe matrix of load vectors
tagDispatch 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.

matrix<NumericT, F2> viennacl::linalg::solve ( const matrix_expression< const matrix_base< NumericT, F1 >, const matrix_base< NumericT, F1 >, op_trans > &  proxy_A,
const matrix_expression< const matrix_base< NumericT, F2 >, const matrix_base< NumericT, F2 >, op_trans > &  proxy_B,
SOLVERTAG  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_AThe transposed system matrix proxy
proxy_BThe transposed matrix of load vectors, where the solution is directly written to
tagDispatch tag
vector<NumericT> viennacl::linalg::solve ( const matrix_expression< const matrix_base< NumericT, F1 >, const matrix_base< NumericT, F1 >, op_trans > &  proxy,
const vector_base< NumericT > &  vec,
SOLVERTAG 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:
proxyThe transposed system matrix proxy
vecThe load vector, where the solution is directly written to
tagDispatch tag
void viennacl::linalg::svd ( viennacl::matrix< SCALARTYPE, row_major, ALIGNMENT > &  A,
viennacl::matrix< SCALARTYPE, row_major, ALIGNMENT > &  QL,
viennacl::matrix< SCALARTYPE, row_major, ALIGNMENT > &  QR 
)

Computes the singular value decomposition of a matrix A. Experimental in 1.3.x.

Parameters:
AThe input matrix. Will be overwritten with a diagonal matrix containing the singular values on return
QLThe left orthogonal matrix
QRThe right orthogonal matrix
viennacl::enable_if< viennacl::is_scalar<S1>::value && viennacl::is_scalar<S2>::value >::type viennacl::linalg::swap ( S1 &  s1,
S2 &  s2 
)

Swaps the contents of two scalars.

Parameters:
s1The first scalar
s2The second scalar
void viennacl::linalg::vector_assign ( vector_base< T > &  vec1,
const T &  alpha,
bool  up_to_internal_size = false 
)

Assign a constant value to a vector (-range/-slice)

Parameters:
vec1The vector to which the value should be assigned
alphaThe value to be assigned
up_to_internal_sizeWhether 'alpha' should be written to padded memory as well. This is used for setting all entries to zero, including padded memory.
void viennacl::linalg::vector_swap ( vector_base< T > &  vec1,
vector_base< T > &  vec2 
)

Swaps the contents of two vectors, data is copied.

Parameters:
vec1The first vector (or -range, or -slice)
vec2The second vector (or -range, or -slice)

Variable Documentation

const char* double_float_conversion_program = "};\n"
const std::string SVD_BIDIAG_PACK_KERNEL = "bidiag_pack"
const std::string SVD_COPY_COL_KERNEL = "copy_col"
const std::string SVD_COPY_ROW_KERNEL = "copy_row"
const std::string SVD_FINAL_ITER_UPDATE_KERNEL = "final_iter_update"
const std::string SVD_GIVENS_NEXT_KERNEL = "givens_next"
const std::string SVD_GIVENS_PREV_KERNEL = "givens_prev"
const std::string SVD_HOUSEHOLDER_UPDATE_A_LEFT_KERNEL = "house_update_A_left"
const std::string SVD_HOUSEHOLDER_UPDATE_A_RIGHT_KERNEL = "house_update_A_right"
const std::string SVD_HOUSEHOLDER_UPDATE_QL_KERNEL = "house_update_QL"
const std::string SVD_HOUSEHOLDER_UPDATE_QR_KERNEL = "house_update_QR"
const std::string SVD_INVERSE_SIGNS_KERNEL = "inverse_signs"
const std::string SVD_MATRIX_TRANSPOSE_KERNEL = "transpose_inplace"
const std::string SVD_UPDATE_QR_COLUMN_KERNEL = "update_qr_column"