sparse_lu.cpp
Go to the documentation of this file.
00001 
00002 // g++ -I.. sparse_lu.cpp -O3 -g0 -I /usr/include/superlu/ -lsuperlu -lgfortran -DSIZE=1000 -DDENSITY=.05 && ./a.out
00003 
00004 #define EIGEN_SUPERLU_SUPPORT
00005 #define EIGEN_UMFPACK_SUPPORT
00006 #include <Eigen/Sparse>
00007 
00008 #define NOGMM
00009 #define NOMTL
00010 
00011 #ifndef SIZE
00012 #define SIZE 10
00013 #endif
00014 
00015 #ifndef DENSITY
00016 #define DENSITY 0.01
00017 #endif
00018 
00019 #ifndef REPEAT
00020 #define REPEAT 1
00021 #endif
00022 
00023 #include "BenchSparseUtil.h"
00024 
00025 #ifndef MINDENSITY
00026 #define MINDENSITY 0.0004
00027 #endif
00028 
00029 #ifndef NBTRIES
00030 #define NBTRIES 10
00031 #endif
00032 
00033 #define BENCH(X) \
00034   timer.reset(); \
00035   for (int _j=0; _j<NBTRIES; ++_j) { \
00036     timer.start(); \
00037     for (int _k=0; _k<REPEAT; ++_k) { \
00038         X  \
00039   } timer.stop(); }
00040 
00041 typedef Matrix<Scalar,Dynamic,1> VectorX;
00042 
00043 #include <Eigen/LU>
00044 
00045 template<int Backend>
00046 void doEigen(const char* name, const EigenSparseMatrix& sm1, const VectorX& b, VectorX& x, int flags = 0)
00047 {
00048   std::cout << name << "..." << std::flush;
00049   BenchTimer timer; timer.start();
00050   SparseLU<EigenSparseMatrix,Backend> lu(sm1, flags);
00051   timer.stop();
00052   if (lu.succeeded())
00053     std::cout << ":\t" << timer.value() << endl;
00054   else
00055   {
00056     std::cout << ":\t FAILED" << endl;
00057     return;
00058   }
00059 
00060   bool ok;
00061   timer.reset(); timer.start();
00062   ok = lu.solve(b,&x);
00063   timer.stop();
00064   if (ok)
00065     std::cout << "  solve:\t" << timer.value() << endl;
00066   else
00067     std::cout << "  solve:\t" << " FAILED" << endl;
00068 
00069   //std::cout << x.transpose() << "\n";
00070 }
00071 
00072 int main(int argc, char *argv[])
00073 {
00074   int rows = SIZE;
00075   int cols = SIZE;
00076   float density = DENSITY;
00077   BenchTimer timer;
00078 
00079   VectorX b = VectorX::Random(cols);
00080   VectorX x = VectorX::Random(cols);
00081 
00082   bool densedone = false;
00083 
00084   //for (float density = DENSITY; density>=MINDENSITY; density*=0.5)
00085 //   float density = 0.5;
00086   {
00087     EigenSparseMatrix sm1(rows, cols);
00088     fillMatrix(density, rows, cols, sm1);
00089 
00090     // dense matrices
00091     #ifdef DENSEMATRIX
00092     if (!densedone)
00093     {
00094       densedone = true;
00095       std::cout << "Eigen Dense\t" << density*100 << "%\n";
00096       DenseMatrix m1(rows,cols);
00097       eiToDense(sm1, m1);
00098 
00099       BenchTimer timer;
00100       timer.start();
00101       FullPivLU<DenseMatrix> lu(m1);
00102       timer.stop();
00103       std::cout << "Eigen/dense:\t" << timer.value() << endl;
00104 
00105       timer.reset();
00106       timer.start();
00107       lu.solve(b,&x);
00108       timer.stop();
00109       std::cout << "  solve:\t" << timer.value() << endl;
00110 //       std::cout << b.transpose() << "\n";
00111 //       std::cout << x.transpose() << "\n";
00112     }
00113     #endif
00114 
00115     #ifdef EIGEN_UMFPACK_SUPPORT
00116     x.setZero();
00117     doEigen<Eigen::UmfPack>("Eigen/UmfPack (auto)", sm1, b, x, 0);
00118     #endif
00119 
00120     #ifdef EIGEN_SUPERLU_SUPPORT
00121     x.setZero();
00122     doEigen<Eigen::SuperLU>("Eigen/SuperLU (nat)", sm1, b, x, Eigen::NaturalOrdering);
00123 //     doEigen<Eigen::SuperLU>("Eigen/SuperLU (MD AT+A)", sm1, b, x, Eigen::MinimumDegree_AT_PLUS_A);
00124 //     doEigen<Eigen::SuperLU>("Eigen/SuperLU (MD ATA)", sm1, b, x, Eigen::MinimumDegree_ATA);
00125     doEigen<Eigen::SuperLU>("Eigen/SuperLU (COLAMD)", sm1, b, x, Eigen::ColApproxMinimumDegree);
00126     #endif
00127 
00128   }
00129 
00130   return 0;
00131 }
00132 


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