1 #ifndef VIENNACL_LINALG_ICHOL_HPP_
2 #define VIENNACL_LINALG_ICHOL_HPP_
54 template<
typename ScalarType>
59 ScalarType * elements = viennacl::linalg::host_based::detail::extract_raw_pointer<ScalarType>(A.
handle());
60 unsigned int const * row_buffer = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(A.
handle1());
61 unsigned int const * col_buffer = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(A.
handle2());
66 unsigned int row_i_begin = row_buffer[i];
67 unsigned int row_i_end = row_buffer[i+1];
71 for (
unsigned int buf_index_aii = row_i_begin; buf_index_aii < row_i_end; ++buf_index_aii)
73 if (col_buffer[buf_index_aii] == i)
75 a_ii = std::sqrt(elements[buf_index_aii]);
76 elements[buf_index_aii] = a_ii;
82 for (
unsigned int buf_index_aii = row_i_begin; buf_index_aii < row_i_end; ++buf_index_aii)
84 if (col_buffer[buf_index_aii] > i)
85 elements[buf_index_aii] /= a_ii;
89 for (
unsigned int buf_index_j = row_i_begin; buf_index_j < row_i_end; ++buf_index_j)
91 unsigned int j = col_buffer[buf_index_j];
95 ScalarType a_ji = elements[buf_index_j];
97 for (
unsigned int buf_index_k = row_i_begin; buf_index_k < row_i_end; ++buf_index_k)
99 unsigned int k = col_buffer[buf_index_k];
103 ScalarType a_ki = elements[buf_index_k];
106 unsigned int row_j_begin = row_buffer[j];
107 unsigned int row_j_end = row_buffer[j+1];
108 for (
unsigned int buf_index_kj = row_j_begin; buf_index_kj < row_j_end; ++buf_index_kj)
110 if (col_buffer[buf_index_kj] == k)
112 elements[buf_index_kj] -= a_ki * a_ji;
126 template <
typename MatrixType>
129 typedef typename MatrixType::value_type ScalarType;
140 template <
typename VectorType>
143 unsigned int const * row_buffer = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(LLT.
handle1());
144 unsigned int const * col_buffer = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(LLT.
handle2());
145 ScalarType
const * elements = viennacl::linalg::host_based::detail::extract_raw_pointer<ScalarType>(LLT.
handle());
148 viennacl::linalg::host_based::detail::csr_trans_inplace_solve<ScalarType>(row_buffer, col_buffer, elements, vec, LLT.
size2(),
lower_tag());
149 viennacl::linalg::host_based::detail::csr_inplace_solve<ScalarType>(row_buffer, col_buffer, elements, vec, LLT.
size2(),
upper_tag());
153 void init(MatrixType
const & mat)
162 ichol0_tag
const & tag_;
171 template <
typename ScalarType,
unsigned int MAT_ALIGNMENT>
206 void init(MatrixType
const & mat)
215 ichol0_tag
const & tag_;
std::size_t vcl_size_t
Definition: forwards.h:58
Incomplete Cholesky preconditioner class with static pattern (ICHOL0), can be supplied to solve()-rou...
Definition: ichol.hpp:127
viennacl::enable_if< viennacl::is_any_sparse_matrix< M1 >::value, matrix_expression< const M1, const M1, op_trans > >::type trans(const M1 &mat)
Returns an expression template class representing a transposed matrix.
Definition: sparse_matrix_operations.hpp:330
A tag for incomplete Cholesky factorization with static pattern (ILU0)
Definition: ichol.hpp:43
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.
Definition: direct_solve.hpp:54
vcl_size_t size1(MatrixType const &mat)
Generic routine for obtaining the number of rows of a matrix (ViennaCL, uBLAS, etc.)
Definition: size.hpp:216
A tag class representing a lower triangular matrix.
Definition: forwards.h:703
This file provides the forward declarations for the main types used within ViennaCL.
const handle_type & handle2() const
Returns the OpenCL handle to the column index array.
Definition: compressed_matrix.hpp:701
result_of::size_type< MatrixType >::type size2(MatrixType const &mat)
Generic routine for obtaining the number of columns of a matrix (ViennaCL, uBLAS, etc...
Definition: size.hpp:245
void apply(VectorType &vec) const
Definition: ichol.hpp:141
const vcl_size_t & size2() const
Returns the number of columns.
Definition: compressed_matrix.hpp:694
Represents a generic 'context' similar to an OpenCL context, but is backend-agnostic and thus also su...
Definition: context.hpp:39
ichol0_precond(MatrixType const &mat, ichol0_tag const &tag)
Definition: ichol.hpp:132
A tag class representing an upper triangular matrix.
Definition: forwards.h:708
Implementation of the compressed_matrix class.
void copy(std::vector< SCALARTYPE > &cpu_vec, circulant_matrix< SCALARTYPE, ALIGNMENT > &gpu_mat)
Copies a circulant matrix from the std::vector to the OpenCL device (either GPU or multi-core CPU) ...
Definition: circulant_matrix.hpp:150
const vcl_size_t & size1() const
Returns the number of rows.
Definition: compressed_matrix.hpp:692
void precondition(viennacl::compressed_matrix< ScalarType > &A, ilu0_tag const &)
Implementation of a ILU-preconditioner with static pattern. Optimized version for CSR matrices...
Definition: ilu0.hpp:79
viennacl::context context(T const &t)
Returns an ID for the currently active memory domain of an object.
Definition: context.hpp:41
Common routines for single-threaded or OpenMP-enabled execution on CPU.
Definition: forwards.h:479
ichol0_precond(MatrixType const &mat, ichol0_tag const &tag)
Definition: ichol.hpp:177
void apply(vector< ScalarType > &vec) const
Definition: ichol.hpp:185
const handle_type & handle() const
Returns the OpenCL handle to the matrix entry array.
Definition: compressed_matrix.hpp:703
const handle_type & handle1() const
Returns the OpenCL handle to the row index array.
Definition: compressed_matrix.hpp:699
void switch_memory_context(T &obj, viennacl::context new_ctx)
Generic convenience routine for migrating data of an object to a new memory domain.
Definition: memory.hpp:624