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