Eigen  3.2.2
 All Classes Namespaces Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
SparseLU.h
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 // Copyright (C) 2012 Gael Guennebaud <gael.guennebaud@inria.fr>
6 //
7 // This Source Code Form is subject to the terms of the Mozilla
8 // Public License v. 2.0. If a copy of the MPL was not distributed
9 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
10 
11 
12 #ifndef EIGEN_SPARSE_LU_H
13 #define EIGEN_SPARSE_LU_H
14 
15 namespace Eigen {
16 
17 template <typename _MatrixType, typename _OrderingType = COLAMDOrdering<typename _MatrixType::Index> > class SparseLU;
18 template <typename MappedSparseMatrixType> struct SparseLUMatrixLReturnType;
19 template <typename MatrixLType, typename MatrixUType> struct SparseLUMatrixUReturnType;
20 
72 template <typename _MatrixType, typename _OrderingType>
73 class SparseLU : public internal::SparseLUImpl<typename _MatrixType::Scalar, typename _MatrixType::Index>
74 {
75  public:
76  typedef _MatrixType MatrixType;
77  typedef _OrderingType OrderingType;
78  typedef typename MatrixType::Scalar Scalar;
79  typedef typename MatrixType::RealScalar RealScalar;
80  typedef typename MatrixType::Index Index;
82  typedef internal::MappedSuperNodalMatrix<Scalar, Index> SCMatrix;
86  typedef internal::SparseLUImpl<Scalar, Index> Base;
87 
88  public:
89  SparseLU():m_isInitialized(true),m_lastError(""),m_Ustore(0,0,0,0,0,0),m_symmetricmode(false),m_diagpivotthresh(1.0),m_detPermR(1)
90  {
91  initperfvalues();
92  }
93  SparseLU(const MatrixType& matrix):m_isInitialized(true),m_lastError(""),m_Ustore(0,0,0,0,0,0),m_symmetricmode(false),m_diagpivotthresh(1.0),m_detPermR(1)
94  {
95  initperfvalues();
96  compute(matrix);
97  }
98 
99  ~SparseLU()
100  {
101  // Free all explicit dynamic pointers
102  }
103 
104  void analyzePattern (const MatrixType& matrix);
105  void factorize (const MatrixType& matrix);
106  void simplicialfactorize(const MatrixType& matrix);
107 
112  void compute (const MatrixType& matrix)
113  {
114  // Analyze
115  analyzePattern(matrix);
116  //Factorize
117  factorize(matrix);
118  }
119 
120  inline Index rows() const { return m_mat.rows(); }
121  inline Index cols() const { return m_mat.cols(); }
123  void isSymmetric(bool sym)
124  {
125  m_symmetricmode = sym;
126  }
127 
134  SparseLUMatrixLReturnType<SCMatrix> matrixL() const
135  {
136  return SparseLUMatrixLReturnType<SCMatrix>(m_Lstore);
137  }
144  SparseLUMatrixUReturnType<SCMatrix,MappedSparseMatrix<Scalar,ColMajor,Index> > matrixU() const
145  {
146  return SparseLUMatrixUReturnType<SCMatrix, MappedSparseMatrix<Scalar,ColMajor,Index> >(m_Lstore, m_Ustore);
147  }
148 
153  inline const PermutationType& rowsPermutation() const
154  {
155  return m_perm_r;
156  }
161  inline const PermutationType& colsPermutation() const
162  {
163  return m_perm_c;
164  }
166  void setPivotThreshold(const RealScalar& thresh)
167  {
168  m_diagpivotthresh = thresh;
169  }
170 
177  template<typename Rhs>
178  inline const internal::solve_retval<SparseLU, Rhs> solve(const MatrixBase<Rhs>& B) const
179  {
180  eigen_assert(m_factorizationIsOk && "SparseLU is not initialized.");
181  eigen_assert(rows()==B.rows()
182  && "SparseLU::solve(): invalid number of rows of the right hand side matrix B");
183  return internal::solve_retval<SparseLU, Rhs>(*this, B.derived());
184  }
185 
190  template<typename Rhs>
191  inline const internal::sparse_solve_retval<SparseLU, Rhs> solve(const SparseMatrixBase<Rhs>& B) const
192  {
193  eigen_assert(m_factorizationIsOk && "SparseLU is not initialized.");
194  eigen_assert(rows()==B.rows()
195  && "SparseLU::solve(): invalid number of rows of the right hand side matrix B");
196  return internal::sparse_solve_retval<SparseLU, Rhs>(*this, B.derived());
197  }
198 
208  {
209  eigen_assert(m_isInitialized && "Decomposition is not initialized.");
210  return m_info;
211  }
212 
216  std::string lastErrorMessage() const
217  {
218  return m_lastError;
219  }
220 
221  template<typename Rhs, typename Dest>
222  bool _solve(const MatrixBase<Rhs> &B, MatrixBase<Dest> &X_base) const
223  {
224  Dest& X(X_base.derived());
225  eigen_assert(m_factorizationIsOk && "The matrix should be factorized first");
226  EIGEN_STATIC_ASSERT((Dest::Flags&RowMajorBit)==0,
227  THIS_METHOD_IS_ONLY_FOR_COLUMN_MAJOR_MATRICES);
228 
229  // Permute the right hand side to form X = Pr*B
230  // on return, X is overwritten by the computed solution
231  X.resize(B.rows(),B.cols());
232 
233  // this ugly const_cast_derived() helps to detect aliasing when applying the permutations
234  for(Index j = 0; j < B.cols(); ++j)
235  X.col(j) = rowsPermutation() * B.const_cast_derived().col(j);
236 
237  //Forward substitution with L
238  this->matrixL().solveInPlace(X);
239  this->matrixU().solveInPlace(X);
240 
241  // Permute back the solution
242  for (Index j = 0; j < B.cols(); ++j)
243  X.col(j) = colsPermutation().inverse() * X.col(j);
244 
245  return true;
246  }
247 
258  Scalar absDeterminant()
259  {
260  eigen_assert(m_factorizationIsOk && "The matrix should be factorized first.");
261  // Initialize with the determinant of the row matrix
262  Scalar det = Scalar(1.);
263  //Note that the diagonal blocks of U are stored in supernodes,
264  // which are available in the L part :)
265  for (Index j = 0; j < this->cols(); ++j)
266  {
267  for (typename SCMatrix::InnerIterator it(m_Lstore, j); it; ++it)
268  {
269  if(it.row() < j) continue;
270  if(it.row() == j)
271  {
272  det *= (std::abs)(it.value());
273  break;
274  }
275  }
276  }
277  return det;
278  }
279 
288  Scalar logAbsDeterminant() const
289  {
290  eigen_assert(m_factorizationIsOk && "The matrix should be factorized first.");
291  Scalar det = Scalar(0.);
292  for (Index j = 0; j < this->cols(); ++j)
293  {
294  for (typename SCMatrix::InnerIterator it(m_Lstore, j); it; ++it)
295  {
296  if(it.row() < j) continue;
297  if(it.row() == j)
298  {
299  det += (std::log)((std::abs)(it.value()));
300  break;
301  }
302  }
303  }
304  return det;
305  }
306 
312  {
313  eigen_assert(m_factorizationIsOk && "The matrix should be factorized first.");
314  return Scalar(m_detPermR);
315  }
316 
317  protected:
318  // Functions
319  void initperfvalues()
320  {
321  m_perfv.panel_size = 1;
322  m_perfv.relax = 1;
323  m_perfv.maxsuper = 128;
324  m_perfv.rowblk = 16;
325  m_perfv.colblk = 8;
326  m_perfv.fillfactor = 20;
327  }
328 
329  // Variables
330  mutable ComputationInfo m_info;
331  bool m_isInitialized;
332  bool m_factorizationIsOk;
333  bool m_analysisIsOk;
334  std::string m_lastError;
335  NCMatrix m_mat; // The input (permuted ) matrix
336  SCMatrix m_Lstore; // The lower triangular matrix (supernodal)
337  MappedSparseMatrix<Scalar,ColMajor,Index> m_Ustore; // The upper triangular matrix
338  PermutationType m_perm_c; // Column permutation
339  PermutationType m_perm_r ; // Row permutation
340  IndexVector m_etree; // Column elimination tree
341 
342  typename Base::GlobalLU_t m_glu;
343 
344  // SparseLU options
345  bool m_symmetricmode;
346  // values for performance
347  internal::perfvalues<Index> m_perfv;
348  RealScalar m_diagpivotthresh; // Specifies the threshold used for a diagonal entry to be an acceptable pivot
349  Index m_nnzL, m_nnzU; // Nonzeros in L and U factors
350  Index m_detPermR; // Determinant of the coefficient matrix
351  private:
352  // Disable copy constructor
353  SparseLU (const SparseLU& );
354 
355 }; // End class SparseLU
356 
357 
358 
359 // Functions needed by the anaysis phase
370 template <typename MatrixType, typename OrderingType>
372 {
373 
374  //TODO It is possible as in SuperLU to compute row and columns scaling vectors to equilibrate the matrix mat.
375 
376  OrderingType ord;
377  ord(mat,m_perm_c);
378 
379  // Apply the permutation to the column of the input matrix
380  //First copy the whole input matrix.
381  m_mat = mat;
382  if (m_perm_c.size()) {
383  m_mat.uncompress(); //NOTE: The effect of this command is only to create the InnerNonzeros pointers. FIXME : This vector is filled but not subsequently used.
384  //Then, permute only the column pointers
385  const Index * outerIndexPtr;
386  if (mat.isCompressed()) outerIndexPtr = mat.outerIndexPtr();
387  else
388  {
389  Index *outerIndexPtr_t = new Index[mat.cols()+1];
390  for(Index i = 0; i <= mat.cols(); i++) outerIndexPtr_t[i] = m_mat.outerIndexPtr()[i];
391  outerIndexPtr = outerIndexPtr_t;
392  }
393  for (Index i = 0; i < mat.cols(); i++)
394  {
395  m_mat.outerIndexPtr()[m_perm_c.indices()(i)] = outerIndexPtr[i];
396  m_mat.innerNonZeroPtr()[m_perm_c.indices()(i)] = outerIndexPtr[i+1] - outerIndexPtr[i];
397  }
398  if(!mat.isCompressed()) delete[] outerIndexPtr;
399  }
400  // Compute the column elimination tree of the permuted matrix
401  IndexVector firstRowElt;
402  internal::coletree(m_mat, m_etree,firstRowElt);
403 
404  // In symmetric mode, do not do postorder here
405  if (!m_symmetricmode) {
406  IndexVector post, iwork;
407  // Post order etree
408  internal::treePostorder(m_mat.cols(), m_etree, post);
409 
410 
411  // Renumber etree in postorder
412  Index m = m_mat.cols();
413  iwork.resize(m+1);
414  for (Index i = 0; i < m; ++i) iwork(post(i)) = post(m_etree(i));
415  m_etree = iwork;
416 
417  // Postmultiply A*Pc by post, i.e reorder the matrix according to the postorder of the etree
418  PermutationType post_perm(m);
419  for (Index i = 0; i < m; i++)
420  post_perm.indices()(i) = post(i);
421 
422  // Combine the two permutations : postorder the permutation for future use
423  if(m_perm_c.size()) {
424  m_perm_c = post_perm * m_perm_c;
425  }
426 
427  } // end postordering
428 
429  m_analysisIsOk = true;
430 }
431 
432 // Functions needed by the numerical factorization phase
433 
434 
453 template <typename MatrixType, typename OrderingType>
454 void SparseLU<MatrixType, OrderingType>::factorize(const MatrixType& matrix)
455 {
456  using internal::emptyIdxLU;
457  eigen_assert(m_analysisIsOk && "analyzePattern() should be called first");
458  eigen_assert((matrix.rows() == matrix.cols()) && "Only for squared matrices");
459 
460  typedef typename IndexVector::Scalar Index;
461 
462 
463  // Apply the column permutation computed in analyzepattern()
464  // m_mat = matrix * m_perm_c.inverse();
465  m_mat = matrix;
466  if (m_perm_c.size())
467  {
468  m_mat.uncompress(); //NOTE: The effect of this command is only to create the InnerNonzeros pointers.
469  //Then, permute only the column pointers
470  const Index * outerIndexPtr;
471  if (matrix.isCompressed()) outerIndexPtr = matrix.outerIndexPtr();
472  else
473  {
474  Index* outerIndexPtr_t = new Index[matrix.cols()+1];
475  for(Index i = 0; i <= matrix.cols(); i++) outerIndexPtr_t[i] = m_mat.outerIndexPtr()[i];
476  outerIndexPtr = outerIndexPtr_t;
477  }
478  for (Index i = 0; i < matrix.cols(); i++)
479  {
480  m_mat.outerIndexPtr()[m_perm_c.indices()(i)] = outerIndexPtr[i];
481  m_mat.innerNonZeroPtr()[m_perm_c.indices()(i)] = outerIndexPtr[i+1] - outerIndexPtr[i];
482  }
483  if(!matrix.isCompressed()) delete[] outerIndexPtr;
484  }
485  else
486  { //FIXME This should not be needed if the empty permutation is handled transparently
487  m_perm_c.resize(matrix.cols());
488  for(Index i = 0; i < matrix.cols(); ++i) m_perm_c.indices()(i) = i;
489  }
490 
491  Index m = m_mat.rows();
492  Index n = m_mat.cols();
493  Index nnz = m_mat.nonZeros();
494  Index maxpanel = m_perfv.panel_size * m;
495  // Allocate working storage common to the factor routines
496  Index lwork = 0;
497  Index info = Base::memInit(m, n, nnz, lwork, m_perfv.fillfactor, m_perfv.panel_size, m_glu);
498  if (info)
499  {
500  m_lastError = "UNABLE TO ALLOCATE WORKING MEMORY\n\n" ;
501  m_factorizationIsOk = false;
502  return ;
503  }
504 
505  // Set up pointers for integer working arrays
506  IndexVector segrep(m); segrep.setZero();
507  IndexVector parent(m); parent.setZero();
508  IndexVector xplore(m); xplore.setZero();
509  IndexVector repfnz(maxpanel);
510  IndexVector panel_lsub(maxpanel);
511  IndexVector xprune(n); xprune.setZero();
512  IndexVector marker(m*internal::LUNoMarker); marker.setZero();
513 
514  repfnz.setConstant(-1);
515  panel_lsub.setConstant(-1);
516 
517  // Set up pointers for scalar working arrays
518  ScalarVector dense;
519  dense.setZero(maxpanel);
520  ScalarVector tempv;
521  tempv.setZero(internal::LUnumTempV(m, m_perfv.panel_size, m_perfv.maxsuper, /*m_perfv.rowblk*/m) );
522 
523  // Compute the inverse of perm_c
524  PermutationType iperm_c(m_perm_c.inverse());
525 
526  // Identify initial relaxed snodes
527  IndexVector relax_end(n);
528  if ( m_symmetricmode == true )
529  Base::heap_relax_snode(n, m_etree, m_perfv.relax, marker, relax_end);
530  else
531  Base::relax_snode(n, m_etree, m_perfv.relax, marker, relax_end);
532 
533 
534  m_perm_r.resize(m);
535  m_perm_r.indices().setConstant(-1);
536  marker.setConstant(-1);
537  m_detPermR = 1; // Record the determinant of the row permutation
538 
539  m_glu.supno(0) = emptyIdxLU; m_glu.xsup.setConstant(0);
540  m_glu.xsup(0) = m_glu.xlsub(0) = m_glu.xusub(0) = m_glu.xlusup(0) = Index(0);
541 
542  // Work on one 'panel' at a time. A panel is one of the following :
543  // (a) a relaxed supernode at the bottom of the etree, or
544  // (b) panel_size contiguous columns, <panel_size> defined by the user
545  Index jcol;
546  IndexVector panel_histo(n);
547  Index pivrow; // Pivotal row number in the original row matrix
548  Index nseg1; // Number of segments in U-column above panel row jcol
549  Index nseg; // Number of segments in each U-column
550  Index irep;
551  Index i, k, jj;
552  for (jcol = 0; jcol < n; )
553  {
554  // Adjust panel size so that a panel won't overlap with the next relaxed snode.
555  Index panel_size = m_perfv.panel_size; // upper bound on panel width
556  for (k = jcol + 1; k < (std::min)(jcol+panel_size, n); k++)
557  {
558  if (relax_end(k) != emptyIdxLU)
559  {
560  panel_size = k - jcol;
561  break;
562  }
563  }
564  if (k == n)
565  panel_size = n - jcol;
566 
567  // Symbolic outer factorization on a panel of columns
568  Base::panel_dfs(m, panel_size, jcol, m_mat, m_perm_r.indices(), nseg1, dense, panel_lsub, segrep, repfnz, xprune, marker, parent, xplore, m_glu);
569 
570  // Numeric sup-panel updates in topological order
571  Base::panel_bmod(m, panel_size, jcol, nseg1, dense, tempv, segrep, repfnz, m_glu);
572 
573  // Sparse LU within the panel, and below the panel diagonal
574  for ( jj = jcol; jj< jcol + panel_size; jj++)
575  {
576  k = (jj - jcol) * m; // Column index for w-wide arrays
577 
578  nseg = nseg1; // begin after all the panel segments
579  //Depth-first-search for the current column
580  VectorBlock<IndexVector> panel_lsubk(panel_lsub, k, m);
581  VectorBlock<IndexVector> repfnz_k(repfnz, k, m);
582  info = Base::column_dfs(m, jj, m_perm_r.indices(), m_perfv.maxsuper, nseg, panel_lsubk, segrep, repfnz_k, xprune, marker, parent, xplore, m_glu);
583  if ( info )
584  {
585  m_lastError = "UNABLE TO EXPAND MEMORY IN COLUMN_DFS() ";
586  m_info = NumericalIssue;
587  m_factorizationIsOk = false;
588  return;
589  }
590  // Numeric updates to this column
591  VectorBlock<ScalarVector> dense_k(dense, k, m);
592  VectorBlock<IndexVector> segrep_k(segrep, nseg1, m-nseg1);
593  info = Base::column_bmod(jj, (nseg - nseg1), dense_k, tempv, segrep_k, repfnz_k, jcol, m_glu);
594  if ( info )
595  {
596  m_lastError = "UNABLE TO EXPAND MEMORY IN COLUMN_BMOD() ";
597  m_info = NumericalIssue;
598  m_factorizationIsOk = false;
599  return;
600  }
601 
602  // Copy the U-segments to ucol(*)
603  info = Base::copy_to_ucol(jj, nseg, segrep, repfnz_k ,m_perm_r.indices(), dense_k, m_glu);
604  if ( info )
605  {
606  m_lastError = "UNABLE TO EXPAND MEMORY IN COPY_TO_UCOL() ";
607  m_info = NumericalIssue;
608  m_factorizationIsOk = false;
609  return;
610  }
611 
612  // Form the L-segment
613  info = Base::pivotL(jj, m_diagpivotthresh, m_perm_r.indices(), iperm_c.indices(), pivrow, m_glu);
614  if ( info )
615  {
616  m_lastError = "THE MATRIX IS STRUCTURALLY SINGULAR ... ZERO COLUMN AT ";
617  std::ostringstream returnInfo;
618  returnInfo << info;
619  m_lastError += returnInfo.str();
620  m_info = NumericalIssue;
621  m_factorizationIsOk = false;
622  return;
623  }
624 
625  // Update the determinant of the row permutation matrix
626  if (pivrow != jj) m_detPermR *= -1;
627 
628  // Prune columns (0:jj-1) using column jj
629  Base::pruneL(jj, m_perm_r.indices(), pivrow, nseg, segrep, repfnz_k, xprune, m_glu);
630 
631  // Reset repfnz for this column
632  for (i = 0; i < nseg; i++)
633  {
634  irep = segrep(i);
635  repfnz_k(irep) = emptyIdxLU;
636  }
637  } // end SparseLU within the panel
638  jcol += panel_size; // Move to the next panel
639  } // end for -- end elimination
640 
641  // Count the number of nonzeros in factors
642  Base::countnz(n, m_nnzL, m_nnzU, m_glu);
643  // Apply permutation to the L subscripts
644  Base::fixupL(n, m_perm_r.indices(), m_glu);
645 
646  // Create supernode matrix L
647  m_Lstore.setInfos(m, n, m_glu.lusup, m_glu.xlusup, m_glu.lsub, m_glu.xlsub, m_glu.supno, m_glu.xsup);
648  // Create the column major upper sparse matrix U;
649  new (&m_Ustore) MappedSparseMatrix<Scalar, ColMajor, Index> ( m, n, m_nnzU, m_glu.xusub.data(), m_glu.usub.data(), m_glu.ucol.data() );
650 
651  m_info = Success;
652  m_factorizationIsOk = true;
653 }
654 
655 template<typename MappedSupernodalType>
656 struct SparseLUMatrixLReturnType : internal::no_assignment_operator
657 {
658  typedef typename MappedSupernodalType::Index Index;
659  typedef typename MappedSupernodalType::Scalar Scalar;
660  SparseLUMatrixLReturnType(const MappedSupernodalType& mapL) : m_mapL(mapL)
661  { }
662  Index rows() { return m_mapL.rows(); }
663  Index cols() { return m_mapL.cols(); }
664  template<typename Dest>
665  void solveInPlace( MatrixBase<Dest> &X) const
666  {
667  m_mapL.solveInPlace(X);
668  }
669  const MappedSupernodalType& m_mapL;
670 };
671 
672 template<typename MatrixLType, typename MatrixUType>
673 struct SparseLUMatrixUReturnType : internal::no_assignment_operator
674 {
675  typedef typename MatrixLType::Index Index;
676  typedef typename MatrixLType::Scalar Scalar;
677  SparseLUMatrixUReturnType(const MatrixLType& mapL, const MatrixUType& mapU)
678  : m_mapL(mapL),m_mapU(mapU)
679  { }
680  Index rows() { return m_mapL.rows(); }
681  Index cols() { return m_mapL.cols(); }
682 
683  template<typename Dest> void solveInPlace(MatrixBase<Dest> &X) const
684  {
685  Index nrhs = X.cols();
686  Index n = X.rows();
687  // Backward solve with U
688  for (Index k = m_mapL.nsuper(); k >= 0; k--)
689  {
690  Index fsupc = m_mapL.supToCol()[k];
691  Index lda = m_mapL.colIndexPtr()[fsupc+1] - m_mapL.colIndexPtr()[fsupc]; // leading dimension
692  Index nsupc = m_mapL.supToCol()[k+1] - fsupc;
693  Index luptr = m_mapL.colIndexPtr()[fsupc];
694 
695  if (nsupc == 1)
696  {
697  for (Index j = 0; j < nrhs; j++)
698  {
699  X(fsupc, j) /= m_mapL.valuePtr()[luptr];
700  }
701  }
702  else
703  {
704  Map<const Matrix<Scalar,Dynamic,Dynamic>, 0, OuterStride<> > A( &(m_mapL.valuePtr()[luptr]), nsupc, nsupc, OuterStride<>(lda) );
705  Map< Matrix<Scalar,Dynamic,Dynamic>, 0, OuterStride<> > U (&(X(fsupc,0)), nsupc, nrhs, OuterStride<>(n) );
706  U = A.template triangularView<Upper>().solve(U);
707  }
708 
709  for (Index j = 0; j < nrhs; ++j)
710  {
711  for (Index jcol = fsupc; jcol < fsupc + nsupc; jcol++)
712  {
713  typename MatrixUType::InnerIterator it(m_mapU, jcol);
714  for ( ; it; ++it)
715  {
716  Index irow = it.index();
717  X(irow, j) -= X(jcol, j) * it.value();
718  }
719  }
720  }
721  } // End For U-solve
722  }
723  const MatrixLType& m_mapL;
724  const MatrixUType& m_mapU;
725 };
726 
727 namespace internal {
728 
729 template<typename _MatrixType, typename Derived, typename Rhs>
730 struct solve_retval<SparseLU<_MatrixType,Derived>, Rhs>
731  : solve_retval_base<SparseLU<_MatrixType,Derived>, Rhs>
732 {
733  typedef SparseLU<_MatrixType,Derived> Dec;
734  EIGEN_MAKE_SOLVE_HELPERS(Dec,Rhs)
735 
736  template<typename Dest> void evalTo(Dest& dst) const
737  {
738  dec()._solve(rhs(),dst);
739  }
740 };
741 
742 template<typename _MatrixType, typename Derived, typename Rhs>
743 struct sparse_solve_retval<SparseLU<_MatrixType,Derived>, Rhs>
744  : sparse_solve_retval_base<SparseLU<_MatrixType,Derived>, Rhs>
745 {
746  typedef SparseLU<_MatrixType,Derived> Dec;
747  EIGEN_MAKE_SPARSE_SOLVE_HELPERS(Dec,Rhs)
748 
749  template<typename Dest> void evalTo(Dest& dst) const
750  {
751  this->defaultEvalTo(dst);
752  }
753 };
754 } // end namespace internal
755 
756 } // End namespace Eigen
757 
758 #endif
SparseLUMatrixLReturnType< SCMatrix > matrixL() const
Definition: SparseLU.h:134
Index rows() const
Definition: SparseMatrix.h:119
void analyzePattern(const MatrixType &matrix)
Definition: SparseLU.h:371
Index cols() const
Definition: SparseMatrix.h:121
const IndicesType & indices() const
Definition: PermutationMatrix.h:358
Transpose< PermutationBase > inverse() const
Definition: PermutationMatrix.h:201
Definition: Constants.h:378
const internal::sparse_solve_retval< SparseLU, Rhs > solve(const SparseMatrixBase< Rhs > &B) const
Definition: SparseLU.h:191
Scalar absDeterminant()
Definition: SparseLU.h:258
void factorize(const MatrixType &matrix)
Definition: SparseLU.h:454
Scalar logAbsDeterminant() const
Definition: SparseLU.h:288
ColXpr col(Index i)
Definition: DenseBase.h:733
const PermutationType & rowsPermutation() const
Definition: SparseLU.h:153
Sparse supernodal LU factorization for general matrices.
Definition: SparseLU.h:17
const PermutationType & colsPermutation() const
Definition: SparseLU.h:161
int coletree(const MatrixType &mat, IndexVector &parent, IndexVector &firstRowElt, typename MatrixType::Index *perm=0)
Definition: SparseColEtree.h:61
Expression of a fixed-size or dynamic-size sub-vector.
Definition: ForwardDeclarations.h:83
Base class of any sparse matrices or sparse expressions.
Definition: SparseMatrixBase.h:26
Derived & derived()
Definition: EigenBase.h:34
const internal::solve_retval< SparseLU, Rhs > solve(const MatrixBase< Rhs > &B) const
Definition: SparseLU.h:178
void isSymmetric(bool sym)
Definition: SparseLU.h:123
void compute(const MatrixType &matrix)
Definition: SparseLU.h:112
Derived & setConstant(Index size, const Scalar &value)
Definition: CwiseNullaryOp.h:348
SparseLUMatrixUReturnType< SCMatrix, MappedSparseMatrix< Scalar, ColMajor, Index > > matrixU() const
Definition: SparseLU.h:144
ComputationInfo info() const
Reports whether previous computation was successful.
Definition: SparseLU.h:207
Definition: Constants.h:376
const unsigned int RowMajorBit
Definition: Constants.h:53
void resize(Index nbRows, Index nbCols)
Definition: PlainObjectBase.h:235
std::string lastErrorMessage() const
Definition: SparseLU.h:216
Index rows() const
Definition: SparseMatrixBase.h:150
void setPivotThreshold(const RealScalar &thresh)
Definition: SparseLU.h:166
ComputationInfo
Definition: Constants.h:374
Base class for all dense matrices, vectors, and expressions.
Definition: MatrixBase.h:48
Scalar signDeterminant()
Definition: SparseLU.h:311
void treePostorder(Index n, IndexVector &parent, IndexVector &post)
Post order a tree.
Definition: SparseColEtree.h:178
Derived & setZero(Index size)
Definition: CwiseNullaryOp.h:515