30 #include <Eigen/Sparse> 
   33 #define VIENNACL_WITH_EIGEN 1 
   53 struct Eigen_dense_matrix
 
   55   typedef typename T::ERROR_NO_EIGEN_TYPE_AVAILABLE   error_type;
 
   59 struct Eigen_dense_matrix<float>
 
   61   typedef Eigen::MatrixXf  type;
 
   65 struct Eigen_dense_matrix<double>
 
   67   typedef Eigen::MatrixXd  type;
 
   75   typedef typename T::ERROR_NO_EIGEN_TYPE_AVAILABLE   error_type;
 
   79 struct Eigen_vector<float>
 
   81   typedef Eigen::VectorXf  type;
 
   85 struct Eigen_vector<double>
 
   87   typedef Eigen::VectorXd  type;
 
  102 template<
typename ScalarType>
 
  109   typedef typename Eigen_dense_matrix<ScalarType>::type  EigenMatrix;
 
  110   typedef typename Eigen_vector<ScalarType>::type        EigenVector;
 
  115   EigenMatrix eigen_densemat(6, 5);
 
  116   EigenMatrix eigen_densemat2(6, 5);
 
  117   eigen_densemat(0,0) = 2.0;   eigen_densemat(0,1) = -1.0;
 
  118   eigen_densemat(1,0) = -1.0;  eigen_densemat(1,1) =  2.0;  eigen_densemat(1,2) = -1.0;
 
  119   eigen_densemat(2,1) = -1.0;  eigen_densemat(2,2) = -1.0;  eigen_densemat(2,3) = -1.0;
 
  120   eigen_densemat(3,2) = -1.0;  eigen_densemat(3,3) =  2.0;  eigen_densemat(3,4) = -1.0;
 
  121                                eigen_densemat(5,4) = -1.0;  eigen_densemat(4,4) = -1.0;
 
  122   Eigen::Map<EigenMatrix> eigen_densemat_map(eigen_densemat.data(), 6, 5); 
 
  127   Eigen::SparseMatrix<ScalarType, Eigen::RowMajor> eigen_sparsemat(6, 5);
 
  128   Eigen::SparseMatrix<ScalarType, Eigen::RowMajor> eigen_sparsemat2(6, 5);
 
  129   eigen_sparsemat.reserve(5*2);
 
  130   eigen_sparsemat.insert(0,0) = 2.0;   eigen_sparsemat.insert(0,1) = -1.0;
 
  131   eigen_sparsemat.insert(1,1) = 2.0;   eigen_sparsemat.insert(1,2) = -1.0;
 
  132   eigen_sparsemat.insert(2,2) = -1.0;  eigen_sparsemat.insert(2,3) = -1.0;
 
  133   eigen_sparsemat.insert(3,3) = 2.0;   eigen_sparsemat.insert(3,4) = -1.0;
 
  134   eigen_sparsemat.insert(5,4) = -1.0;
 
  140   EigenVector eigen_rhs(5);
 
  141   Eigen::Map<EigenVector> eigen_rhs_map(eigen_rhs.data(), 5);
 
  142   EigenVector eigen_result(6);
 
  143   EigenVector eigen_temp(6);
 
  164   viennacl::copy(&(eigen_rhs[0]), &(eigen_rhs[0]) + 5, vcl_rhs.begin());  
 
  171   std::cout << 
"VCL sparsematrix dimensions: " << vcl_sparsemat.size1() << 
", " << vcl_sparsemat.size2() << std::endl;
 
  181   eigen_result = eigen_densemat * eigen_rhs;
 
  184   std::cout << 
"Difference for dense matrix-vector product: " << (eigen_result - eigen_temp).norm() << std::endl;
 
  185   std::cout << 
"Difference for dense matrix-vector product (Eigen->ViennaCL->Eigen): " 
  186             << (eigen_densemat2 * eigen_rhs - eigen_temp).norm() << std::endl;
 
  191   eigen_result = eigen_sparsemat * eigen_rhs;
 
  194   std::cout << 
"Difference for sparse matrix-vector product: " << (eigen_result - eigen_temp).norm() << std::endl;
 
  195   std::cout << 
"Difference for sparse matrix-vector product (Eigen->ViennaCL->Eigen): " 
  196             << (eigen_sparsemat2 * eigen_rhs - eigen_temp).norm() << std::endl;
 
  203 int main(
int, 
char *[])
 
  205   std::cout << 
"----------------------------------------------" << std::endl;
 
  206   std::cout << 
"## Single precision" << std::endl;
 
  207   std::cout << 
"----------------------------------------------" << std::endl;
 
  208   run_tutorial<float>();
 
  210 #ifdef VIENNACL_HAVE_OPENCL 
  214     std::cout << 
"----------------------------------------------" << std::endl;
 
  215     std::cout << 
"## Double precision" << std::endl;
 
  216     std::cout << 
"----------------------------------------------" << std::endl;
 
  217     run_tutorial<double>();
 
  223   std::cout << std::endl;
 
  224   std::cout << 
"!!!! TUTORIAL COMPLETED SUCCESSFULLY !!!!" << std::endl;
 
  225   std::cout << std::endl;
 
Generic interface for matrix-vector and matrix-matrix products. See viennacl/linalg/vector_operations...
Implementation of the dense matrix class. 
viennacl::ocl::device const & current_device()
Convenience function for returning the active device in the current context. 
VectorT prod(std::vector< std::vector< T, A1 >, A2 > const &matrix, VectorT const &vector)
Implementation of the compressed_matrix class. 
bool double_support() const 
ViennaCL convenience function: Returns true if the device supports double precision. 
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) ...
A sparse square matrix in compressed sparse rows format.