1 #ifndef VIENNACL_LINALG_DETAIL_SPAI_SPAI_HPP
2 #define VIENNACL_LINALG_DETAIL_SPAI_SPAI_HPP
44 #include "boost/numeric/ublas/vector.hpp"
45 #include "boost/numeric/ublas/matrix.hpp"
46 #include "boost/numeric/ublas/matrix_proxy.hpp"
47 #include "boost/numeric/ublas/vector_proxy.hpp"
48 #include "boost/numeric/ublas/storage.hpp"
49 #include "boost/numeric/ublas/io.hpp"
50 #include "boost/numeric/ublas/lu.hpp"
51 #include "boost/numeric/ublas/triangular.hpp"
52 #include "boost/numeric/ublas/matrix_expression.hpp"
68 #define VIENNACL_SPAI_K_b 20
80 template<
typename SparseVectorType>
82 for(
typename SparseVectorType::const_iterator vec_it = v.begin(); vec_it!= v.end(); ++vec_it){
83 std::cout<<
"[ "<<vec_it->first<<
" ]:"<<vec_it->second<<std::endl;
86 template<
typename DenseMatrixType>
88 for(
int i = 0; i < m.size2(); ++i){
89 for(
int j = 0; j < m.size1(); ++j){
90 std::cout<<m(j, i)<<
" ";
101 template<
typename SparseVectorType,
typename ScalarType>
103 for(
typename SparseVectorType::const_iterator v_it = v.begin(); v_it != v.end(); ++v_it){
104 res_v[v_it->first] += b*v_it->second;
114 template<
typename SparseVectorType,
typename ScalarType>
116 const unsigned int ind, SparseVectorType& res){
117 for(
typename SparseVectorType::const_iterator v_it = v.begin(); v_it != v.end(); ++v_it){
120 res[ind] -=
static_cast<ScalarType
>(1);
129 template<
typename SparseVectorType>
130 void build_index_set(
const std::vector<SparseVectorType>& A_v_c,
const SparseVectorType& v, std::vector<unsigned int>& J,
131 std::vector<unsigned int>& I){
142 template<
typename SparseMatrixType,
typename DenseMatrixType>
143 void initProjectSubMatrix(
const SparseMatrixType& A_in,
const std::vector<unsigned int>& J, std::vector<unsigned int>& I,
144 DenseMatrixType& A_out)
146 A_out.resize(I.size(), J.size(),
false);
150 A_out(i,j) = A_in(I[i],J[j]);
165 template<
typename SparseMatrixType,
typename DenseMatrixType,
typename SparseVectorType,
typename VectorType>
167 const std::vector<SparseVectorType>& A_v_c,
168 const std::vector<SparseVectorType>& M_v,
169 std::vector<std::vector<unsigned int> >& g_I,
170 std::vector<std::vector<unsigned int> >& g_J,
171 std::vector<DenseMatrixType>& g_A_I_J,
172 std::vector<VectorType>& g_b_v){
173 #ifdef VIENNACL_WITH_OPENMP
174 #pragma omp parallel for
176 for (
long i = 0; i < static_cast<long>(M_v.size()); ++i){
191 template<
typename SparseVectorType>
193 const std::vector<SparseVectorType>& M_v,
194 std::vector<std::vector<unsigned int> >& g_J,
195 std::vector<std::vector<unsigned int> >& g_I)
197 #ifdef VIENNACL_WITH_OPENMP
198 #pragma omp parallel for
200 for (
long i = 0; i < static_cast<long>(M_v.size()); ++i){
216 template<
typename ScalarType,
unsigned int MAT_ALIGNMENT,
typename SparseVectorType>
218 const std::vector<SparseVectorType>& A_v_c,
219 const std::vector<SparseVectorType>& M_v,
220 std::vector<cl_uint> g_is_update,
221 std::vector<std::vector<unsigned int> >& g_I,
222 std::vector<std::vector<unsigned int> >& g_J,
230 block_assembly(A, g_J, g_I, g_A_I_J, g_is_update, is_empty_block);
231 block_qr<ScalarType>(g_I, g_J, g_A_I_J, g_bv, g_is_update, ctx);
245 template<
typename ScalarType,
typename SparseVectorType>
247 unsigned int start_m_ind,
248 const std::vector<unsigned int> & J,
249 SparseVectorType & m)
251 unsigned int cnt = 0;
253 m[J[i]] = m_in[start_m_ind + cnt++];
272 template<
typename SparseVectorType,
typename ScalarType>
274 std::vector<SparseVectorType> & M_v,
275 std::vector<std::vector<unsigned int> >& g_I,
276 std::vector<std::vector<unsigned int> > & g_J,
279 std::vector<SparseVectorType> & g_res,
280 std::vector<cl_uint> & g_is_update,
285 unsigned int y_sz, m_sz;
286 std::vector<cl_uint> y_inds(M_v.size() + 1,
static_cast<cl_uint
>(0));
287 std::vector<cl_uint> m_inds(M_v.size() + 1,
static_cast<cl_uint
>(0));
292 std::vector<ScalarType> y_v(y_sz, static_cast<ScalarType>(0));
298 y_v[y_inds[i] + j] =
static_cast<ScalarType
>(1.0);
303 std::vector<ScalarType> m_v(m_sz, static_cast<cl_uint>(0));
309 static_cast<unsigned int>(
sizeof(ScalarType)*y_v.size()),
312 static_cast<unsigned int>(
sizeof(ScalarType)*m_v.size()),
315 static_cast<unsigned int>(
sizeof(cl_uint)*(g_I.size() + 1)),
318 static_cast<unsigned int>(
sizeof(cl_uint)*(g_is_update.size())),
326 g_A_I_J_vcl.
handle1(), g_is_update_vcl,
328 static_cast<unsigned int>(M_v.size())));
332 sizeof(ScalarType)*(m_v.size()),
333 &(m_v[0]), 0, NULL, NULL);
337 for(
long i = 0; i < static_cast<long>(M_v.size()); ++i)
344 compute_spai_residual<SparseVectorType, ScalarType>(A_v_c, M_v[i],
static_cast<unsigned int>(i), g_res[i]);
345 ScalarType res_norm = 0;
368 template<
typename SparseVectorType,
typename DenseMatrixType,
typename VectorType>
370 std::vector<DenseMatrixType>& g_R,
371 std::vector<VectorType>& g_b_v,
372 std::vector<std::vector<unsigned int> >& g_I,
373 std::vector<std::vector<unsigned int> >& g_J,
374 std::vector<SparseVectorType>& g_res,
375 std::vector<bool>& g_is_update,
376 std::vector<SparseVectorType>& M_v,
378 typedef typename DenseMatrixType::value_type ScalarType;
380 #ifdef VIENNACL_WITH_OPENMP
381 #pragma omp parallel for
383 for (
long i = 0; i < static_cast<long>(M_v.size()); ++i){
385 VectorType y = boost::numeric::ublas::zero_vector<ScalarType>(g_I[i].size());
387 projectI<VectorType, ScalarType>(g_I[i], y,
static_cast<unsigned int>(tag.
getBegInd() + i));
389 VectorType m_new = boost::numeric::ublas::zero_vector<ScalarType>(g_R[i].size2());
393 compute_spai_residual<SparseVectorType, ScalarType>(A_v_c, M_v[i],
static_cast<unsigned int>(tag.
getBegInd() + i), g_res[i]);
394 ScalarType res_norm = 0;
405 template<
typename VectorType>
408 for(
unsigned int i = 0; i < parallel_is_update.size(); ++i){
409 if(parallel_is_update[i])
421 template<
typename SparseMatrixType,
typename SparseVectorType>
423 for(
typename SparseMatrixType::const_iterator1 row_it = M_in.begin1(); row_it!= M_in.end1(); ++row_it){
425 for(
typename SparseMatrixType::const_iterator2 col_it = row_it.begin(); col_it != row_it.end(); ++col_it){
426 M_v[
static_cast<unsigned int>(col_it.index2())][static_cast<unsigned int>(col_it.index1())] = *col_it;
433 template<
typename SparseMatrixType,
typename SparseVectorType>
435 for(
typename SparseMatrixType::const_iterator1 row_it = M_in.begin1(); row_it!= M_in.end1(); ++row_it){
436 for(
typename SparseMatrixType::const_iterator2 col_it = row_it.begin(); col_it != row_it.end(); ++col_it){
437 M_v[
static_cast<unsigned int>(col_it.index1())][static_cast<unsigned int>(col_it.index2())] = *col_it;
445 template <
typename SizeType>
446 void write_set_to_array(
const std::vector<std::vector<SizeType> >& ind_set, std::vector<cl_uint>& a){
449 for(
vcl_size_t i = 0; i < ind_set.size(); ++i){
450 for(
vcl_size_t j = 0; j < ind_set[i].size(); ++j){
451 a[cnt++] =
static_cast<cl_uint
>(ind_set[i][j]);
467 template<
typename ScalarType,
unsigned int MAT_ALIGNMENT>
469 const std::vector<std::vector<unsigned int> >& g_I,
471 std::vector<cl_uint>& g_is_update,
472 bool& is_empty_block)
475 unsigned int sz_I, sz_J, sz_blocks;
476 std::vector<cl_uint> matrix_dims(g_I.size()*2,
static_cast<cl_uint
>(0));
477 std::vector<cl_uint> i_ind(g_I.size() + 1,
static_cast<cl_uint
>(0));
478 std::vector<cl_uint> j_ind(g_I.size() + 1,
static_cast<cl_uint
>(0));
479 std::vector<cl_uint> blocks_ind(g_I.size() + 1,
static_cast<cl_uint
>(0));
486 std::vector<cl_uint> I_set(sz_I, static_cast<cl_uint>(0));
488 std::vector<cl_uint> J_set(sz_J, static_cast<cl_uint>(0));
494 if (I_set.size() > 0 && J_set.size() > 0)
499 std::vector<ScalarType> con_A_I_J(sz_blocks, static_cast<ScalarType>(0));
505 static_cast<unsigned int>(
sizeof(ScalarType)*(sz_blocks)),
511 static_cast<unsigned int>(
sizeof(cl_uint)*2*static_cast<cl_uint>(g_I.size())),
517 static_cast<unsigned int>(
sizeof(cl_uint)*(g_I.size() + 1)),
523 static_cast<unsigned int>(
sizeof(cl_uint)*sz_I),
529 static_cast<unsigned int>(
sizeof(cl_uint)*sz_J),
535 static_cast<unsigned int>(
sizeof(cl_uint)*(g_I.size() + 1)),
541 static_cast<unsigned int>(
sizeof(cl_uint)*(g_I.size() + 1)),
546 static_cast<unsigned int>(
sizeof(cl_uint)*g_is_update.size()),
558 static_cast<unsigned int>(g_I.size())));
559 is_empty_block =
false;
562 is_empty_block =
true;
572 template<
typename SparseMatrixType,
typename SparseVectorType>
578 for(
unsigned int i = 0; i < M_v.size(); ++i){
579 for(
typename SparseVectorType::const_iterator vec_it = M_v[i].begin(); vec_it!=M_v[i].end(); ++vec_it){
580 M(vec_it->first, i) = vec_it->second;
586 for(
unsigned int i = 0; i < M_v.size(); ++i){
587 for(
typename SparseVectorType::const_iterator vec_it = M_v[i].begin(); vec_it!=M_v[i].end(); ++vec_it){
588 M(i, vec_it->first) = vec_it->second;
598 template<
typename MatrixType>
600 typedef typename MatrixType::value_type ScalarType;
601 std::vector<std::map<vcl_size_t, ScalarType> > temp_A(A_in.size2());
602 A.resize(A_in.size2(), A_in.size1(),
false);
604 for (
typename MatrixType::const_iterator1 row_it = A_in.begin1();
605 row_it != A_in.end1();
608 for (
typename MatrixType::const_iterator2 col_it = row_it.begin();
609 col_it != row_it.end();
612 temp_A[col_it.index2()][col_it.index1()] = *col_it;
618 for (
typename std::map<vcl_size_t, ScalarType>::const_iterator it = temp_A[i].begin();
619 it != temp_A[i].end();
621 A(i, it->first) = it->second;
641 template <
typename MatrixType>
643 typedef typename MatrixType::value_type ScalarType;
644 typedef typename boost::numeric::ublas::vector<ScalarType> VectorType;
646 typedef typename boost::numeric::ublas::matrix<ScalarType> DenseMatrixType;
648 unsigned int cur_iter = 0;
651 std::vector<SparseVectorType> A_v_c(M.size2());
652 std::vector<SparseVectorType> M_v(M.size2());
658 go_on = (tag.
getEndInd() <
static_cast<long>(M.size2()));
662 std::vector<bool> g_is_update(l_sz,
true);
669 std::vector<SparseVectorType> l_M_v(l_sz);
680 std::vector<DenseMatrixType> g_A_I_J(l_sz);
682 std::vector<VectorType> g_b_v(l_sz);
684 std::vector<SparseVectorType> g_res(l_sz);
686 std::vector<std::vector<unsigned int> > g_I(l_sz);
688 std::vector<std::vector<unsigned int> > g_J(l_sz);
692 if(cur_iter == 0)
block_set_up(A, A_v_c, l_M_v, g_I, g_J, g_A_I_J, g_b_v);
693 else block_update(A, A_v_c, g_res, g_is_update, g_I, g_J, g_b_v, g_A_I_J, tag);
709 M.resize(M.size1(), M.size2(),
false);
722 template <
typename ScalarType,
unsigned int MAT_ALIGNMENT>
724 const boost::numeric::ublas::compressed_matrix<ScalarType>& cpu_A,
725 boost::numeric::ublas::compressed_matrix<ScalarType>& cpu_M,
731 unsigned int cur_iter = 0;
732 std::vector<cl_uint> g_is_update(cpu_M.size2(),
static_cast<cl_uint
>(1));
735 std::vector<SparseVectorType> A_v_c(cpu_M.size2());
736 std::vector<SparseVectorType> M_v(cpu_M.size2());
739 std::vector<SparseVectorType> g_res(cpu_M.size2());
740 std::vector<std::vector<unsigned int> > g_I(cpu_M.size2());
741 std::vector<std::vector<unsigned int> > g_J(cpu_M.size2());
752 block_set_up(A, A_v_c, M_v, g_is_update, g_I, g_J, g_A_I_J_vcl, g_bv_vcl);
754 block_update(A, A_v_c, g_is_update, g_res, g_J, g_I, g_A_I_J_vcl, g_bv_vcl, tag);
759 least_square_solve<SparseVectorType, ScalarType>(A_v_c, M_v, g_I, g_J, g_A_I_J_vcl, g_bv_vcl, g_res, g_is_update, tag,
viennacl::traits::context(A));
761 if(tag.getIsStatic())
break;
764 cpu_M.resize(cpu_M.size1(), cpu_M.size2(),
false);
767 M.
resize(static_cast<unsigned int>(cpu_M.size1()), static_cast<unsigned int>(cpu_M.size2()));
viennacl::ocl::kernel & get_kernel(std::string const &program_name, std::string const &kernel_name)
Convenience function for retrieving the kernel of a program directly from the context.
Definition: context.hpp:470
Main kernel class for generating OpenCL kernels for the sparse approximate inverse preconditioners...
Definition: spai.hpp:570
viennacl::ocl::handle< cl_mem > & handle1()
Return handle to start indices.
Definition: block_vector.hpp:60
std::size_t vcl_size_t
Definition: forwards.h:58
void least_square_solve(std::vector< SparseVectorType > &A_v_c, std::vector< SparseVectorType > &M_v, std::vector< std::vector< unsigned int > > &g_I, std::vector< std::vector< unsigned int > > &g_J, block_matrix &g_A_I_J_vcl, block_vector &g_bv_vcl, std::vector< SparseVectorType > &g_res, std::vector< cl_uint > &g_is_update, const spai_tag &tag, viennacl::context ctx)
Solution of Least square problem on GPU.
Definition: spai.hpp:273
Represents a contigious matrices on GPU.
Definition: block_matrix.hpp:49
bool is_all_update(VectorType ¶llel_is_update)
Definition: spai.hpp:406
void resize(vcl_size_t new_size1, vcl_size_t new_size2, bool preserve=true)
Resize the matrix.
Definition: compressed_matrix.hpp:618
void sparse_norm_2(const SparseVectorType &v, ScalarType &norm)
Computation of Euclidean norm for sparse vector.
Definition: spai-dynamic.hpp:130
void index_set_up(const std::vector< SparseVectorType > &A_v_c, const std::vector< SparseVectorType > &M_v, std::vector< std::vector< unsigned int > > &g_J, std::vector< std::vector< unsigned int > > &g_I)
Setting up index set of columns and rows for all columns.
Definition: spai.hpp:192
Implementations of dense matrix related operations including matrix-vector products.
viennacl::ocl::handle< cl_mem > & handle()
Returns a handle to the elements.
Definition: block_matrix.hpp:56
bool getIsStatic() const
Definition: spai_tag.hpp:94
Generic interface for matrix-vector and matrix-matrix products. See viennacl/linalg/vector_operations...
Represents an OpenCL kernel within ViennaCL.
Definition: kernel.hpp:59
Implementation of the dense matrix class.
void vectorize_column_matrix(const SparseMatrixType &M_in, std::vector< SparseVectorType > &M_v)
Solution of Least square problem on CPU.
Definition: spai.hpp:422
void backwardSolve(const MatrixType &R, const VectorType &y, VectorType &x)
Solution of linear:R*x=y system by backward substitution.
Definition: spai-static.hpp:98
void block_set_up(const SparseMatrixType &A, const std::vector< SparseVectorType > &A_v_c, const std::vector< SparseVectorType > &M_v, std::vector< std::vector< unsigned int > > &g_I, std::vector< std::vector< unsigned int > > &g_J, std::vector< DenseMatrixType > &g_A_I_J, std::vector< VectorType > &g_b_v)
Setting up blocks and QR factorizing them on CPU.
Definition: spai.hpp:166
OpenCL kernel file for sparse approximate inverse operations.
void block_update(const SparseMatrixType &A, const std::vector< SparseVectorType > &A_v_c, std::vector< SparseVectorType > &g_res, std::vector< bool > &g_is_update, std::vector< std::vector< unsigned int > > &g_I, std::vector< std::vector< unsigned int > > &g_J, std::vector< VectorType > &g_b_v, std::vector< DenseMatrixType > &g_A_I_J, spai_tag const &tag)
CPU-based dynamic update for SPAI preconditioner.
Definition: spai-dynamic.hpp:258
Manages an OpenCL context and provides the respective convenience functions for creating buffers...
Definition: context.hpp:51
viennacl::ocl::handle< cl_mem > & handle2()
Returns a handle to the start indices of matrix.
Definition: block_matrix.hpp:64
void projectRows(const std::vector< SparseVectorType > &A_v_c, const std::vector< unsigned int > &J, std::vector< unsigned int > &I)
Row projection for matrix A(:,J) -> A(I,J), building index set of non-zero rows.
Definition: spai-static.hpp:169
viennacl::ocl::handle< cl_mem > & handle()
Return handle to the elements.
Definition: block_vector.hpp:56
void compute_spai_residual(const std::vector< SparseVectorType > &A_v_c, const SparseVectorType &v, const unsigned int ind, SparseVectorType &res)
Computation of residual res = A*v - e.
Definition: spai.hpp:115
void insert_sparse_columns(const std::vector< SparseVectorType > &M_v, SparseMatrixType &M, bool is_right)
Insertion of vectorized matrix column into original sparse matrix.
Definition: spai.hpp:573
long getEndInd() const
Definition: spai_tag.hpp:103
viennacl::ocl::handle< cl_command_queue > const & handle() const
Definition: command_queue.hpp:81
A tag for SPAI Contains values for the algorithm. Must be passed to spai_precond constructor.
Definition: spai_tag.hpp:63
const handle_type & handle2() const
Returns the OpenCL handle to the column index array.
Definition: compressed_matrix.hpp:701
Generic interface for the computation of inner products. See viennacl/linalg/vector_operations.hpp for implementations.
Implementation of a helper sparse vector class for SPAI. Experimental.
void computeSPAI(const MatrixType &A, MatrixType &M, spai_tag &tag)
Construction of SPAI preconditioner on CPU.
Definition: spai.hpp:642
void print_matrix(DenseMatrixType &m)
Definition: spai.hpp:87
void setBegInd(long beg_ind)
Definition: spai_tag.hpp:130
Represents a generic 'context' similar to an OpenCL context, but is backend-agnostic and thus also su...
Definition: context.hpp:39
#define VIENNACL_ERR_CHECK(err)
Definition: error.hpp:655
Implementation of a bunch of (small) matrices on GPU. Experimental.
long getBegInd() const
Definition: spai_tag.hpp:100
viennacl::ocl::command_queue & get_queue()
Definition: context.hpp:249
void enqueue(KernelType &k, viennacl::ocl::command_queue const &queue)
Enqueues a kernel in the provided queue.
Definition: enqueue.hpp:48
viennacl::ocl::handle< cl_mem > create_memory(cl_mem_flags flags, unsigned int size, void *ptr=NULL) const
Creates a memory buffer within the context.
Definition: context.hpp:199
void initProjectSubMatrix(const SparseMatrixType &A_in, const std::vector< unsigned int > &J, std::vector< unsigned int > &I, DenseMatrixType &A_out)
Initializes a dense matrix from a sparse one.
Definition: spai.hpp:143
const OCL_TYPE & get() const
Definition: handle.hpp:189
void write_set_to_array(const std::vector< std::vector< SizeType > > &ind_set, std::vector< cl_uint > &a)
Definition: spai.hpp:446
Implementations of incomplete factorization preconditioners. Convenience header file.
Implementation of a static SPAI. Experimental.
Implementation of the compressed_matrix class.
Implementation of a bunch of vectors on GPU. Experimental.
Implementations of operations using sparse matrices.
void block_assembly(const viennacl::compressed_matrix< ScalarType, MAT_ALIGNMENT > &A, const std::vector< std::vector< unsigned int > > &g_J, const std::vector< std::vector< unsigned int > > &g_I, block_matrix &g_A_I_J_vcl, std::vector< cl_uint > &g_is_update, bool &is_empty_block)
Assembly of blocks on GPU by a gived set of row indices: g_I and column indices: g_J.
Definition: spai.hpp:468
void print_sparse_vector(const SparseVectorType &v)
Definition: spai.hpp:81
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
void get_size(const std::vector< std::vector< SizeType > > &inds, SizeType &size)
Computes size of particular container of index set.
Definition: qr.hpp:129
void setEndInd(long end_ind)
Definition: spai_tag.hpp:132
Represents sparse vector based on std::map<unsigned int, ScalarType>
Definition: sparse_vector.hpp:52
Provides a QR factorization using a block-based approach.
Implementation of the spai tag holding SPAI configuration parameters. Experimental.
void init_start_inds(const std::vector< std::vector< SizeType > > &inds, std::vector< cl_uint > &start_inds)
Initializes start indices of particular index set.
Definition: qr.hpp:141
void apply_q_trans_vec(const MatrixType &R, const VectorType &b_v, VectorType &y)
Recovery Q from matrix R and vector of betas b_v.
Definition: qr.hpp:330
void fanOutVector(const VectorType &m_in, const std::vector< unsigned int > &J, SparseVectorType &m)
Projects solution of LS problem onto original column m.
Definition: spai-static.hpp:86
#define VIENNACL_SPAI_K_b
Definition: spai.hpp:68
void compute_blocks_size(const std::vector< std::vector< unsigned int > > &g_I, const std::vector< std::vector< unsigned int > > &g_J, unsigned int &sz, std::vector< cl_uint > &blocks_ind, std::vector< cl_uint > &matrix_dims)
**************************************** BLOCK FUNCTIONS ************************************// ...
Definition: qr.hpp:112
Implementations of the OpenCL backend, where all contexts are stored in.
bool getIsRight() const
Definition: spai_tag.hpp:97
void build_index_set(const std::vector< SparseVectorType > &A_v_c, const SparseVectorType &v, std::vector< unsigned int > &J, std::vector< unsigned int > &I)
Setting up index set of columns and rows for certain column.
Definition: spai.hpp:130
void custom_fan_out(const std::vector< ScalarType > &m_in, unsigned int start_m_ind, const std::vector< unsigned int > &J, SparseVectorType &m)
Elicitation of sparse vector m for particular column from m_in - contigious vector for all columns...
Definition: spai.hpp:246
unsigned int getIterationLimit() const
Definition: spai_tag.hpp:91
viennacl::context context(T const &t)
Returns an ID for the currently active memory domain of an object.
Definition: context.hpp:41
void buildColumnIndexSet(const SparseVectorType &v, std::vector< unsigned int > &J)
Builds index set of projected columns for current column of preconditioner.
Definition: spai-static.hpp:132
Represents a contigious vector on GPU.
Definition: block_vector.hpp:49
void vectorize_row_matrix(const SparseMatrixType &M_in, std::vector< SparseVectorType > &M_v)
Definition: spai.hpp:434
void single_qr(MatrixType &R, VectorType &b_v)
Inplace QR factorization via Householder reflections c.f. Gene H. Golub, Charles F. Van Loan "Matrix Computations" 3rd edition p.224.
Definition: qr.hpp:257
size_type global_work_size(int index=0) const
Returns the global work size at the respective dimension.
Definition: kernel.hpp:759
double getResidualNormThreshold() const
Definition: spai_tag.hpp:85
static void init(viennacl::ocl::context &ctx)
Definition: spai.hpp:577
viennacl::ocl::context const & context() const
Definition: handle.hpp:191
viennacl::ocl::handle< cl_mem > & handle1()
Returns a handle to the matrix dimensions.
Definition: block_matrix.hpp:60
size_type local_work_size(int index=0) const
Returns the local work size at the respective dimension.
Definition: kernel.hpp:750
void sparse_transpose(const MatrixType &A_in, MatrixType &A)
Transposition of sparse matrix.
Definition: spai.hpp:599
Implementation of the ViennaCL scalar class.
Implementation of a dynamic SPAI. Provides the routines for automatic pattern updates Experimental...
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 add_sparse_vectors(const SparseVectorType &v, const ScalarType b, SparseVectorType &res_v)
Add two sparse vectors res_v = b*v.
Definition: spai.hpp:102