1 #ifndef VIENNACL_COMPRESSED_MATRIX_HPP_ 
    2 #define VIENNACL_COMPRESSED_MATRIX_HPP_ 
   36 #ifdef VIENNACL_WITH_UBLAS 
   37 #include <boost/numeric/ublas/matrix_sparse.hpp> 
   49   template<
typename CPUMatrixT, 
typename NumericT, 
unsigned int AlignmentV>
 
   59     std::vector<NumericT> elements(nonzeros);
 
   64     for (
typename CPUMatrixT::const_iterator1 row_it = cpu_matrix.begin1();
 
   65          row_it != cpu_matrix.end1();
 
   68       row_buffer.set(row_index, data_index);
 
   71       for (
typename CPUMatrixT::const_iterator2 col_it = row_it.begin();
 
   72            col_it != row_it.end();
 
   75         col_buffer.set(data_index, col_it.index2());
 
   76         elements[data_index] = *col_it;
 
   79       data_index = viennacl::tools::align_to_multiple<vcl_size_t>(data_index, AlignmentV); 
 
   81     row_buffer.set(row_index, data_index);
 
   83     gpu_matrix.
set(row_buffer.get(),
 
  111 template<
typename CPUMatrixT, 
typename NumericT, 
unsigned int AlignmentV>
 
  112 void copy(
const CPUMatrixT & cpu_matrix,
 
  115   if ( cpu_matrix.size1() > 0 && cpu_matrix.size2() > 0 )
 
  119     for (
typename CPUMatrixT::const_iterator1 row_it = cpu_matrix.begin1();
 
  120          row_it != cpu_matrix.end1();
 
  124       for (
typename CPUMatrixT::const_iterator2 col_it = row_it.begin();
 
  125            col_it != row_it.end();
 
  130       num_entries += viennacl::tools::align_to_multiple<vcl_size_t>(entries_per_row, AlignmentV);
 
  133     if (num_entries == 0) 
 
  148 template<
typename SizeT, 
typename NumericT, 
unsigned int AlignmentV>
 
  149 void copy(
const std::vector< std::map<SizeT, NumericT> > & cpu_matrix,
 
  154   for (
vcl_size_t i=0; i<cpu_matrix.size(); ++i)
 
  156     if (cpu_matrix[i].
size() > 0)
 
  157       nonzeros += ((cpu_matrix[i].
size() - 1) / AlignmentV + 1) * AlignmentV;
 
  158     if (cpu_matrix[i].
size() > 0)
 
  159       max_col = std::max<vcl_size_t>(max_col, (cpu_matrix[i].rbegin())->first);
 
  167 #ifdef VIENNACL_WITH_UBLAS 
  172 template<
typename ScalarType, 
typename F, vcl_
size_t IB, 
typename IA, 
typename TA>
 
  173 void copy(
const boost::numeric::ublas::compressed_matrix<ScalarType, F, IB, IA, TA> & ublas_matrix,
 
  181   for (
vcl_size_t i=0; i<=ublas_matrix.size1(); ++i)
 
  182     row_buffer.
set(i, ublas_matrix.index1_data()[i]);
 
  185   for (
vcl_size_t i=0; i<ublas_matrix.nnz(); ++i)
 
  186     col_buffer.
set(i, ublas_matrix.index2_data()[i]);
 
  188   gpu_matrix.
set(row_buffer.get(),
 
  190                  &(ublas_matrix.value_data()[0]),
 
  191       ublas_matrix.size1(),
 
  192       ublas_matrix.size2(),
 
  198 #ifdef VIENNACL_WITH_ARMADILLO 
  204 template<
typename NumericT, 
unsigned int AlignmentV>
 
  205 void copy(arma::SpMat<NumericT> 
const & arma_matrix,
 
  208   assert( (vcl_matrix.
size1() == 0 || 
static_cast<vcl_size_t>(arma_matrix.n_rows) == vcl_matrix.
size1()) && 
bool(
"Size mismatch") );
 
  209   assert( (vcl_matrix.
size2() == 0 || 
static_cast<vcl_size_t>(arma_matrix.n_cols) == vcl_matrix.
size2()) && 
bool(
"Size mismatch") );
 
  216   for (
vcl_size_t col=0; col < static_cast<vcl_size_t>(arma_matrix.n_cols); ++col)
 
  220     for (
vcl_size_t i = col_begin; i < col_end; ++i)
 
  222       unsigned int row = arma_matrix.row_indices[i];
 
  223       row_buffer.set(row, row_buffer[row] + 1);
 
  228   unsigned int offset = 0;
 
  229   for (
vcl_size_t i=0; i<row_buffer.size(); ++i)
 
  231     unsigned int tmp = row_buffer[i];
 
  232     row_buffer.set(i, offset);
 
  237   std::vector<unsigned int> row_offsets(arma_matrix.n_rows);
 
  238   for (
vcl_size_t col=0; col < static_cast<vcl_size_t>(arma_matrix.n_cols); ++col)
 
  242     for (
vcl_size_t i = col_begin; i < col_end; ++i)
 
  244       unsigned int row = arma_matrix.row_indices[i];
 
  245       col_buffer.set(row_buffer[row] + row_offsets[row], col);
 
  246       value_buffer.set(row_buffer[row] + row_offsets[row], arma_matrix.values[i]);
 
  247       row_offsets[
row] += 1;
 
  251   vcl_matrix.
set(row_buffer.get(), col_buffer.get(), 
reinterpret_cast<NumericT*
>(value_buffer.get()),
 
  252                  arma_matrix.n_rows, arma_matrix.n_cols, arma_matrix.n_nonzero);
 
  256 #ifdef VIENNACL_WITH_EIGEN 
  261 template<
typename NumericT, 
int flags, 
unsigned int AlignmentV>
 
  262 void copy(
const Eigen::SparseMatrix<NumericT, flags> & eigen_matrix,
 
  263           compressed_matrix<NumericT, AlignmentV> & gpu_matrix)
 
  265   assert( (gpu_matrix.size1() == 0 || 
static_cast<vcl_size_t>(eigen_matrix.rows()) == gpu_matrix.size1()) && 
bool(
"Size mismatch") );
 
  266   assert( (gpu_matrix.size2() == 0 || 
static_cast<vcl_size_t>(eigen_matrix.cols()) == gpu_matrix.size2()) && 
bool(
"Size mismatch") );
 
  268   std::vector< std::map<unsigned int, NumericT> >  stl_matrix(eigen_matrix.rows());
 
  270   for (
int k=0; k < eigen_matrix.outerSize(); ++k)
 
  271     for (
typename Eigen::SparseMatrix<NumericT, flags>::InnerIterator it(eigen_matrix, k); it; ++it)
 
  272       stl_matrix[it.row()][it.col()] = it.value();
 
  274   copy(tools::const_sparse_matrix_adapter<NumericT>(stl_matrix, eigen_matrix.rows(), eigen_matrix.cols()), gpu_matrix);
 
  279 #ifdef VIENNACL_WITH_MTL4 
  284 template<
typename NumericT, 
unsigned int AlignmentV>
 
  285 void copy(
const mtl::compressed2D<NumericT> & cpu_matrix,
 
  286           compressed_matrix<NumericT, AlignmentV> & gpu_matrix)
 
  288   assert( (gpu_matrix.size1() == 0 || 
static_cast<vcl_size_t>(cpu_matrix.num_rows()) == gpu_matrix.size1()) && 
bool(
"Size mismatch") );
 
  289   assert( (gpu_matrix.size2() == 0 || 
static_cast<vcl_size_t>(cpu_matrix.num_cols()) == gpu_matrix.size2()) && 
bool(
"Size mismatch") );
 
  291   typedef mtl::compressed2D<NumericT>  MatrixType;
 
  293   std::vector< std::map<unsigned int, NumericT> >  stl_matrix(cpu_matrix.num_rows());
 
  295   using mtl::traits::range_generator;
 
  299   typedef typename min<range_generator<mtl::tag::row, MatrixType>,
 
  300       range_generator<mtl::tag::col, MatrixType> >::type   range_type;
 
  304   typedef typename range_type::type                               c_type;
 
  306   typedef typename mtl::traits::range_generator<mtl::tag::nz, c_type>::type ic_type;
 
  309   typename mtl::traits::row<MatrixType>::type                              
row(cpu_matrix);
 
  310   typename mtl::traits::col<MatrixType>::type                              col(cpu_matrix);
 
  311   typename mtl::traits::const_value<MatrixType>::type                      value(cpu_matrix);
 
  314   for (c_type cursor(my_range.begin(cpu_matrix)), cend(my_range.end(cpu_matrix)); cursor != cend; ++cursor)
 
  315     for (ic_type icursor(mtl::begin<mtl::tag::nz>(cursor)), icend(mtl::end<mtl::tag::nz>(cursor)); icursor != icend; ++icursor)
 
  316       stl_matrix[
row(*icursor)][col(*icursor)] = value(*icursor);
 
  318   copy(tools::const_sparse_matrix_adapter<NumericT>(stl_matrix, cpu_matrix.num_rows(), cpu_matrix.num_cols()), gpu_matrix);
 
  340 template<
typename CPUMatrixT, 
typename NumericT, 
unsigned int AlignmentV>
 
  342           CPUMatrixT & cpu_matrix )
 
  347   if ( gpu_matrix.
size1() > 0 && gpu_matrix.
size2() > 0 )
 
  352     std::vector<NumericT> elements(gpu_matrix.
nnz());
 
  364       while (data_index < row_buffer[row])
 
  366         if (col_buffer[data_index] >= gpu_matrix.
size2())
 
  368           std::cerr << 
"ViennaCL encountered invalid data at colbuffer[" << data_index << 
"]: " << col_buffer[data_index] << std::endl;
 
  372         if (std::fabs(elements[data_index]) > 
static_cast<NumericT>(0))
 
  373           cpu_matrix(row-1, static_cast<vcl_size_t>(col_buffer[data_index])) = elements[data_index];
 
  386 template<
typename NumericT, 
unsigned int AlignmentV>
 
  388           std::vector< std::map<unsigned int, NumericT> > & cpu_matrix)
 
  390   assert( (cpu_matrix.size() == gpu_matrix.
size1()) && 
bool(
"Size mismatch") );
 
  393   copy(gpu_matrix, temp);
 
  396 #ifdef VIENNACL_WITH_UBLAS 
  401 template<
typename ScalarType, 
unsigned int AlignmentV, 
typename F, vcl_
size_t IB, 
typename IA, 
typename TA>
 
  403           boost::numeric::ublas::compressed_matrix<ScalarType> & ublas_matrix)
 
  414   ublas_matrix.clear();
 
  415   ublas_matrix.reserve(gpu_matrix.
nnz());
 
  417   ublas_matrix.set_filled(gpu_matrix.
size1() + 1, gpu_matrix.
nnz());
 
  419   for (
vcl_size_t i=0; i<ublas_matrix.size1() + 1; ++i)
 
  420     ublas_matrix.index1_data()[i] = row_buffer[i];
 
  422   for (
vcl_size_t i=0; i<ublas_matrix.nnz(); ++i)
 
  423     ublas_matrix.index2_data()[i] = col_buffer[i];
 
  430 #ifdef VIENNACL_WITH_ARMADILLO 
  436 template<
typename NumericT, 
unsigned int AlignmentV>
 
  438           arma::SpMat<NumericT> & arma_matrix)
 
  440   assert( (static_cast<vcl_size_t>(arma_matrix.n_rows) == vcl_matrix.
size1()) && 
bool(
"Size mismatch") );
 
  441   assert( (static_cast<vcl_size_t>(arma_matrix.n_cols) == vcl_matrix.
size2()) && 
bool(
"Size mismatch") );
 
  443   if ( vcl_matrix.
size1() > 0 && vcl_matrix.
size2() > 0 )
 
  458       while (data_index < row_buffer[row])
 
  460         assert(col_buffer[data_index] < vcl_matrix.
size2() && bool(
"ViennaCL encountered invalid data at col_buffer"));
 
  461         if (elements[data_index] != static_cast<NumericT>(0.0))
 
  462           arma_matrix(row-1, col_buffer[data_index]) = elements[data_index];
 
  470 #ifdef VIENNACL_WITH_EIGEN 
  472 template<
typename NumericT, 
int flags, 
unsigned int AlignmentV>
 
  473 void copy(compressed_matrix<NumericT, AlignmentV> & gpu_matrix,
 
  474           Eigen::SparseMatrix<NumericT, flags> & eigen_matrix)
 
  476   assert( (static_cast<vcl_size_t>(eigen_matrix.rows()) == gpu_matrix.size1()) && 
bool(
"Size mismatch") );
 
  477   assert( (static_cast<vcl_size_t>(eigen_matrix.cols()) == gpu_matrix.size2()) && 
bool(
"Size mismatch") );
 
  479   if ( gpu_matrix.size1() > 0 && gpu_matrix.size2() > 0 )
 
  484     std::vector<NumericT> elements(gpu_matrix.nnz());
 
  490     eigen_matrix.setZero();
 
  494       while (data_index < row_buffer[row])
 
  496         assert(col_buffer[data_index] < gpu_matrix.size2() && bool(
"ViennaCL encountered invalid data at col_buffer"));
 
  497         if (elements[data_index] != static_cast<NumericT>(0.0))
 
  498           eigen_matrix.insert(row-1, col_buffer[data_index]) = elements[data_index];
 
  508 #ifdef VIENNACL_WITH_MTL4 
  510 template<
typename NumericT, 
unsigned int AlignmentV>
 
  511 void copy(compressed_matrix<NumericT, AlignmentV> & gpu_matrix,
 
  512           mtl::compressed2D<NumericT> & mtl4_matrix)
 
  514   assert( (static_cast<vcl_size_t>(mtl4_matrix.num_rows()) == gpu_matrix.size1()) && 
bool(
"Size mismatch") );
 
  515   assert( (static_cast<vcl_size_t>(mtl4_matrix.num_cols()) == gpu_matrix.size2()) && 
bool(
"Size mismatch") );
 
  517   if ( gpu_matrix.size1() > 0 && gpu_matrix.size2() > 0 )
 
  523     std::vector<NumericT> elements(gpu_matrix.nnz());
 
  532     mtl::matrix::inserter< mtl::compressed2D<NumericT> >  ins(mtl4_matrix);
 
  536       while (data_index < row_buffer[row])
 
  538         assert(col_buffer[data_index] < gpu_matrix.size2() && bool(
"ViennaCL encountered invalid data at col_buffer"));
 
  539         if (elements[data_index] != static_cast<NumericT>(0.0))
 
  540           ins(row-1, col_buffer[data_index]) << 
typename mtl::Collection< mtl::compressed2D<NumericT> >::value_type(elements[data_index]);
 
  558 template<
class NumericT, 
unsigned int AlignmentV >
 
  577     : rows_(rows), cols_(cols), nonzeros_(nonzeros), row_block_num_(0)
 
  584 #ifdef VIENNACL_WITH_OPENCL 
  587       row_buffer_.opencl_handle().context(ctx.opencl_context());
 
  588       col_buffer_.opencl_handle().context(ctx.opencl_context());
 
  589       elements_.opencl_handle().context(ctx.opencl_context());
 
  590       row_blocks_.opencl_handle().context(ctx.opencl_context());
 
  613     : rows_(rows), cols_(cols), nonzeros_(0), row_block_num_(0)
 
  620 #ifdef VIENNACL_WITH_OPENCL 
  623       row_buffer_.opencl_handle().context(ctx.opencl_context());
 
  624       col_buffer_.opencl_handle().context(ctx.opencl_context());
 
  625       elements_.opencl_handle().context(ctx.opencl_context());
 
  626       row_blocks_.opencl_handle().context(ctx.opencl_context());
 
  648 #ifdef VIENNACL_WITH_OPENCL 
  651       row_buffer_.opencl_handle().context(ctx.opencl_context());
 
  652       col_buffer_.opencl_handle().context(ctx.opencl_context());
 
  653       elements_.opencl_handle().context(ctx.opencl_context());
 
  654       row_blocks_.opencl_handle().context(ctx.opencl_context());
 
  660 #ifdef VIENNACL_WITH_OPENCL 
  670   explicit compressed_matrix(cl_mem mem_row_buffer, cl_mem mem_col_buffer, cl_mem mem_elements,
 
  672     rows_(rows), cols_(cols), nonzeros_(nonzeros), row_block_num_(0)
 
  675     row_buffer_.opencl_handle() = mem_row_buffer;
 
  676     row_buffer_.opencl_handle().inc();             
 
  677     row_buffer_.
raw_size(
sizeof(cl_uint) * (rows + 1));
 
  680     col_buffer_.opencl_handle() = mem_col_buffer;
 
  681     col_buffer_.opencl_handle().inc();             
 
  682     col_buffer_.
raw_size(
sizeof(cl_uint) * nonzeros);
 
  685     elements_.opencl_handle() = mem_elements;
 
  686     elements_.opencl_handle().inc();               
 
  696     : rows_(0), cols_(0), nonzeros_(0), row_block_num_(0)
 
  705 #ifdef VIENNACL_WITH_OPENCL 
  708       row_buffer_.opencl_handle().context(ctx.opencl_context());
 
  709       col_buffer_.opencl_handle().context(ctx.opencl_context());
 
  710       elements_.opencl_handle().context(ctx.opencl_context());
 
  711       row_blocks_.opencl_handle().context(ctx.opencl_context());
 
  722     assert( (rows_ == 0 || rows_ == other.
size1()) && 
bool(
"Size mismatch") );
 
  723     assert( (cols_ == 0 || cols_ == other.
size2()) && 
bool(
"Size mismatch") );
 
  725     rows_ = other.
size1();
 
  726     cols_ = other.
size2();
 
  727     nonzeros_ = other.
nnz();
 
  728     row_block_num_ = other.row_block_num_;
 
  730     viennacl::backend::typesafe_memory_copy<unsigned int>(other.row_buffer_, row_buffer_);
 
  731     viennacl::backend::typesafe_memory_copy<unsigned int>(other.col_buffer_, col_buffer_);
 
  732     viennacl::backend::typesafe_memory_copy<unsigned int>(other.row_blocks_, row_blocks_);
 
  733     viennacl::backend::typesafe_memory_copy<NumericT>(other.elements_, elements_);
 
  741     assert( (rows_ == 0 || rows_ == proxy.
lhs().size1()) && 
bool(
"Size mismatch") );
 
  742     assert( (cols_ == 0 || cols_ == proxy.
rhs().size2()) && 
bool(
"Size mismatch") );
 
  764   void set(
const void * row_jumper,
 
  765            const void * col_buffer,
 
  771     assert( (rows > 0)     && 
bool(
"Error in compressed_matrix::set(): Number of rows must be larger than zero!"));
 
  772     assert( (cols > 0)     && 
bool(
"Error in compressed_matrix::set(): Number of columns must be larger than zero!"));
 
  773     assert( (nonzeros > 0) && 
bool(
"Error in compressed_matrix::set(): Number of nonzeros must be larger than zero!"));
 
  785     nonzeros_ = nonzeros;
 
  796     if (new_nonzeros > nonzeros_)
 
  819       nonzeros_ = new_nonzeros;
 
  831     assert(new_size1 > 0 && new_size2 > 0 && 
bool(
"Cannot resize to zero size!"));
 
  833     if (new_size1 != rows_ || new_size2 != cols_)
 
  845         std::vector<std::map<unsigned int, NumericT> > stl_sparse_matrix;
 
  848           stl_sparse_matrix.resize(rows_);
 
  851           stl_sparse_matrix.resize(new_size1);
 
  852           stl_sparse_matrix[0][0] = 0;      
 
  855         stl_sparse_matrix.resize(new_size1);
 
  858         if (new_size2 < cols_ && rows_ > 0)
 
  860           for (
vcl_size_t i=0; i<stl_sparse_matrix.size(); ++i)
 
  862             std::list<unsigned int> to_delete;
 
  863             for (
typename std::map<unsigned int, NumericT>::iterator it = stl_sparse_matrix[i].begin();
 
  864                  it != stl_sparse_matrix[i].end();
 
  867               if (it->first >= new_size2)
 
  868                 to_delete.push_back(it->first);
 
  871             for (std::list<unsigned int>::iterator it = to_delete.begin(); it != to_delete.end(); ++it)
 
  872               stl_sparse_matrix[i].erase(*it);
 
  892     std::vector<NumericT> host_elements(1);
 
  904     assert( (i < rows_) && (j < cols_) && 
bool(
"compressed_matrix access out of bounds!"));
 
  909     if (index < nonzeros_)
 
  913     std::vector< std::map<unsigned int, NumericT> > cpu_backup(rows_);
 
  916     cpu_backup[i][
static_cast<unsigned int>(j)] = 0.0;
 
  919     index = element_index(i, j);
 
  921     assert(index < nonzeros_);
 
  959     viennacl::backend::switch_memory_context<unsigned int>(row_buffer_, new_ctx);
 
  960     viennacl::backend::switch_memory_context<unsigned int>(col_buffer_, new_ctx);
 
  961     viennacl::backend::switch_memory_context<unsigned int>(row_blocks_, new_ctx);
 
  962     viennacl::backend::switch_memory_context<NumericT>(elements_, new_ctx);
 
  982     viennacl::backend::memory_read(col_buffer_, col_indices.element_size()*row_indices[0], row_indices.element_size()*col_indices.size(), col_indices.get());
 
  984     for (
vcl_size_t k=0; k<col_indices.size(); ++k)
 
  986       if (col_indices[k] == j)
 
  987         return row_indices[0] + k;
 
 1011     row_blocks.set(0, 0);
 
 1015       num_entries_in_current_batch += entries_in_row;
 
 1017       if (num_entries_in_current_batch > shared_mem_size)
 
 1019         vcl_size_t rows_in_batch = i - row_blocks[row_block_num_];
 
 1020         if (rows_in_batch > 0) 
 
 1021           row_blocks.set(++row_block_num_, i--);
 
 1023           row_blocks.set(++row_block_num_, i+1);
 
 1024         num_entries_in_current_batch = 0;
 
 1027     if (num_entries_in_current_batch > 0)
 
 1028       row_blocks.set(++row_block_num_, rows_);
 
 1030     if (row_block_num_ > 0) 
 
 1032                                        row_blocks.element_size() * (row_block_num_ + 1),
 
 1058 template<
typename NumericT, 
unsigned int AlignmentV>
 
 1059 std::ostream & operator<<(std::ostream & os, compressed_matrix<NumericT, AlignmentV> 
const & A)
 
 1061   std::vector<std::map<unsigned int, NumericT> > tmp(A.size1());
 
 1063   os << 
"compressed_matrix of size (" << A.size1() << 
", " << A.size2() << 
") with " << A.nnz() << 
" nonzeros:" << std::endl;
 
 1067     for (
typename std::map<unsigned int, NumericT>::const_iterator it = tmp[i].begin(); it != tmp[i].end(); ++it)
 
 1068       os << 
"  (" << i << 
", " << it->first << 
")\t" << it->second << std::endl;
 
 1084   template<
typename T, 
unsigned int A>
 
 1085   struct op_executor<vector_base<T>, 
op_assign, vector_expression<const compressed_matrix<T, A>, const vector_base<T>, op_prod> >
 
 1087     static void apply(vector_base<T> & lhs, vector_expression<
const compressed_matrix<T, A>, 
const vector_base<T>, op_prod> 
const & rhs)
 
 1101   template<
typename T, 
unsigned int A>
 
 1102   struct op_executor<vector_base<T>, op_inplace_add, vector_expression<const compressed_matrix<T, A>, const vector_base<T>, op_prod> >
 
 1104     static void apply(vector_base<T> & lhs, vector_expression<
const compressed_matrix<T, A>, 
const vector_base<T>, op_prod> 
const & rhs)
 
 1118   template<
typename T, 
unsigned int A>
 
 1119   struct op_executor<vector_base<T>, op_inplace_sub, vector_expression<const compressed_matrix<T, A>, const vector_base<T>, op_prod> >
 
 1121     static void apply(vector_base<T> & lhs, vector_expression<
const compressed_matrix<T, A>, 
const vector_base<T>, op_prod> 
const & rhs)
 
 1137   template<
typename T, 
unsigned int A, 
typename LHS, 
typename RHS, 
typename OP>
 
 1138   struct op_executor<vector_base<T>, 
op_assign, vector_expression<const compressed_matrix<T, A>, const vector_expression<const LHS, const RHS, OP>, op_prod> >
 
 1140     static void apply(vector_base<T> & lhs, vector_expression<
const compressed_matrix<T, A>, 
const vector_expression<const LHS, const RHS, OP>, op_prod> 
const & rhs)
 
 1148   template<
typename T, 
unsigned int A, 
typename LHS, 
typename RHS, 
typename OP>
 
 1149   struct op_executor<vector_base<T>, op_inplace_add, vector_expression<const compressed_matrix<T, A>, vector_expression<const LHS, const RHS, OP>, op_prod> >
 
 1151     static void apply(vector_base<T> & lhs, vector_expression<
const compressed_matrix<T, A>, vector_expression<const LHS, const RHS, OP>, op_prod> 
const & rhs)
 
 1161   template<
typename T, 
unsigned int A, 
typename LHS, 
typename RHS, 
typename OP>
 
 1162   struct op_executor<vector_base<T>, op_inplace_sub, vector_expression<const compressed_matrix<T, A>, const vector_expression<const LHS, const RHS, OP>, op_prod> >
 
 1164     static void apply(vector_base<T> & lhs, vector_expression<
const compressed_matrix<T, A>, 
const vector_expression<const LHS, const RHS, OP>, op_prod> 
const & rhs)
 
const vcl_size_t & size2() const 
Returns the number of columns. 
Helper class implementing an array on the host. Default case: No conversion necessary. 
vcl_size_t element_size() const 
This class represents a single scalar value on the GPU and behaves mostly like a built-in scalar type...
const vcl_size_t & size1() const 
Returns the number of rows. 
void switch_memory_context(viennacl::context new_ctx)
Switches the memory context of the matrix. 
vcl_size_t size1(MatrixType const &mat)
Generic routine for obtaining the number of rows of a matrix (ViennaCL, uBLAS, etc.) 
A proxy class for entries in a vector. 
Expression template class for representing a tree of expressions which ultimately result in a matrix...
This file provides the forward declarations for the main types used within ViennaCL. 
void memory_read(mem_handle const &src_buffer, vcl_size_t src_offset, vcl_size_t bytes_to_read, void *ptr, bool async=false)
Reads data from a buffer back to main RAM. 
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. 
const vcl_size_t & nnz() const 
Returns the number of nonzero entries. 
vcl_size_t element_size(memory_types)
Represents a generic 'context' similar to an OpenCL context, but is backend-agnostic and thus also su...
void generate_row_block_information()
Builds the row block information needed for fast sparse matrix-vector multiplications. 
vcl_size_t size(VectorType const &vec)
Generic routine for obtaining the size of a vector (ViennaCL, uBLAS, etc.) 
entry_proxy< NumericT > operator()(vcl_size_t i, vcl_size_t j)
Returns a reference to the (i,j)-th entry of the sparse matrix. If (i,j) does not exist (zero)...
handle_type & handle3()
Returns the OpenCL handle to the row block array. 
handle_type & handle()
Returns the OpenCL handle to the matrix entry array. 
Implementations of operations using sparse matrices. 
const handle_type & handle2() const 
Returns the OpenCL handle to the column index array. 
void clear()
Resets all entries in the matrix back to zero without changing the matrix size. Resets the sparsity p...
scalar< typename viennacl::tools::CHECK_SCALAR_TEMPLATE_ARGUMENT< NumericT >::ResultType > value_type
compressed_matrix(vcl_size_t rows, vcl_size_t cols, viennacl::context ctx)
Construction of a compressed matrix with the supplied number of rows and columns. If the number of no...
vector_expression< const matrix_base< NumericT, F >, const unsigned int, op_row > row(const matrix_base< NumericT, F > &A, unsigned int i)
viennacl::memory_types memory_type() const 
void copy_impl(const CPUMatrixT &cpu_matrix, compressed_compressed_matrix< NumericT > &gpu_matrix, vcl_size_t nonzero_rows, vcl_size_t nonzeros)
viennacl::memory_types memory_context() const 
Returns the current memory context to determine whether the matrix is set up for OpenMP, OpenCL, or CUDA. 
RHS & rhs() const 
Get right hand side operand. 
void reserve(vcl_size_t new_nonzeros, bool preserve=true)
Allocate memory for the supplied number of nonzeros in the matrix. Old values are preserved...
const handle_type & handle3() const 
Returns the OpenCL handle to the row block array. 
handle_type & handle1()
Returns the OpenCL handle to the row index array. 
compressed_matrix(vcl_size_t rows, vcl_size_t cols, vcl_size_t nonzeros=0, viennacl::context ctx=viennacl::context())
Construction of a compressed matrix with the supplied number of rows and columns. If the number of no...
void switch_active_handle_id(memory_types new_id)
Switches the currently active handle. If no support for that backend is provided, an exception is thr...
void memory_copy(mem_handle const &src_buffer, mem_handle &dst_buffer, vcl_size_t src_offset, vcl_size_t dst_offset, vcl_size_t bytes_to_copy)
Copies 'bytes_to_copy' bytes from address 'src_buffer + src_offset' to memory starting at address 'ds...
compressed_matrix(viennacl::context ctx)
Creates an empty compressed_matrix, but sets the respective context information. 
viennacl::context context(T const &t)
Returns an ID for the currently active memory domain of an object. 
The vector type with operator-overloads and proxy classes is defined here. Linear algebra operations ...
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) ...
void set(const void *row_jumper, const void *col_buffer, const NumericT *elements, vcl_size_t rows, vcl_size_t cols, vcl_size_t nonzeros)
Sets the row, column and value arrays of the compressed matrix. 
void set(vcl_size_t index, U value)
viennacl::backend::mem_handle handle_type
Main abstraction class for multiple memory domains. Represents a buffer in either main RAM...
vcl_size_t raw_size() const 
Returns the number of bytes of the currently active buffer. 
A sparse square matrix in compressed sparse rows format. 
void memory_create(mem_handle &handle, vcl_size_t size_in_bytes, viennacl::context const &ctx, const void *host_ptr=NULL)
Creates an array of the specified size. If the second argument is provided, the buffer is initialized...
T min(const T &lhs, const T &rhs)
Minimum. 
void prod_impl(const matrix_base< NumericT > &mat, const vector_base< NumericT > &vec, vector_base< NumericT > &result)
Carries out matrix-vector multiplication. 
compressed_matrix(matrix_expression< const compressed_matrix, const compressed_matrix, op_prod > const &proxy)
Assignment a compressed matrix from the product of two compressed_matrix objects (C = A * B)...
viennacl::backend::mem_handle & handle(T &obj)
Returns the generic memory handle of an object. Non-const version. 
LHS & lhs() const 
Get left hand side operand. 
vcl_size_t raw_size() const 
const vcl_size_t & blocks1() const 
Returns the internal number of row blocks for an adaptive SpMV. 
A proxy class for a single element of a vector or matrix. This proxy should not be noticed by end-use...
compressed_matrix()
Default construction of a compressed matrix. No memory is allocated. 
void resize(vcl_size_t new_size1, vcl_size_t new_size2, bool preserve=true)
Resize the matrix. 
void memory_shallow_copy(mem_handle const &src_buffer, mem_handle &dst_buffer)
A 'shallow' copy operation from an initialized buffer to an uninitialized buffer. The uninitialized b...
compressed_matrix & operator=(compressed_matrix const &other)
Assignment a compressed matrix from possibly another memory domain. 
handle_type & handle2()
Returns the OpenCL handle to the column index array. 
compressed_matrix & operator=(matrix_expression< const compressed_matrix, const compressed_matrix, op_prod > const &proxy)
Assignment a compressed matrix from the product of two compressed_matrix objects (C = A * B)...
memory_types get_active_handle_id() const 
Returns an ID for the currently active memory buffer. Other memory buffers might contain old or no da...