spbenchsolver.h
Go to the documentation of this file.
1 // This file is part of Eigen, a lightweight C++ template library
2 // for linear algebra.
3 //
4 // Copyright (C) 2012 Désiré Nuentsa-Wakam <desire.nuentsa_wakam@inria.fr>
5 //
6 // This Source Code Form is subject to the terms of the Mozilla
7 // Public License v. 2.0. If a copy of the MPL was not distributed
8 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
9 
10 
11 #include <iostream>
12 #include <fstream>
13 #include <Eigen/SparseCore>
14 #include <bench/BenchTimer.h>
15 #include <cstdlib>
16 #include <string>
17 #include <Eigen/Cholesky>
18 #include <Eigen/Jacobi>
19 #include <Eigen/Householder>
20 #include <Eigen/IterativeLinearSolvers>
21 #include <unsupported/Eigen/IterativeSolvers>
22 #include <Eigen/LU>
23 #include <unsupported/Eigen/SparseExtra>
24 #include <Eigen/SparseLU>
25 
26 #include "spbenchstyle.h"
27 
28 #ifdef EIGEN_METIS_SUPPORT
29 #include <Eigen/MetisSupport>
30 #endif
31 
32 #ifdef EIGEN_CHOLMOD_SUPPORT
33 #include <Eigen/CholmodSupport>
34 #endif
35 
36 #ifdef EIGEN_UMFPACK_SUPPORT
37 #include <Eigen/UmfPackSupport>
38 #endif
39 
40 #ifdef EIGEN_KLU_SUPPORT
41 #include <Eigen/KLUSupport>
42 #endif
43 
44 #ifdef EIGEN_PARDISO_SUPPORT
45 #include <Eigen/PardisoSupport>
46 #endif
47 
48 #ifdef EIGEN_SUPERLU_SUPPORT
49 #include <Eigen/SuperLUSupport>
50 #endif
51 
52 #ifdef EIGEN_PASTIX_SUPPORT
53 #include <Eigen/PaStiXSupport>
54 #endif
55 
56 // CONSTANTS
57 #define EIGEN_UMFPACK 10
58 #define EIGEN_KLU 11
59 #define EIGEN_SUPERLU 20
60 #define EIGEN_PASTIX 30
61 #define EIGEN_PARDISO 40
62 #define EIGEN_SPARSELU_COLAMD 50
63 #define EIGEN_SPARSELU_METIS 51
64 #define EIGEN_BICGSTAB 60
65 #define EIGEN_BICGSTAB_ILUT 61
66 #define EIGEN_GMRES 70
67 #define EIGEN_GMRES_ILUT 71
68 #define EIGEN_SIMPLICIAL_LDLT 80
69 #define EIGEN_CHOLMOD_LDLT 90
70 #define EIGEN_PASTIX_LDLT 100
71 #define EIGEN_PARDISO_LDLT 110
72 #define EIGEN_SIMPLICIAL_LLT 120
73 #define EIGEN_CHOLMOD_SUPERNODAL_LLT 130
74 #define EIGEN_CHOLMOD_SIMPLICIAL_LLT 140
75 #define EIGEN_PASTIX_LLT 150
76 #define EIGEN_PARDISO_LLT 160
77 #define EIGEN_CG 170
78 #define EIGEN_CG_PRECOND 180
79 
80 using namespace Eigen;
81 using namespace std;
82 
83 
84 // Global variables for input parameters
85 int MaximumIters; // Maximum number of iterations
86 double RelErr; // Relative error of the computed solution
87 double best_time_val; // Current best time overall solvers
88 int best_time_id; // id of the best solver for the current system
89 
90 template<typename T> inline typename NumTraits<T>::Real test_precision() { return NumTraits<T>::dummy_precision(); }
91 template<> inline float test_precision<float>() { return 1e-3f; }
92 template<> inline double test_precision<double>() { return 1e-6; }
93 template<> inline float test_precision<std::complex<float> >() { return test_precision<float>(); }
94 template<> inline double test_precision<std::complex<double> >() { return test_precision<double>(); }
95 
96 void printStatheader(std::ofstream& out)
97 {
98  // Print XML header
99  // NOTE It would have been much easier to write these XML documents using external libraries like tinyXML or Xerces-C++.
100 
101  out << "<?xml version='1.0' encoding='UTF-8'?> \n";
102  out << "<?xml-stylesheet type='text/xsl' href='#stylesheet' ?> \n";
103  out << "<!DOCTYPE BENCH [\n<!ATTLIST xsl:stylesheet\n id\t ID #REQUIRED>\n]>";
104  out << "\n\n<!-- Generated by the Eigen library -->\n";
105 
106  out << "\n<BENCH> \n" ; //root XML element
107  // Print the xsl style section
108  printBenchStyle(out);
109  // List all available solvers
110  out << " <AVAILSOLVER> \n";
111 #ifdef EIGEN_UMFPACK_SUPPORT
112  out <<" <SOLVER ID='" << EIGEN_UMFPACK << "'>\n";
113  out << " <TYPE> LU </TYPE> \n";
114  out << " <PACKAGE> UMFPACK </PACKAGE> \n";
115  out << " </SOLVER> \n";
116 #endif
117 #ifdef EIGEN_KLU_SUPPORT
118  out <<" <SOLVER ID='" << EIGEN_KLU << "'>\n";
119  out << " <TYPE> LU </TYPE> \n";
120  out << " <PACKAGE> KLU </PACKAGE> \n";
121  out << " </SOLVER> \n";
122 #endif
123 #ifdef EIGEN_SUPERLU_SUPPORT
124  out <<" <SOLVER ID='" << EIGEN_SUPERLU << "'>\n";
125  out << " <TYPE> LU </TYPE> \n";
126  out << " <PACKAGE> SUPERLU </PACKAGE> \n";
127  out << " </SOLVER> \n";
128 #endif
129 #ifdef EIGEN_CHOLMOD_SUPPORT
130  out <<" <SOLVER ID='" << EIGEN_CHOLMOD_SIMPLICIAL_LLT << "'>\n";
131  out << " <TYPE> LLT SP</TYPE> \n";
132  out << " <PACKAGE> CHOLMOD </PACKAGE> \n";
133  out << " </SOLVER> \n";
134 
135  out <<" <SOLVER ID='" << EIGEN_CHOLMOD_SUPERNODAL_LLT << "'>\n";
136  out << " <TYPE> LLT</TYPE> \n";
137  out << " <PACKAGE> CHOLMOD </PACKAGE> \n";
138  out << " </SOLVER> \n";
139 
140  out <<" <SOLVER ID='" << EIGEN_CHOLMOD_LDLT << "'>\n";
141  out << " <TYPE> LDLT </TYPE> \n";
142  out << " <PACKAGE> CHOLMOD </PACKAGE> \n";
143  out << " </SOLVER> \n";
144 #endif
145 #ifdef EIGEN_PARDISO_SUPPORT
146  out <<" <SOLVER ID='" << EIGEN_PARDISO << "'>\n";
147  out << " <TYPE> LU </TYPE> \n";
148  out << " <PACKAGE> PARDISO </PACKAGE> \n";
149  out << " </SOLVER> \n";
150 
151  out <<" <SOLVER ID='" << EIGEN_PARDISO_LLT << "'>\n";
152  out << " <TYPE> LLT </TYPE> \n";
153  out << " <PACKAGE> PARDISO </PACKAGE> \n";
154  out << " </SOLVER> \n";
155 
156  out <<" <SOLVER ID='" << EIGEN_PARDISO_LDLT << "'>\n";
157  out << " <TYPE> LDLT </TYPE> \n";
158  out << " <PACKAGE> PARDISO </PACKAGE> \n";
159  out << " </SOLVER> \n";
160 #endif
161 #ifdef EIGEN_PASTIX_SUPPORT
162  out <<" <SOLVER ID='" << EIGEN_PASTIX << "'>\n";
163  out << " <TYPE> LU </TYPE> \n";
164  out << " <PACKAGE> PASTIX </PACKAGE> \n";
165  out << " </SOLVER> \n";
166 
167  out <<" <SOLVER ID='" << EIGEN_PASTIX_LLT << "'>\n";
168  out << " <TYPE> LLT </TYPE> \n";
169  out << " <PACKAGE> PASTIX </PACKAGE> \n";
170  out << " </SOLVER> \n";
171 
172  out <<" <SOLVER ID='" << EIGEN_PASTIX_LDLT << "'>\n";
173  out << " <TYPE> LDLT </TYPE> \n";
174  out << " <PACKAGE> PASTIX </PACKAGE> \n";
175  out << " </SOLVER> \n";
176 #endif
177 
178  out <<" <SOLVER ID='" << EIGEN_BICGSTAB << "'>\n";
179  out << " <TYPE> BICGSTAB </TYPE> \n";
180  out << " <PACKAGE> EIGEN </PACKAGE> \n";
181  out << " </SOLVER> \n";
182 
183  out <<" <SOLVER ID='" << EIGEN_BICGSTAB_ILUT << "'>\n";
184  out << " <TYPE> BICGSTAB_ILUT </TYPE> \n";
185  out << " <PACKAGE> EIGEN </PACKAGE> \n";
186  out << " </SOLVER> \n";
187 
188  out <<" <SOLVER ID='" << EIGEN_GMRES_ILUT << "'>\n";
189  out << " <TYPE> GMRES_ILUT </TYPE> \n";
190  out << " <PACKAGE> EIGEN </PACKAGE> \n";
191  out << " </SOLVER> \n";
192 
193  out <<" <SOLVER ID='" << EIGEN_SIMPLICIAL_LDLT << "'>\n";
194  out << " <TYPE> LDLT </TYPE> \n";
195  out << " <PACKAGE> EIGEN </PACKAGE> \n";
196  out << " </SOLVER> \n";
197 
198  out <<" <SOLVER ID='" << EIGEN_SIMPLICIAL_LLT << "'>\n";
199  out << " <TYPE> LLT </TYPE> \n";
200  out << " <PACKAGE> EIGEN </PACKAGE> \n";
201  out << " </SOLVER> \n";
202 
203  out <<" <SOLVER ID='" << EIGEN_CG << "'>\n";
204  out << " <TYPE> CG </TYPE> \n";
205  out << " <PACKAGE> EIGEN </PACKAGE> \n";
206  out << " </SOLVER> \n";
207 
208  out <<" <SOLVER ID='" << EIGEN_SPARSELU_COLAMD << "'>\n";
209  out << " <TYPE> LU_COLAMD </TYPE> \n";
210  out << " <PACKAGE> EIGEN </PACKAGE> \n";
211  out << " </SOLVER> \n";
212 
213 #ifdef EIGEN_METIS_SUPPORT
214  out <<" <SOLVER ID='" << EIGEN_SPARSELU_METIS << "'>\n";
215  out << " <TYPE> LU_METIS </TYPE> \n";
216  out << " <PACKAGE> EIGEN </PACKAGE> \n";
217  out << " </SOLVER> \n";
218 #endif
219  out << " </AVAILSOLVER> \n";
220 
221 }
222 
223 
224 template<typename Solver, typename Scalar>
225 void call_solver(Solver &solver, const int solver_id, const typename Solver::MatrixType& A, const Matrix<Scalar, Dynamic, 1>& b, const Matrix<Scalar, Dynamic, 1>& refX,std::ofstream& statbuf)
226 {
227 
228  double total_time;
229  double compute_time;
230  double solve_time;
231  double rel_error;
233  BenchTimer timer;
234  timer.reset();
235  timer.start();
236  solver.compute(A);
237  if (solver.info() != Success)
238  {
239  std::cerr << "Solver failed ... \n";
240  return;
241  }
242  timer.stop();
243  compute_time = timer.value();
244  statbuf << " <TIME>\n";
245  statbuf << " <COMPUTE> " << timer.value() << "</COMPUTE>\n";
246  std::cout<< "COMPUTE TIME : " << timer.value() <<std::endl;
247 
248  timer.reset();
249  timer.start();
250  x = solver.solve(b);
251  if (solver.info() == NumericalIssue)
252  {
253  std::cerr << "Solver failed ... \n";
254  return;
255  }
256  timer.stop();
257  solve_time = timer.value();
258  statbuf << " <SOLVE> " << timer.value() << "</SOLVE>\n";
259  std::cout<< "SOLVE TIME : " << timer.value() <<std::endl;
260 
261  total_time = solve_time + compute_time;
262  statbuf << " <TOTAL> " << total_time << "</TOTAL>\n";
263  std::cout<< "TOTAL TIME : " << total_time <<std::endl;
264  statbuf << " </TIME>\n";
265 
266  // Verify the relative error
267  if(refX.size() != 0)
268  rel_error = (refX - x).norm()/refX.norm();
269  else
270  {
271  // Compute the relative residual norm
273  temp = A * x;
274  rel_error = (b-temp).norm()/b.norm();
275  }
276  statbuf << " <ERROR> " << rel_error << "</ERROR>\n";
277  std::cout<< "REL. ERROR : " << rel_error << "\n\n" ;
278  if ( rel_error <= RelErr )
279  {
280  // check the best time if convergence
281  if(!best_time_val || (best_time_val > total_time))
282  {
283  best_time_val = total_time;
284  best_time_id = solver_id;
285  }
286  }
287 }
288 
289 template<typename Solver, typename Scalar>
290 void call_directsolver(Solver& solver, const int solver_id, const typename Solver::MatrixType& A, const Matrix<Scalar, Dynamic, 1>& b, const Matrix<Scalar, Dynamic, 1>& refX, std::string& statFile)
291 {
292  std::ofstream statbuf(statFile.c_str(), std::ios::app);
293  statbuf << " <SOLVER_STAT ID='" << solver_id <<"'>\n";
294  call_solver(solver, solver_id, A, b, refX,statbuf);
295  statbuf << " </SOLVER_STAT>\n";
296  statbuf.close();
297 }
298 
299 template<typename Solver, typename Scalar>
300 void call_itersolver(Solver &solver, const int solver_id, const typename Solver::MatrixType& A, const Matrix<Scalar, Dynamic, 1>& b, const Matrix<Scalar, Dynamic, 1>& refX, std::string& statFile)
301 {
302  solver.setTolerance(RelErr);
303  solver.setMaxIterations(MaximumIters);
304 
305  std::ofstream statbuf(statFile.c_str(), std::ios::app);
306  statbuf << " <SOLVER_STAT ID='" << solver_id <<"'>\n";
307  call_solver(solver, solver_id, A, b, refX,statbuf);
308  statbuf << " <ITER> "<< solver.iterations() << "</ITER>\n";
309  statbuf << " </SOLVER_STAT>\n";
310  std::cout << "ITERATIONS : " << solver.iterations() <<"\n\n\n";
311 
312 }
313 
314 
315 template <typename Scalar>
316 void SelectSolvers(const SparseMatrix<Scalar>&A, unsigned int sym, Matrix<Scalar, Dynamic, 1>& b, const Matrix<Scalar, Dynamic, 1>& refX, std::string& statFile)
317 {
319  // First, deal with Nonsymmetric and symmetric matrices
320  best_time_id = 0;
321  best_time_val = 0.0;
322  //UMFPACK
323  #ifdef EIGEN_UMFPACK_SUPPORT
324  {
325  cout << "Solving with UMFPACK LU ... \n";
327  call_directsolver(solver, EIGEN_UMFPACK, A, b, refX,statFile);
328  }
329  #endif
330  //KLU
331  #ifdef EIGEN_KLU_SUPPORT
332  {
333  cout << "Solving with KLU LU ... \n";
335  call_directsolver(solver, EIGEN_KLU, A, b, refX,statFile);
336  }
337  #endif
338  //SuperLU
339  #ifdef EIGEN_SUPERLU_SUPPORT
340  {
341  cout << "\nSolving with SUPERLU ... \n";
343  call_directsolver(solver, EIGEN_SUPERLU, A, b, refX,statFile);
344  }
345  #endif
346 
347  // PaStix LU
348  #ifdef EIGEN_PASTIX_SUPPORT
349  {
350  cout << "\nSolving with PASTIX LU ... \n";
352  call_directsolver(solver, EIGEN_PASTIX, A, b, refX,statFile) ;
353  }
354  #endif
355 
356  //PARDISO LU
357  #ifdef EIGEN_PARDISO_SUPPORT
358  {
359  cout << "\nSolving with PARDISO LU ... \n";
361  call_directsolver(solver, EIGEN_PARDISO, A, b, refX,statFile);
362  }
363  #endif
364 
365  // Eigen SparseLU METIS
366  cout << "\n Solving with Sparse LU AND COLAMD ... \n";
368  call_directsolver(solver, EIGEN_SPARSELU_COLAMD, A, b, refX, statFile);
369  // Eigen SparseLU METIS
370  #ifdef EIGEN_METIS_SUPPORT
371  {
372  cout << "\n Solving with Sparse LU AND METIS ... \n";
374  call_directsolver(solver, EIGEN_SPARSELU_METIS, A, b, refX, statFile);
375  }
376  #endif
377 
378  //BiCGSTAB
379  {
380  cout << "\nSolving with BiCGSTAB ... \n";
382  call_itersolver(solver, EIGEN_BICGSTAB, A, b, refX,statFile);
383  }
384  //BiCGSTAB+ILUT
385  {
386  cout << "\nSolving with BiCGSTAB and ILUT ... \n";
388  call_itersolver(solver, EIGEN_BICGSTAB_ILUT, A, b, refX,statFile);
389  }
390 
391 
392  //GMRES
393 // {
394 // cout << "\nSolving with GMRES ... \n";
395 // GMRES<SpMat> solver;
396 // call_itersolver(solver, EIGEN_GMRES, A, b, refX,statFile);
397 // }
398  //GMRES+ILUT
399  {
400  cout << "\nSolving with GMRES and ILUT ... \n";
402  call_itersolver(solver, EIGEN_GMRES_ILUT, A, b, refX,statFile);
403  }
404 
405  // Hermitian and not necessarily positive-definites
406  if (sym != NonSymmetric)
407  {
408  // Internal Cholesky
409  {
410  cout << "\nSolving with Simplicial LDLT ... \n";
412  call_directsolver(solver, EIGEN_SIMPLICIAL_LDLT, A, b, refX,statFile);
413  }
414 
415  // CHOLMOD
416  #ifdef EIGEN_CHOLMOD_SUPPORT
417  {
418  cout << "\nSolving with CHOLMOD LDLT ... \n";
420  solver.setMode(CholmodLDLt);
421  call_directsolver(solver,EIGEN_CHOLMOD_LDLT, A, b, refX,statFile);
422  }
423  #endif
424 
425  //PASTIX LLT
426  #ifdef EIGEN_PASTIX_SUPPORT
427  {
428  cout << "\nSolving with PASTIX LDLT ... \n";
430  call_directsolver(solver,EIGEN_PASTIX_LDLT, A, b, refX,statFile);
431  }
432  #endif
433 
434  //PARDISO LLT
435  #ifdef EIGEN_PARDISO_SUPPORT
436  {
437  cout << "\nSolving with PARDISO LDLT ... \n";
439  call_directsolver(solver,EIGEN_PARDISO_LDLT, A, b, refX,statFile);
440  }
441  #endif
442  }
443 
444  // Now, symmetric POSITIVE DEFINITE matrices
445  if (sym == SPD)
446  {
447 
448  //Internal Sparse Cholesky
449  {
450  cout << "\nSolving with SIMPLICIAL LLT ... \n";
452  call_directsolver(solver,EIGEN_SIMPLICIAL_LLT, A, b, refX,statFile);
453  }
454 
455  // CHOLMOD
456  #ifdef EIGEN_CHOLMOD_SUPPORT
457  {
458  // CholMOD SuperNodal LLT
459  cout << "\nSolving with CHOLMOD LLT (Supernodal)... \n";
462  call_directsolver(solver,EIGEN_CHOLMOD_SUPERNODAL_LLT, A, b, refX,statFile);
463  // CholMod Simplicial LLT
464  cout << "\nSolving with CHOLMOD LLT (Simplicial) ... \n";
466  call_directsolver(solver,EIGEN_CHOLMOD_SIMPLICIAL_LLT, A, b, refX,statFile);
467  }
468  #endif
469 
470  //PASTIX LLT
471  #ifdef EIGEN_PASTIX_SUPPORT
472  {
473  cout << "\nSolving with PASTIX LLT ... \n";
475  call_directsolver(solver,EIGEN_PASTIX_LLT, A, b, refX,statFile);
476  }
477  #endif
478 
479  //PARDISO LLT
480  #ifdef EIGEN_PARDISO_SUPPORT
481  {
482  cout << "\nSolving with PARDISO LLT ... \n";
484  call_directsolver(solver,EIGEN_PARDISO_LLT, A, b, refX,statFile);
485  }
486  #endif
487 
488  // Internal CG
489  {
490  cout << "\nSolving with CG ... \n";
492  call_itersolver(solver,EIGEN_CG, A, b, refX,statFile);
493  }
494  //CG+IdentityPreconditioner
495 // {
496 // cout << "\nSolving with CG and IdentityPreconditioner ... \n";
497 // ConjugateGradient<SpMat, Lower, IdentityPreconditioner> solver;
498 // call_itersolver(solver,EIGEN_CG_PRECOND, A, b, refX,statFile);
499 // }
500  } // End SPD matrices
501 }
502 
503 /* Browse all the matrices available in the specified folder
504  * and solve the associated linear system.
505  * The results of each solve are printed in the standard output
506  * and optionally in the provided html file
507  */
508 template <typename Scalar>
509 void Browse_Matrices(const string folder, bool statFileExists, std::string& statFile, int maxiters, double tol)
510 {
511  MaximumIters = maxiters; // Maximum number of iterations, global variable
512  RelErr = tol; //Relative residual error as stopping criterion for iterative solvers
513  MatrixMarketIterator<Scalar> it(folder);
514  for ( ; it; ++it)
515  {
516  //print the infos for this linear system
517  if(statFileExists)
518  {
519  std::ofstream statbuf(statFile.c_str(), std::ios::app);
520  statbuf << "<LINEARSYSTEM> \n";
521  statbuf << " <MATRIX> \n";
522  statbuf << " <NAME> " << it.matname() << " </NAME>\n";
523  statbuf << " <SIZE> " << it.matrix().rows() << " </SIZE>\n";
524  statbuf << " <ENTRIES> " << it.matrix().nonZeros() << "</ENTRIES>\n";
525  if (it.sym()!=NonSymmetric)
526  {
527  statbuf << " <SYMMETRY> Symmetric </SYMMETRY>\n" ;
528  if (it.sym() == SPD)
529  statbuf << " <POSDEF> YES </POSDEF>\n";
530  else
531  statbuf << " <POSDEF> NO </POSDEF>\n";
532 
533  }
534  else
535  {
536  statbuf << " <SYMMETRY> NonSymmetric </SYMMETRY>\n" ;
537  statbuf << " <POSDEF> NO </POSDEF>\n";
538  }
539  statbuf << " </MATRIX> \n";
540  statbuf.close();
541  }
542 
543  cout<< "\n\n===================================================== \n";
544  cout<< " ====== SOLVING WITH MATRIX " << it.matname() << " ====\n";
545  cout<< " =================================================== \n\n";
547  if(it.hasrefX()) refX = it.refX();
548  // Call all suitable solvers for this linear system
549  SelectSolvers<Scalar>(it.matrix(), it.sym(), it.rhs(), refX, statFile);
550 
551  if(statFileExists)
552  {
553  std::ofstream statbuf(statFile.c_str(), std::ios::app);
554  statbuf << " <BEST_SOLVER ID='"<< best_time_id
555  << "'></BEST_SOLVER>\n";
556  statbuf << " </LINEARSYSTEM> \n";
557  statbuf.close();
558  }
559  }
560 }
561 
562 bool get_options(int argc, char **args, string option, string* value=0)
563 {
564  int idx = 1, found=false;
565  while (idx<argc && !found){
566  if (option.compare(args[idx]) == 0){
567  found = true;
568  if(value) *value = args[idx+1];
569  }
570  idx+=2;
571  }
572  return found;
573 }
#define EIGEN_SIMPLICIAL_LDLT
Definition: spbenchsolver.h:68
void setMode(CholmodMode mode)
void printStatheader(std::ofstream &out)
Definition: spbenchsolver.h:96
void SelectSolvers(const SparseMatrix< Scalar > &A, unsigned int sym, Matrix< Scalar, Dynamic, 1 > &b, const Matrix< Scalar, Dynamic, 1 > &refX, std::string &statFile)
#define EIGEN_PARDISO
Definition: spbenchsolver.h:61
A sparse LU factorization and solver based on UmfPack.
Scalar * b
Definition: benchVecAdd.cpp:17
double test_precision< double >()
Definition: main.h:417
void Browse_Matrices(const string folder, bool statFileExists, std::string &statFile, int maxiters, double tol)
#define EIGEN_CHOLMOD_LDLT
Definition: spbenchsolver.h:69
#define EIGEN_BICGSTAB_ILUT
Definition: spbenchsolver.h:65
A sparse direct LU factorization and solver based on the SuperLU library.
Index rows() const
Definition: SparseMatrix.h:138
std::ofstream out("Result.txt")
bool get_options(int argc, char **args, string option, string *value=0)
void printBenchStyle(std::ofstream &out)
Definition: spbenchstyle.h:13
Definition: pytypes.h:2012
double value(int TIMER=CPU_TIMER) const
Definition: BenchTimer.h:104
Namespace containing all symbols from the Eigen library.
Definition: jet.h:637
A conjugate gradient solver for sparse (or dense) self-adjoint problems.
Definition: BFloat16.h:88
MatrixXf MatrixType
BiCGSTAB< SparseMatrix< double > > solver
Holds information about the various numeric (i.e. scalar) types allowed by Eigen. ...
Definition: NumTraits.h:232
#define EIGEN_GMRES_ILUT
Definition: spbenchsolver.h:67
void call_itersolver(Solver &solver, const int solver_id, const typename Solver::MatrixType &A, const Matrix< Scalar, Dynamic, 1 > &b, const Matrix< Scalar, Dynamic, 1 > &refX, std::string &statFile)
A direct sparse LDLT Cholesky factorizations without square root.
Sparse supernodal LU factorization for general matrices.
Definition: SparseLU.h:17
#define EIGEN_SPARSELU_METIS
Definition: spbenchsolver.h:63
A sparse direct Cholesky (LLT) factorization and solver based on the PARDISO library.
float test_precision< float >()
Definition: main.h:416
#define EIGEN_PARDISO_LLT
Definition: spbenchsolver.h:76
A sparse direct supernodal Cholesky (LLT) factorization and solver based on the PaStiX library...
Definition: PaStiXSupport.h:33
int best_time_id
Definition: spbenchsolver.h:88
int MaximumIters
Definition: spbenchsolver.h:85
double best_time_val
Definition: spbenchsolver.h:87
#define EIGEN_PASTIX
Definition: spbenchsolver.h:60
Point2(* f)(const Point3 &, OptionalJacobian< 2, 3 >)
Array< double, 1, 3 > e(1./3., 0.5, 2.)
#define EIGEN_SUPERLU
Definition: spbenchsolver.h:59
Iterator to browse matrices from a specified folder.
A sparse direct LU factorization and solver based on the PARDISO library.
A general Cholesky factorization and solver based on Cholmod.
#define EIGEN_KLU
Definition: spbenchsolver.h:58
#define EIGEN_UMFPACK
Definition: spbenchsolver.h:57
Eigen::SparseMatrix< double > SpMat
double RelErr
Definition: spbenchsolver.h:86
#define EIGEN_SIMPLICIAL_LLT
Definition: spbenchsolver.h:72
NumTraits< T >::Real test_precision()
Definition: main.h:415
A GMRES solver for sparse square problems.
Definition: GMRES.h:221
#define EIGEN_PASTIX_LLT
Definition: spbenchsolver.h:75
#define EIGEN_CHOLMOD_SIMPLICIAL_LLT
Definition: spbenchsolver.h:74
A sparse direct Cholesky (LDLT) factorization and solver based on the PARDISO library.
#define EIGEN_BICGSTAB
Definition: spbenchsolver.h:64
Interface to the PaStix solver.
Definition: PaStiXSupport.h:31
#define EIGEN_SPARSELU_COLAMD
Definition: spbenchsolver.h:62
A sparse direct supernodal Cholesky (LLT) factorization and solver based on the PaStiX library...
Definition: PaStiXSupport.h:32
const G double tol
Definition: Group.h:86
A bi conjugate gradient stabilized solver for sparse square problems.
Definition: BiCGSTAB.h:113
#define EIGEN_PASTIX_LDLT
Definition: spbenchsolver.h:70
A direct sparse LLT Cholesky factorizations.
#define EIGEN_CHOLMOD_SUPERNODAL_LLT
Definition: spbenchsolver.h:73
set noclip points set clip one set noclip two set bar set border lt lw set xdata set ydata set zdata set x2data set y2data set boxwidth set dummy x
#define EIGEN_PARDISO_LDLT
Definition: spbenchsolver.h:71
void call_directsolver(Solver &solver, const int solver_id, const typename Solver::MatrixType &A, const Matrix< Scalar, Dynamic, 1 > &b, const Matrix< Scalar, Dynamic, 1 > &refX, std::string &statFile)
static BenchTimer timer
#define EIGEN_CG
Definition: spbenchsolver.h:77
void call_solver(Solver &solver, const int solver_id, const typename Solver::MatrixType &A, const Matrix< Scalar, Dynamic, 1 > &b, const Matrix< Scalar, Dynamic, 1 > &refX, std::ofstream &statbuf)


gtsam
Author(s):
autogenerated on Tue Jul 4 2023 02:36:18