Home | History | Annotate | Download | only in bench 
      1 
      2 #include <Eigen/Sparse>
      3 #include <bench/BenchTimer.h>
      4 #include <set>
      5 
      6 using namespace std;
      7 using namespace Eigen;
      8 using namespace Eigen;
      9 
     10 #ifndef SIZE
     11 #define SIZE 1024
     12 #endif
     13 
     14 #ifndef DENSITY
     15 #define DENSITY 0.01
     16 #endif
     17 
     18 #ifndef SCALAR
     19 #define SCALAR double
     20 #endif
     21 
     22 typedef SCALAR Scalar;
     23 typedef Matrix<Scalar,Dynamic,Dynamic> DenseMatrix;
     24 typedef Matrix<Scalar,Dynamic,1> DenseVector;
     25 typedef SparseMatrix<Scalar> EigenSparseMatrix;
     26 
     27 void fillMatrix(float density, int rows, int cols,  EigenSparseMatrix& dst)
     28 {
     29   dst.reserve(double(rows)*cols*density);
     30   for(int j = 0; j < cols; j++)
     31   {
     32     for(int i = 0; i < rows; i++)
     33     {
     34       Scalar v = (internal::random<float>(0,1) < density) ? internal::random<Scalar>() : 0;
     35       if (v!=0)
     36         dst.insert(i,j) = v;
     37     }
     38   }
     39   dst.finalize();
     40 }
     41 
     42 void fillMatrix2(int nnzPerCol, int rows, int cols,  EigenSparseMatrix& dst)
     43 {
     44 //   std::cout << "alloc " << nnzPerCol*cols << "\n";
     45   dst.reserve(nnzPerCol*cols);
     46   for(int j = 0; j < cols; j++)
     47   {
     48     std::set<int> aux;
     49     for(int i = 0; i < nnzPerCol; i++)
     50     {
     51       int k = internal::random<int>(0,rows-1);
     52       while (aux.find(k)!=aux.end())
     53         k = internal::random<int>(0,rows-1);
     54       aux.insert(k);
     55 
     56       dst.insert(k,j) = internal::random<Scalar>();
     57     }
     58   }
     59   dst.finalize();
     60 }
     61 
     62 void eiToDense(const EigenSparseMatrix& src, DenseMatrix& dst)
     63 {
     64   dst.setZero();
     65   for (int j=0; j<src.cols(); ++j)
     66     for (EigenSparseMatrix::InnerIterator it(src.derived(), j); it; ++it)
     67       dst(it.index(),j) = it.value();
     68 }
     69 
     70 #ifndef NOGMM
     71 #include "gmm/gmm.h"
     72 typedef gmm::csc_matrix<Scalar> GmmSparse;
     73 typedef gmm::col_matrix< gmm::wsvector<Scalar> > GmmDynSparse;
     74 void eiToGmm(const EigenSparseMatrix& src, GmmSparse& dst)
     75 {
     76   GmmDynSparse tmp(src.rows(), src.cols());
     77   for (int j=0; j<src.cols(); ++j)
     78     for (EigenSparseMatrix::InnerIterator it(src.derived(), j); it; ++it)
     79       tmp(it.index(),j) = it.value();
     80   gmm::copy(tmp, dst);
     81 }
     82 #endif
     83 
     84 #ifndef NOMTL
     85 #include <boost/numeric/mtl/mtl.hpp>
     86 typedef mtl::compressed2D<Scalar, mtl::matrix::parameters<mtl::tag::col_major> > MtlSparse;
     87 typedef mtl::compressed2D<Scalar, mtl::matrix::parameters<mtl::tag::row_major> > MtlSparseRowMajor;
     88 void eiToMtl(const EigenSparseMatrix& src, MtlSparse& dst)
     89 {
     90   mtl::matrix::inserter<MtlSparse> ins(dst);
     91   for (int j=0; j<src.cols(); ++j)
     92     for (EigenSparseMatrix::InnerIterator it(src.derived(), j); it; ++it)
     93       ins[it.index()][j] = it.value();
     94 }
     95 #endif
     96 
     97 #ifdef CSPARSE
     98 extern "C" {
     99 #include "cs.h"
    100 }
    101 void eiToCSparse(const EigenSparseMatrix& src, cs* &dst)
    102 {
    103   cs* aux = cs_spalloc (0, 0, 1, 1, 1);
    104   for (int j=0; j<src.cols(); ++j)
    105     for (EigenSparseMatrix::InnerIterator it(src.derived(), j); it; ++it)
    106       if (!cs_entry(aux, it.index(), j, it.value()))
    107       {
    108         std::cout << "cs_entry error\n";
    109         exit(2);
    110       }
    111    dst = cs_compress(aux);
    112 //    cs_spfree(aux);
    113 }
    114 #endif // CSPARSE
    115 
    116 #ifndef NOUBLAS
    117 #include <boost/numeric/ublas/vector.hpp>
    118 #include <boost/numeric/ublas/matrix.hpp>
    119 #include <boost/numeric/ublas/io.hpp>
    120 #include <boost/numeric/ublas/triangular.hpp>
    121 #include <boost/numeric/ublas/vector_sparse.hpp>
    122 #include <boost/numeric/ublas/matrix_sparse.hpp>
    123 #include <boost/numeric/ublas/vector_of_vector.hpp>
    124 #include <boost/numeric/ublas/operation.hpp>
    125 
    126 typedef boost::numeric::ublas::compressed_matrix<Scalar,boost::numeric::ublas::column_major> UBlasSparse;
    127 
    128 void eiToUblas(const EigenSparseMatrix& src, UBlasSparse& dst)
    129 {
    130   dst.resize(src.rows(), src.cols(), false);
    131   for (int j=0; j<src.cols(); ++j)
    132     for (EigenSparseMatrix::InnerIterator it(src.derived(), j); it; ++it)
    133       dst(it.index(),j) = it.value();
    134 }
    135 
    136 template <typename EigenType, typename UblasType>
    137 void eiToUblasVec(const EigenType& src, UblasType& dst)
    138 {
    139   dst.resize(src.size());
    140   for (int j=0; j<src.size(); ++j)
    141       dst[j] = src.coeff(j);
    142 }
    143 #endif
    144 
    145 #ifdef OSKI
    146 extern "C" {
    147 #include <oski/oski.h>
    148 }
    149 #endif
    150