sparse_dense_product.cpp
Go to the documentation of this file.
00001 
00002 //g++ -O3 -g0 -DNDEBUG  sparse_product.cpp -I.. -I/home/gael/Coding/LinearAlgebra/mtl4/ -DDENSITY=0.005 -DSIZE=10000 && ./a.out
00003 //g++ -O3 -g0 -DNDEBUG  sparse_product.cpp -I.. -I/home/gael/Coding/LinearAlgebra/mtl4/ -DDENSITY=0.05 -DSIZE=2000 && ./a.out
00004 // -DNOGMM -DNOMTL -DCSPARSE
00005 // -I /home/gael/Coding/LinearAlgebra/CSparse/Include/ /home/gael/Coding/LinearAlgebra/CSparse/Lib/libcsparse.a
00006 #ifndef SIZE
00007 #define SIZE 10000
00008 #endif
00009 
00010 #ifndef DENSITY
00011 #define DENSITY 0.01
00012 #endif
00013 
00014 #ifndef REPEAT
00015 #define REPEAT 1
00016 #endif
00017 
00018 #include "BenchSparseUtil.h"
00019 
00020 #ifndef MINDENSITY
00021 #define MINDENSITY 0.0004
00022 #endif
00023 
00024 #ifndef NBTRIES
00025 #define NBTRIES 10
00026 #endif
00027 
00028 #define BENCH(X) \
00029   timer.reset(); \
00030   for (int _j=0; _j<NBTRIES; ++_j) { \
00031     timer.start(); \
00032     for (int _k=0; _k<REPEAT; ++_k) { \
00033         X  \
00034   } timer.stop(); }
00035 
00036 
00037 #ifdef CSPARSE
00038 cs* cs_sorted_multiply(const cs* a, const cs* b)
00039 {
00040   cs* A = cs_transpose (a, 1) ;
00041   cs* B = cs_transpose (b, 1) ;
00042   cs* D = cs_multiply (B,A) ;   /* D = B'*A' */
00043   cs_spfree (A) ;
00044   cs_spfree (B) ;
00045   cs_dropzeros (D) ;      /* drop zeros from D */
00046   cs* C = cs_transpose (D, 1) ;   /* C = D', so that C is sorted */
00047   cs_spfree (D) ;
00048   return C;
00049 }
00050 #endif
00051 
00052 int main(int argc, char *argv[])
00053 {
00054   int rows = SIZE;
00055   int cols = SIZE;
00056   float density = DENSITY;
00057 
00058   EigenSparseMatrix sm1(rows,cols);
00059   DenseVector v1(cols), v2(cols);
00060   v1.setRandom();
00061 
00062   BenchTimer timer;
00063   for (float density = DENSITY; density>=MINDENSITY; density*=0.5)
00064   {
00065     fillMatrix(density, rows, cols, sm1);
00066 
00067     // dense matrices
00068     #ifdef DENSEMATRIX
00069     {
00070       std::cout << "Eigen Dense\t" << density*100 << "%\n";
00071       DenseMatrix m1(rows,cols);
00072       eiToDense(sm1, m1);
00073 
00074       timer.reset();
00075       timer.start();
00076       for (int k=0; k<REPEAT; ++k)
00077         v2 = m1 * v1;
00078       timer.stop();
00079       std::cout << "   a * v:\t" << timer.value() << endl;
00080 
00081       timer.reset();
00082       timer.start();
00083       for (int k=0; k<REPEAT; ++k)
00084         v2 = m1.transpose() * v1;
00085       timer.stop();
00086       std::cout << "   a' * v:\t" << timer.value() << endl;
00087     }
00088     #endif
00089 
00090     // eigen sparse matrices
00091     {
00092       std::cout << "Eigen sparse\t" << sm1.nonZeros()/float(sm1.rows()*sm1.cols())*100 << "%\n";
00093 
00094       BENCH(asm("#myc"); v2 = sm1 * v1; asm("#myd");)
00095       std::cout << "   a * v:\t" << timer.value() << endl;
00096 
00097 
00098       BENCH( { asm("#mya"); v2 = sm1.transpose() * v1; asm("#myb"); })
00099 
00100       std::cout << "   a' * v:\t" << timer.value() << endl;
00101     }
00102 
00103 //     {
00104 //       DynamicSparseMatrix<Scalar> m1(sm1);
00105 //       std::cout << "Eigen dyn-sparse\t" << m1.nonZeros()/float(m1.rows()*m1.cols())*100 << "%\n";
00106 //
00107 //       BENCH(for (int k=0; k<REPEAT; ++k) v2 = m1 * v1;)
00108 //       std::cout << "   a * v:\t" << timer.value() << endl;
00109 //
00110 //       BENCH(for (int k=0; k<REPEAT; ++k) v2 = m1.transpose() * v1;)
00111 //       std::cout << "   a' * v:\t" << timer.value() << endl;
00112 //     }
00113 
00114     // GMM++
00115     #ifndef NOGMM
00116     {
00117       std::cout << "GMM++ sparse\t" << density*100 << "%\n";
00118       //GmmDynSparse  gmmT3(rows,cols);
00119       GmmSparse m1(rows,cols);
00120       eiToGmm(sm1, m1);
00121 
00122       std::vector<Scalar> gmmV1(cols), gmmV2(cols);
00123       Map<Matrix<Scalar,Dynamic,1> >(&gmmV1[0], cols) = v1;
00124       Map<Matrix<Scalar,Dynamic,1> >(&gmmV2[0], cols) = v2;
00125 
00126       BENCH( asm("#myx"); gmm::mult(m1, gmmV1, gmmV2); asm("#myy"); )
00127       std::cout << "   a * v:\t" << timer.value() << endl;
00128 
00129       BENCH( gmm::mult(gmm::transposed(m1), gmmV1, gmmV2); )
00130       std::cout << "   a' * v:\t" << timer.value() << endl;
00131     }
00132     #endif
00133 
00134     // MTL4
00135     #ifndef NOMTL
00136     {
00137       std::cout << "MTL4\t" << density*100 << "%\n";
00138       MtlSparse m1(rows,cols);
00139       eiToMtl(sm1, m1);
00140       mtl::dense_vector<Scalar> mtlV1(cols, 1.0);
00141       mtl::dense_vector<Scalar> mtlV2(cols, 1.0);
00142 
00143       timer.reset();
00144       timer.start();
00145       for (int k=0; k<REPEAT; ++k)
00146         mtlV2 = m1 * mtlV1;
00147       timer.stop();
00148       std::cout << "   a * v:\t" << timer.value() << endl;
00149 
00150       timer.reset();
00151       timer.start();
00152       for (int k=0; k<REPEAT; ++k)
00153         mtlV2 = trans(m1) * mtlV1;
00154       timer.stop();
00155       std::cout << "   a' * v:\t" << timer.value() << endl;
00156     }
00157     #endif
00158 
00159     std::cout << "\n\n";
00160   }
00161 
00162   return 0;
00163 }
00164 


libicr
Author(s): Robert Krug
autogenerated on Mon Jan 6 2014 11:33:32