SparseQR.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-2013 Desire Nuentsa <desire.nuentsa_wakam@inria.fr>
5 // Copyright (C) 2012-2013 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 #ifndef EIGEN_SPARSE_QR_H
12 #define EIGEN_SPARSE_QR_H
13 
14 namespace Eigen {
15 
16 template<typename MatrixType, typename OrderingType> class SparseQR;
17 template<typename SparseQRType> struct SparseQRMatrixQReturnType;
18 template<typename SparseQRType> struct SparseQRMatrixQTransposeReturnType;
19 template<typename SparseQRType, typename Derived> struct SparseQR_QProduct;
20 namespace internal {
21  template <typename SparseQRType> struct traits<SparseQRMatrixQReturnType<SparseQRType> >
22  {
23  typedef typename SparseQRType::MatrixType ReturnType;
24  typedef typename ReturnType::Index Index;
25  typedef typename ReturnType::StorageKind StorageKind;
26  };
27  template <typename SparseQRType> struct traits<SparseQRMatrixQTransposeReturnType<SparseQRType> >
28  {
29  typedef typename SparseQRType::MatrixType ReturnType;
30  };
31  template <typename SparseQRType, typename Derived> struct traits<SparseQR_QProduct<SparseQRType, Derived> >
32  {
33  typedef typename Derived::PlainObject ReturnType;
34  };
35 } // End namespace internal
36 
63 template<typename _MatrixType, typename _OrderingType>
64 class SparseQR
65 {
66  public:
67  typedef _MatrixType MatrixType;
68  typedef _OrderingType OrderingType;
69  typedef typename MatrixType::Scalar Scalar;
70  typedef typename MatrixType::RealScalar RealScalar;
71  typedef typename MatrixType::Index Index;
76  public:
77  SparseQR () : m_isInitialized(false), m_analysisIsok(false), m_lastError(""), m_useDefaultThreshold(true),m_isQSorted(false)
78  { }
79 
80  SparseQR(const MatrixType& mat) : m_isInitialized(false), m_analysisIsok(false), m_lastError(""), m_useDefaultThreshold(true),m_isQSorted(false)
81  {
82  compute(mat);
83  }
84  void compute(const MatrixType& mat)
85  {
86  analyzePattern(mat);
87  factorize(mat);
88  }
89  void analyzePattern(const MatrixType& mat);
90  void factorize(const MatrixType& mat);
91 
94  inline Index rows() const { return m_pmat.rows(); }
95 
98  inline Index cols() const { return m_pmat.cols();}
99 
102  const QRMatrixType& matrixR() const { return m_R; }
103 
108  Index rank() const
109  {
110  eigen_assert(m_isInitialized && "The factorization should be called first, use compute()");
111  return m_nonzeropivots;
112  }
113 
133  { return SparseQRMatrixQReturnType<SparseQR>(*this); }
134 
138  const PermutationType& colsPermutation() const
139  {
140  eigen_assert(m_isInitialized && "Decomposition is not initialized.");
141  return m_outputPerm_c;
142  }
143 
147  std::string lastErrorMessage() const { return m_lastError; }
148 
150  template<typename Rhs, typename Dest>
151  bool _solve(const MatrixBase<Rhs> &B, MatrixBase<Dest> &dest) const
152  {
153  eigen_assert(m_isInitialized && "The factorization should be called first, use compute()");
154  eigen_assert(this->rows() == B.rows() && "SparseQR::solve() : invalid number of rows in the right hand side matrix");
155 
156  Index rank = this->rank();
157 
158  // Compute Q^T * b;
159  typename Dest::PlainObject y, b;
160  y = this->matrixQ().transpose() * B;
161  b = y;
162 
163  // Solve with the triangular matrix R
164  y.topRows(rank) = this->matrixR().topLeftCorner(rank, rank).template triangularView<Upper>().solve(b.topRows(rank));
165  y.bottomRows(y.size()-rank).setZero();
166 
167  // Apply the column permutation
168  if (m_perm_c.size()) dest.topRows(cols()) = colsPermutation() * y.topRows(cols());
169  else dest = y.topRows(cols());
170 
171  m_info = Success;
172  return true;
173  }
174 
175 
181  void setPivotThreshold(const RealScalar& threshold)
182  {
183  m_useDefaultThreshold = false;
184  m_threshold = threshold;
185  }
186 
191  template<typename Rhs>
193  {
194  eigen_assert(m_isInitialized && "The factorization should be called first, use compute()");
195  eigen_assert(this->rows() == B.rows() && "SparseQR::solve() : invalid number of rows in the right hand side matrix");
196  return internal::solve_retval<SparseQR, Rhs>(*this, B.derived());
197  }
198  template<typename Rhs>
200  {
201  eigen_assert(m_isInitialized && "The factorization should be called first, use compute()");
202  eigen_assert(this->rows() == B.rows() && "SparseQR::solve() : invalid number of rows in the right hand side matrix");
204  }
205 
215  {
216  eigen_assert(m_isInitialized && "Decomposition is not initialized.");
217  return m_info;
218  }
219 
220  protected:
221  inline void sort_matrix_Q()
222  {
223  if(this->m_isQSorted) return;
224  // The matrix Q is sorted during the transposition
226  this->m_Q = mQrm;
227  this->m_isQSorted = true;
228  }
229 
230 
231  protected:
236  std::string m_lastError;
237  QRMatrixType m_pmat; // Temporary matrix
238  QRMatrixType m_R; // The triangular factor matrix
239  QRMatrixType m_Q; // The orthogonal reflectors
240  ScalarVector m_hcoeffs; // The Householder coefficients
241  PermutationType m_perm_c; // Fill-reducing Column permutation
242  PermutationType m_pivotperm; // The permutation for rank revealing
243  PermutationType m_outputPerm_c; // The final column permutation
244  RealScalar m_threshold; // Threshold to determine null Householder reflections
245  bool m_useDefaultThreshold; // Use default threshold
246  Index m_nonzeropivots; // Number of non zero pivots found
247  IndexVector m_etree; // Column elimination tree
248  IndexVector m_firstRowElt; // First element in each row
249  bool m_isQSorted; // whether Q is sorted or not
250 
251  template <typename, typename > friend struct SparseQR_QProduct;
252  template <typename > friend struct SparseQRMatrixQReturnType;
253 
254 };
255 
263 template <typename MatrixType, typename OrderingType>
265 {
266  // Compute the column fill reducing ordering
267  OrderingType ord;
268  ord(mat, m_perm_c);
269  Index n = mat.cols();
270  Index m = mat.rows();
271 
272  if (!m_perm_c.size())
273  {
274  m_perm_c.resize(n);
275  m_perm_c.indices().setLinSpaced(n, 0,n-1);
276  }
277 
278  // Compute the column elimination tree of the permuted matrix
279  m_outputPerm_c = m_perm_c.inverse();
280  internal::coletree(mat, m_etree, m_firstRowElt, m_outputPerm_c.indices().data());
281 
282  m_R.resize(n, n);
283  m_Q.resize(m, n);
284 
285  // Allocate space for nonzero elements : rough estimation
286  m_R.reserve(2*mat.nonZeros()); //FIXME Get a more accurate estimation through symbolic factorization with the etree
287  m_Q.reserve(2*mat.nonZeros());
288  m_hcoeffs.resize(n);
289  m_analysisIsok = true;
290 }
291 
299 template <typename MatrixType, typename OrderingType>
301 {
302  using std::abs;
303  using std::max;
304 
305  eigen_assert(m_analysisIsok && "analyzePattern() should be called before this step");
306  Index m = mat.rows();
307  Index n = mat.cols();
308  IndexVector mark(m); mark.setConstant(-1); // Record the visited nodes
309  IndexVector Ridx(n), Qidx(m); // Store temporarily the row indexes for the current column of R and Q
310  Index nzcolR, nzcolQ; // Number of nonzero for the current column of R and Q
311  ScalarVector tval(m); // The dense vector used to compute the current column
312  bool found_diag;
313 
314  m_pmat = mat;
315  m_pmat.uncompress(); // To have the innerNonZeroPtr allocated
316  // Apply the fill-in reducing permutation lazily:
317  for (int i = 0; i < n; i++)
318  {
319  Index p = m_perm_c.size() ? m_perm_c.indices()(i) : i;
320  m_pmat.outerIndexPtr()[p] = mat.outerIndexPtr()[i];
321  m_pmat.innerNonZeroPtr()[p] = mat.outerIndexPtr()[i+1] - mat.outerIndexPtr()[i];
322  }
323 
324  /* Compute the default threshold, see :
325  * Tim Davis, "Algorithm 915, SuiteSparseQR: Multifrontal Multithreaded Rank-Revealing
326  * Sparse QR Factorization, ACM Trans. on Math. Soft. 38(1), 2011, Page 8:3
327  */
328  if(m_useDefaultThreshold)
329  {
330  RealScalar max2Norm = 0.0;
331  for (int j = 0; j < n; j++) max2Norm = (max)(max2Norm, m_pmat.col(j).norm());
332  m_threshold = 20 * (m + n) * max2Norm * NumTraits<RealScalar>::epsilon();
333  }
334 
335  // Initialize the numerical permutation
336  m_pivotperm.setIdentity(n);
337 
338  Index nonzeroCol = 0; // Record the number of valid pivots
339 
340  // Left looking rank-revealing QR factorization: compute a column of R and Q at a time
341  for (Index col = 0; col < n; ++col)
342  {
343  mark.setConstant(-1);
344  m_R.startVec(col);
345  m_Q.startVec(col);
346  mark(nonzeroCol) = col;
347  Qidx(0) = nonzeroCol;
348  nzcolR = 0; nzcolQ = 1;
349  found_diag = false;
350  tval.setZero();
351 
352  // Symbolic factorization: find the nonzero locations of the column k of the factors R and Q, i.e.,
353  // all the nodes (with indexes lower than rank) reachable through the column elimination tree (etree) rooted at node k.
354  // Note: if the diagonal entry does not exist, then its contribution must be explicitly added,
355  // thus the trick with found_diag that permits to do one more iteration on the diagonal element if this one has not been found.
356  for (typename MatrixType::InnerIterator itp(m_pmat, col); itp || !found_diag; ++itp)
357  {
358  Index curIdx = nonzeroCol ;
359  if(itp) curIdx = itp.row();
360  if(curIdx == nonzeroCol) found_diag = true;
361 
362  // Get the nonzeros indexes of the current column of R
363  Index st = m_firstRowElt(curIdx); // The traversal of the etree starts here
364  if (st < 0 )
365  {
366  m_lastError = "Empty row found during numerical factorization";
367  m_info = InvalidInput;
368  return;
369  }
370 
371  // Traverse the etree
372  Index bi = nzcolR;
373  for (; mark(st) != col; st = m_etree(st))
374  {
375  Ridx(nzcolR) = st; // Add this row to the list,
376  mark(st) = col; // and mark this row as visited
377  nzcolR++;
378  }
379 
380  // Reverse the list to get the topological ordering
381  Index nt = nzcolR-bi;
382  for(Index i = 0; i < nt/2; i++) std::swap(Ridx(bi+i), Ridx(nzcolR-i-1));
383 
384  // Copy the current (curIdx,pcol) value of the input matrix
385  if(itp) tval(curIdx) = itp.value();
386  else tval(curIdx) = Scalar(0);
387 
388  // Compute the pattern of Q(:,k)
389  if(curIdx > nonzeroCol && mark(curIdx) != col )
390  {
391  Qidx(nzcolQ) = curIdx; // Add this row to the pattern of Q,
392  mark(curIdx) = col; // and mark it as visited
393  nzcolQ++;
394  }
395  }
396 
397  // Browse all the indexes of R(:,col) in reverse order
398  for (Index i = nzcolR-1; i >= 0; i--)
399  {
400  Index curIdx = m_pivotperm.indices()(Ridx(i));
401 
402  // Apply the curIdx-th householder vector to the current column (temporarily stored into tval)
403  Scalar tdot(0);
404 
405  // First compute q' * tval
406  tdot = m_Q.col(curIdx).dot(tval);
407 
408  tdot *= m_hcoeffs(curIdx);
409 
410  // Then update tval = tval - q * tau
411  // FIXME: tval -= tdot * m_Q.col(curIdx) should amount to the same (need to check/add support for efficient "dense ?= sparse")
412  for (typename QRMatrixType::InnerIterator itq(m_Q, curIdx); itq; ++itq)
413  tval(itq.row()) -= itq.value() * tdot;
414 
415  // Detect fill-in for the current column of Q
416  if(m_etree(Ridx(i)) == nonzeroCol)
417  {
418  for (typename QRMatrixType::InnerIterator itq(m_Q, curIdx); itq; ++itq)
419  {
420  Index iQ = itq.row();
421  if (mark(iQ) != col)
422  {
423  Qidx(nzcolQ++) = iQ; // Add this row to the pattern of Q,
424  mark(iQ) = col; // and mark it as visited
425  }
426  }
427  }
428  } // End update current column
429 
430  // Compute the Householder reflection that eliminate the current column
431  // FIXME this step should call the Householder module.
432  Scalar tau;
433  RealScalar beta;
434  Scalar c0 = nzcolQ ? tval(Qidx(0)) : Scalar(0);
435 
436  // First, the squared norm of Q((col+1):m, col)
437  RealScalar sqrNorm = 0.;
438  for (Index itq = 1; itq < nzcolQ; ++itq) sqrNorm += numext::abs2(tval(Qidx(itq)));
439 
440  if(sqrNorm == RealScalar(0) && numext::imag(c0) == RealScalar(0))
441  {
442  tau = RealScalar(0);
443  beta = numext::real(c0);
444  tval(Qidx(0)) = 1;
445  }
446  else
447  {
448  beta = std::sqrt(numext::abs2(c0) + sqrNorm);
449  if(numext::real(c0) >= RealScalar(0))
450  beta = -beta;
451  tval(Qidx(0)) = 1;
452  for (Index itq = 1; itq < nzcolQ; ++itq)
453  tval(Qidx(itq)) /= (c0 - beta);
454  tau = numext::conj((beta-c0) / beta);
455 
456  }
457 
458  // Insert values in R
459  for (Index i = nzcolR-1; i >= 0; i--)
460  {
461  Index curIdx = Ridx(i);
462  if(curIdx < nonzeroCol)
463  {
464  m_R.insertBackByOuterInnerUnordered(col, curIdx) = tval(curIdx);
465  tval(curIdx) = Scalar(0.);
466  }
467  }
468 
469  if(abs(beta) >= m_threshold)
470  {
471  m_R.insertBackByOuterInner(col, nonzeroCol) = beta;
472  nonzeroCol++;
473  // The householder coefficient
474  m_hcoeffs(col) = tau;
475  // Record the householder reflections
476  for (Index itq = 0; itq < nzcolQ; ++itq)
477  {
478  Index iQ = Qidx(itq);
479  m_Q.insertBackByOuterInnerUnordered(col,iQ) = tval(iQ);
480  tval(iQ) = Scalar(0.);
481  }
482  }
483  else
484  {
485  // Zero pivot found: move implicitly this column to the end
486  m_hcoeffs(col) = Scalar(0);
487  for (Index j = nonzeroCol; j < n-1; j++)
488  std::swap(m_pivotperm.indices()(j), m_pivotperm.indices()[j+1]);
489 
490  // Recompute the column elimination tree
491  internal::coletree(m_pmat, m_etree, m_firstRowElt, m_pivotperm.indices().data());
492  }
493  }
494 
495  // Finalize the column pointers of the sparse matrices R and Q
496  m_Q.finalize();
497  m_Q.makeCompressed();
498  m_R.finalize();
499  m_R.makeCompressed();
500  m_isQSorted = false;
501 
502  m_nonzeropivots = nonzeroCol;
503 
504  if(nonzeroCol<n)
505  {
506  // Permute the triangular factor to put the 'dead' columns to the end
507  MatrixType tempR(m_R);
508  m_R = tempR * m_pivotperm;
509 
510  // Update the column permutation
511  m_outputPerm_c = m_outputPerm_c * m_pivotperm;
512  }
513 
514  m_isInitialized = true;
515  m_factorizationIsok = true;
516  m_info = Success;
517 }
518 
519 namespace internal {
520 
521 template<typename _MatrixType, typename OrderingType, typename Rhs>
522 struct solve_retval<SparseQR<_MatrixType,OrderingType>, Rhs>
523  : solve_retval_base<SparseQR<_MatrixType,OrderingType>, Rhs>
524 {
526  EIGEN_MAKE_SOLVE_HELPERS(Dec,Rhs)
527 
528  template<typename Dest> void evalTo(Dest& dst) const
529  {
530  dec()._solve(rhs(),dst);
531  }
532 };
533 template<typename _MatrixType, typename OrderingType, typename Rhs>
534 struct sparse_solve_retval<SparseQR<_MatrixType, OrderingType>, Rhs>
535  : sparse_solve_retval_base<SparseQR<_MatrixType, OrderingType>, Rhs>
536 {
539 
540  template<typename Dest> void evalTo(Dest& dst) const
541  {
542  this->defaultEvalTo(dst);
543  }
544 };
545 } // end namespace internal
546 
547 template <typename SparseQRType, typename Derived>
548 struct SparseQR_QProduct : ReturnByValue<SparseQR_QProduct<SparseQRType, Derived> >
549 {
550  typedef typename SparseQRType::QRMatrixType MatrixType;
551  typedef typename SparseQRType::Scalar Scalar;
552  typedef typename SparseQRType::Index Index;
553  // Get the references
554  SparseQR_QProduct(const SparseQRType& qr, const Derived& other, bool transpose) :
555  m_qr(qr),m_other(other),m_transpose(transpose) {}
556  inline Index rows() const { return m_transpose ? m_qr.rows() : m_qr.cols(); }
557  inline Index cols() const { return m_other.cols(); }
558 
559  // Assign to a vector
560  template<typename DesType>
561  void evalTo(DesType& res) const
562  {
563  Index n = m_qr.cols();
564  res = m_other;
565  if (m_transpose)
566  {
567  eigen_assert(m_qr.m_Q.rows() == m_other.rows() && "Non conforming object sizes");
568  //Compute res = Q' * other column by column
569  for(Index j = 0; j < res.cols(); j++){
570  for (Index k = 0; k < n; k++)
571  {
572  Scalar tau = Scalar(0);
573  tau = m_qr.m_Q.col(k).dot(res.col(j));
574  tau = tau * m_qr.m_hcoeffs(k);
575  res.col(j) -= tau * m_qr.m_Q.col(k);
576  }
577  }
578  }
579  else
580  {
581  eigen_assert(m_qr.m_Q.cols() == m_other.rows() && "Non conforming object sizes");
582  // Compute res = Q' * other column by column
583  for(Index j = 0; j < res.cols(); j++)
584  {
585  for (Index k = n-1; k >=0; k--)
586  {
587  Scalar tau = Scalar(0);
588  tau = m_qr.m_Q.col(k).dot(res.col(j));
589  tau = tau * m_qr.m_hcoeffs(k);
590  res.col(j) -= tau * m_qr.m_Q.col(k);
591  }
592  }
593  }
594  }
595 
596  const SparseQRType& m_qr;
597  const Derived& m_other;
599 };
600 
601 template<typename SparseQRType>
602 struct SparseQRMatrixQReturnType : public EigenBase<SparseQRMatrixQReturnType<SparseQRType> >
603 {
604  typedef typename SparseQRType::Index Index;
605  typedef typename SparseQRType::Scalar Scalar;
607  SparseQRMatrixQReturnType(const SparseQRType& qr) : m_qr(qr) {}
608  template<typename Derived>
610  {
611  return SparseQR_QProduct<SparseQRType,Derived>(m_qr,other.derived(),false);
612  }
614  {
616  }
617  inline Index rows() const { return m_qr.rows(); }
618  inline Index cols() const { return m_qr.cols(); }
619  // To use for operations with the transpose of Q
621  {
623  }
624  template<typename Dest> void evalTo(MatrixBase<Dest>& dest) const
625  {
626  dest.derived() = m_qr.matrixQ() * Dest::Identity(m_qr.rows(), m_qr.rows());
627  }
628  template<typename Dest> void evalTo(SparseMatrixBase<Dest>& dest) const
629  {
630  Dest idMat(m_qr.rows(), m_qr.rows());
631  idMat.setIdentity();
632  // Sort the sparse householder reflectors if needed
633  const_cast<SparseQRType *>(&m_qr)->sort_matrix_Q();
634  dest.derived() = SparseQR_QProduct<SparseQRType, Dest>(m_qr, idMat, false);
635  }
636 
637  const SparseQRType& m_qr;
638 };
639 
640 template<typename SparseQRType>
642 {
643  SparseQRMatrixQTransposeReturnType(const SparseQRType& qr) : m_qr(qr) {}
644  template<typename Derived>
646  {
647  return SparseQR_QProduct<SparseQRType,Derived>(m_qr,other.derived(), true);
648  }
649  const SparseQRType& m_qr;
650 };
651 
652 } // end namespace Eigen
653 
654 #endif
_OrderingType OrderingType
Definition: SparseQR.h:68
Matrix< Scalar, Dynamic, 1 > ScalarVector
Definition: SparseQR.h:74
const AutoDiffScalar< DerType > & conj(const AutoDiffScalar< DerType > &x)
IntermediateState sqrt(const Expression &arg)
QRMatrixType m_Q
Definition: SparseQR.h:239
_MatrixType MatrixType
Definition: SparseQR.h:67
ScalarVector m_hcoeffs
Definition: SparseQR.h:240
RowsBlockXpr topRows(Index n)
Definition: DenseBase.h:381
Index cols() const
Definition: SparseQR.h:557
Index rows() const
Definition: SparseQR.h:556
PermutationType m_perm_c
Definition: SparseQR.h:241
SparseQRMatrixQReturnType(const SparseQRType &qr)
Definition: SparseQR.h:607
const internal::sparse_solve_retval< SparseQR, Rhs > solve(const SparseMatrixBase< Rhs > &B) const
Definition: SparseQR.h:199
void evalTo(MatrixBase< Dest > &dest) const
Definition: SparseQR.h:624
const SparseQRType & m_qr
Definition: SparseQR.h:596
QRMatrixType m_pmat
Definition: SparseQR.h:237
Index rank() const
Definition: SparseQR.h:108
SparseQR_QProduct(const SparseQRType &qr, const Derived &other, bool transpose)
Definition: SparseQR.h:554
bool _solve(const MatrixBase< Rhs > &B, MatrixBase< Dest > &dest) const
Definition: SparseQR.h:151
iterative scaling algorithm to equilibrate rows and column norms in matrices
Definition: matrix.hpp:471
Holds information about the various numeric (i.e. scalar) types allowed by Eigen. ...
Definition: NumTraits.h:88
SparseMatrix< Scalar, ColMajor, Index > QRMatrixType
Definition: SparseQR.h:72
MatrixType::Index Index
Definition: SparseQR.h:71
SparseQRType::Scalar Scalar
Definition: SparseQR.h:605
EIGEN_STRONG_INLINE const CwiseUnaryOp< internal::scalar_abs2_op< Scalar >, const Derived > abs2() const
ComputationInfo info() const
Reports whether previous computation was successful.
Definition: SparseQR.h:214
const ImagReturnType imag() const
const QRMatrixType & matrixR() const
Definition: SparseQR.h:102
SparseQRType::QRMatrixType MatrixType
Definition: SparseQR.h:550
void factorize(const MatrixType &mat)
Performs the numerical QR factorization of the input matrix.
Definition: SparseQR.h:300
const Derived & m_other
Definition: SparseQR.h:597
bool m_isInitialized
Definition: SparseQR.h:232
SparseQR_QProduct< SparseQRType, Derived > operator*(const MatrixBase< Derived > &other)
Definition: SparseQR.h:645
int coletree(const MatrixType &mat, IndexVector &parent, IndexVector &firstRowElt, typename MatrixType::Index *perm=0)
EIGEN_STRONG_INLINE const CwiseUnaryOp< internal::scalar_abs_op< Scalar >, const Derived > abs() const
SparseQRType::Index Index
Definition: SparseQR.h:604
RealReturnType real() const
Sparse left-looking rank-revealing QR factorization.
Definition: SparseQR.h:16
IndexVector m_firstRowElt
Definition: SparseQR.h:248
Base class of any sparse matrices or sparse expressions.
PermutationMatrix< Dynamic, Dynamic, Index > PermutationType
Definition: SparseQR.h:75
Matrix< Scalar, Dynamic, Dynamic > DenseMatrix
Definition: SparseQR.h:606
Index cols() const
Definition: SparseQR.h:98
PermutationType m_pivotperm
Definition: SparseQR.h:242
SparseQRMatrixQTransposeReturnType< SparseQRType > adjoint() const
Definition: SparseQR.h:613
ComputationInfo m_info
Definition: SparseQR.h:235
SparseQRType::Index Index
Definition: SparseQR.h:552
QRMatrixType m_R
Definition: SparseQR.h:238
Derived & setZero(Index size)
#define EIGEN_MAKE_SPARSE_SOLVE_HELPERS(DecompositionType, Rhs)
Definition: SparseSolve.h:71
MatrixType::Scalar Scalar
Definition: SparseQR.h:69
Derived & setConstant(Index size, const Scalar &value)
void evalTo(DesType &res) const
Definition: SparseQR.h:561
bool m_isQSorted
Definition: SparseQR.h:249
void rhs(const real_t *x, real_t *f)
void evalTo(SparseMatrixBase< Dest > &dest) const
Definition: SparseQR.h:628
void analyzePattern(const MatrixType &mat)
Preprocessing step of a QR factorization.
Definition: SparseQR.h:264
SparseQRType::Scalar Scalar
Definition: SparseQR.h:551
PermutationType m_outputPerm_c
Definition: SparseQR.h:243
Index rows() const
Definition: SparseQR.h:94
const Derived & derived() const
void sort_matrix_Q()
Definition: SparseQR.h:221
const PermutationType & colsPermutation() const
Definition: SparseQR.h:138
const internal::solve_retval< SparseQR, Rhs > solve(const MatrixBase< Rhs > &B) const
Definition: SparseQR.h:192
SparseQRMatrixQReturnType< SparseQR > matrixQ() const
Definition: SparseQR.h:132
const SparseQRType & m_qr
Definition: SparseQR.h:637
SparseQRMatrixQTransposeReturnType(const SparseQRType &qr)
Definition: SparseQR.h:643
#define EIGEN_MAKE_SOLVE_HELPERS(DecompositionType, Rhs)
Definition: Solve.h:61
std::string m_lastError
Definition: SparseQR.h:236
Matrix< Index, Dynamic, 1 > IndexVector
Definition: SparseQR.h:73
ColXpr col(Index i)
Definition: BlockMethods.h:708
MatrixType::RealScalar RealScalar
Definition: SparseQR.h:70
#define eigen_assert(x)
bool m_analysisIsok
Definition: SparseQR.h:233
SparseQR(const MatrixType &mat)
Definition: SparseQR.h:80
Base class for all dense matrices, vectors, and expressions.
Definition: MatrixBase.h:48
IndexVector m_etree
Definition: SparseQR.h:247
RealScalar m_threshold
Definition: SparseQR.h:244
SparseQRMatrixQTransposeReturnType< SparseQRType > transpose() const
Definition: SparseQR.h:620
Index m_nonzeropivots
Definition: SparseQR.h:246
void setPivotThreshold(const RealScalar &threshold)
Definition: SparseQR.h:181
std::string lastErrorMessage() const
Definition: SparseQR.h:147
bool m_useDefaultThreshold
Definition: SparseQR.h:245
SparseQR_QProduct< SparseQRType, Derived > operator*(const MatrixBase< Derived > &other)
Definition: SparseQR.h:609
bool m_factorizationIsok
Definition: SparseQR.h:234
void compute(const MatrixType &mat)
Definition: SparseQR.h:84


acado
Author(s): Milan Vukov, Rien Quirynen
autogenerated on Mon Jun 10 2019 12:35:09