1 #ifndef VIENNACL_LINALG_ICHOL_HPP_ 
    2 #define VIENNACL_LINALG_ICHOL_HPP_ 
   54 template<
typename NumericT>
 
   59   NumericT           * elements   = viennacl::linalg::host_based::detail::extract_raw_pointer<NumericT>(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       NumericT 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         NumericT 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 MatrixT>
 
  129   typedef typename MatrixT::value_type      NumericType;
 
  140   template<
typename VectorT>
 
  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     NumericType  
const * elements   = viennacl::linalg::host_based::detail::extract_raw_pointer<NumericType>(LLT.
handle());
 
  148     viennacl::linalg::host_based::detail::csr_trans_inplace_solve<NumericType>(row_buffer, col_buffer, elements, vec, LLT.
size2(), 
lower_tag());
 
  149     viennacl::linalg::host_based::detail::csr_inplace_solve<NumericType>(row_buffer, col_buffer, elements, vec, LLT.
size2(), 
upper_tag());
 
  153   void init(MatrixT 
const & mat)
 
  162   ichol0_tag 
const & tag_;
 
  171 template<
typename NumericT, 
unsigned int AlignmentV>
 
  206   void init(MatrixType 
const & mat)
 
  215   ichol0_tag 
const & tag_;
 
const vcl_size_t & size2() const 
Returns the number of columns. 
void inplace_solve(const matrix_base< NumericT > &A, matrix_base< NumericT > &B, SolverTagT)
Direct inplace solver for triangular systems with multiple right hand sides, i.e. A \ B (MATLAB notat...
Incomplete Cholesky preconditioner class with static pattern (ICHOL0), can be supplied to solve()-rou...
A tag for incomplete Cholesky factorization with static pattern (ILU0) 
const vcl_size_t & size1() const 
Returns the number of rows. 
vcl_size_t size1(MatrixType const &mat)
Generic routine for obtaining the number of rows of a matrix (ViennaCL, uBLAS, etc.) 
A tag class representing a lower triangular matrix. 
void precondition(viennacl::compressed_matrix< NumericT > &A, ilu0_tag const &)
Implementation of a ILU-preconditioner with static pattern. Optimized version for CSR matrices...
This file provides the forward declarations for the main types used within ViennaCL. 
void apply(vector< NumericT > &vec) const 
result_of::size_type< MatrixType >::type size2(MatrixType const &mat)
Generic routine for obtaining the number of columns of a matrix (ViennaCL, uBLAS, etc...
const handle_type & handle() const 
Returns the OpenCL handle to the matrix entry array. 
const handle_type & handle1() const 
Returns the OpenCL handle to the row index array. 
Represents a generic 'context' similar to an OpenCL context, but is backend-agnostic and thus also su...
A tag class representing an upper triangular matrix. 
Implementation of the compressed_matrix class. 
const handle_type & handle2() const 
Returns the OpenCL handle to the column index array. 
ichol0_precond(MatrixT const &mat, ichol0_tag const &tag)
Common routines for single-threaded or OpenMP-enabled execution on CPU. 
void trans(const matrix_expression< const matrix_base< NumericT, SizeT, DistanceT >, const matrix_base< NumericT, SizeT, DistanceT >, op_trans > &proxy, matrix_base< NumericT > &temp_trans)
void apply(VectorT &vec) const 
viennacl::context context(T const &t)
Returns an ID for the currently active memory domain of an object. 
void copy(std::vector< NumericT > &cpu_vec, circulant_matrix< NumericT, AlignmentV > &gpu_mat)
Copies a circulant matrix from the std::vector to the OpenCL device (either GPU or multi-core CPU) ...
ichol0_precond(MatrixType const &mat, ichol0_tag const &tag)
void switch_memory_context(T &obj, viennacl::context new_ctx)
Generic convenience routine for migrating data of an object to a new memory domain.