UmfPackSupport.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) 2008-2011 Gael Guennebaud <gael.guennebaud@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 #ifndef EIGEN_UMFPACKSUPPORT_H
11 #define EIGEN_UMFPACKSUPPORT_H
12 
13 // for compatibility with super old version of umfpack,
14 // not sure this is really needed, but this is harmless.
15 #ifndef SuiteSparse_long
16 #ifdef UF_long
17 #define SuiteSparse_long UF_long
18 #else
19 #error neither SuiteSparse_long nor UF_long are defined
20 #endif
21 #endif
22 
23 namespace Eigen {
24 
25 /* TODO extract L, extract U, compute det, etc... */
26 
27 // generic double/complex<double> wrapper functions:
28 
29 
30  // Defaults
31 inline void umfpack_defaults(double control[UMFPACK_CONTROL], double, int)
32 { umfpack_di_defaults(control); }
33 
34 inline void umfpack_defaults(double control[UMFPACK_CONTROL], std::complex<double>, int)
35 { umfpack_zi_defaults(control); }
36 
37 inline void umfpack_defaults(double control[UMFPACK_CONTROL], double, SuiteSparse_long)
38 { umfpack_dl_defaults(control); }
39 
40 inline void umfpack_defaults(double control[UMFPACK_CONTROL], std::complex<double>, SuiteSparse_long)
41 { umfpack_zl_defaults(control); }
42 
43 // Report info
44 inline void umfpack_report_info(double control[UMFPACK_CONTROL], double info[UMFPACK_INFO], double, int)
45 { umfpack_di_report_info(control, info);}
46 
47 inline void umfpack_report_info(double control[UMFPACK_CONTROL], double info[UMFPACK_INFO], std::complex<double>, int)
48 { umfpack_zi_report_info(control, info);}
49 
50 inline void umfpack_report_info(double control[UMFPACK_CONTROL], double info[UMFPACK_INFO], double, SuiteSparse_long)
51 { umfpack_dl_report_info(control, info);}
52 
53 inline void umfpack_report_info(double control[UMFPACK_CONTROL], double info[UMFPACK_INFO], std::complex<double>, SuiteSparse_long)
54 { umfpack_zl_report_info(control, info);}
55 
56 // Report status
57 inline void umfpack_report_status(double control[UMFPACK_CONTROL], int status, double, int)
58 { umfpack_di_report_status(control, status);}
59 
60 inline void umfpack_report_status(double control[UMFPACK_CONTROL], int status, std::complex<double>, int)
61 { umfpack_zi_report_status(control, status);}
62 
63 inline void umfpack_report_status(double control[UMFPACK_CONTROL], int status, double, SuiteSparse_long)
64 { umfpack_dl_report_status(control, status);}
65 
66 inline void umfpack_report_status(double control[UMFPACK_CONTROL], int status, std::complex<double>, SuiteSparse_long)
67 { umfpack_zl_report_status(control, status);}
68 
69 // report control
70 inline void umfpack_report_control(double control[UMFPACK_CONTROL], double, int)
71 { umfpack_di_report_control(control);}
72 
73 inline void umfpack_report_control(double control[UMFPACK_CONTROL], std::complex<double>, int)
74 { umfpack_zi_report_control(control);}
75 
76 inline void umfpack_report_control(double control[UMFPACK_CONTROL], double, SuiteSparse_long)
77 { umfpack_dl_report_control(control);}
78 
79 inline void umfpack_report_control(double control[UMFPACK_CONTROL], std::complex<double>, SuiteSparse_long)
80 { umfpack_zl_report_control(control);}
81 
82 // Free numeric
83 inline void umfpack_free_numeric(void **Numeric, double, int)
84 { umfpack_di_free_numeric(Numeric); *Numeric = 0; }
85 
86 inline void umfpack_free_numeric(void **Numeric, std::complex<double>, int)
87 { umfpack_zi_free_numeric(Numeric); *Numeric = 0; }
88 
89 inline void umfpack_free_numeric(void **Numeric, double, SuiteSparse_long)
90 { umfpack_dl_free_numeric(Numeric); *Numeric = 0; }
91 
92 inline void umfpack_free_numeric(void **Numeric, std::complex<double>, SuiteSparse_long)
93 { umfpack_zl_free_numeric(Numeric); *Numeric = 0; }
94 
95 // Free symbolic
96 inline void umfpack_free_symbolic(void **Symbolic, double, int)
97 { umfpack_di_free_symbolic(Symbolic); *Symbolic = 0; }
98 
99 inline void umfpack_free_symbolic(void **Symbolic, std::complex<double>, int)
100 { umfpack_zi_free_symbolic(Symbolic); *Symbolic = 0; }
101 
102 inline void umfpack_free_symbolic(void **Symbolic, double, SuiteSparse_long)
103 { umfpack_dl_free_symbolic(Symbolic); *Symbolic = 0; }
104 
105 inline void umfpack_free_symbolic(void **Symbolic, std::complex<double>, SuiteSparse_long)
106 { umfpack_zl_free_symbolic(Symbolic); *Symbolic = 0; }
107 
108 // Symbolic
109 inline int umfpack_symbolic(int n_row,int n_col,
110  const int Ap[], const int Ai[], const double Ax[], void **Symbolic,
111  const double Control [UMFPACK_CONTROL], double Info [UMFPACK_INFO])
112 {
113  return umfpack_di_symbolic(n_row,n_col,Ap,Ai,Ax,Symbolic,Control,Info);
114 }
115 
116 inline int umfpack_symbolic(int n_row,int n_col,
117  const int Ap[], const int Ai[], const std::complex<double> Ax[], void **Symbolic,
118  const double Control [UMFPACK_CONTROL], double Info [UMFPACK_INFO])
119 {
120  return umfpack_zi_symbolic(n_row,n_col,Ap,Ai,&numext::real_ref(Ax[0]),0,Symbolic,Control,Info);
121 }
123  const SuiteSparse_long Ap[], const SuiteSparse_long Ai[], const double Ax[], void **Symbolic,
124  const double Control [UMFPACK_CONTROL], double Info [UMFPACK_INFO])
125 {
126  return umfpack_dl_symbolic(n_row,n_col,Ap,Ai,Ax,Symbolic,Control,Info);
127 }
128 
130  const SuiteSparse_long Ap[], const SuiteSparse_long Ai[], const std::complex<double> Ax[], void **Symbolic,
131  const double Control [UMFPACK_CONTROL], double Info [UMFPACK_INFO])
132 {
133  return umfpack_zl_symbolic(n_row,n_col,Ap,Ai,&numext::real_ref(Ax[0]),0,Symbolic,Control,Info);
134 }
135 
136 // Numeric
137 inline int umfpack_numeric( const int Ap[], const int Ai[], const double Ax[],
138  void *Symbolic, void **Numeric,
139  const double Control[UMFPACK_CONTROL],double Info [UMFPACK_INFO])
140 {
141  return umfpack_di_numeric(Ap,Ai,Ax,Symbolic,Numeric,Control,Info);
142 }
143 
144 inline int umfpack_numeric( const int Ap[], const int Ai[], const std::complex<double> Ax[],
145  void *Symbolic, void **Numeric,
146  const double Control[UMFPACK_CONTROL],double Info [UMFPACK_INFO])
147 {
148  return umfpack_zi_numeric(Ap,Ai,&numext::real_ref(Ax[0]),0,Symbolic,Numeric,Control,Info);
149 }
150 inline SuiteSparse_long umfpack_numeric(const SuiteSparse_long Ap[], const SuiteSparse_long Ai[], const double Ax[],
151  void *Symbolic, void **Numeric,
152  const double Control[UMFPACK_CONTROL],double Info [UMFPACK_INFO])
153 {
154  return umfpack_dl_numeric(Ap,Ai,Ax,Symbolic,Numeric,Control,Info);
155 }
156 
157 inline SuiteSparse_long umfpack_numeric(const SuiteSparse_long Ap[], const SuiteSparse_long Ai[], const std::complex<double> Ax[],
158  void *Symbolic, void **Numeric,
159  const double Control[UMFPACK_CONTROL],double Info [UMFPACK_INFO])
160 {
161  return umfpack_zl_numeric(Ap,Ai,&numext::real_ref(Ax[0]),0,Symbolic,Numeric,Control,Info);
162 }
163 
164 // solve
165 inline int umfpack_solve( int sys, const int Ap[], const int Ai[], const double Ax[],
166  double X[], const double B[], void *Numeric,
167  const double Control[UMFPACK_CONTROL], double Info[UMFPACK_INFO])
168 {
169  return umfpack_di_solve(sys,Ap,Ai,Ax,X,B,Numeric,Control,Info);
170 }
171 
172 inline int umfpack_solve( int sys, const int Ap[], const int Ai[], const std::complex<double> Ax[],
173  std::complex<double> X[], const std::complex<double> B[], void *Numeric,
174  const double Control[UMFPACK_CONTROL], double Info[UMFPACK_INFO])
175 {
176  return umfpack_zi_solve(sys,Ap,Ai,&numext::real_ref(Ax[0]),0,&numext::real_ref(X[0]),0,&numext::real_ref(B[0]),0,Numeric,Control,Info);
177 }
178 
179 inline SuiteSparse_long umfpack_solve(int sys, const SuiteSparse_long Ap[], const SuiteSparse_long Ai[], const double Ax[],
180  double X[], const double B[], void *Numeric,
181  const double Control[UMFPACK_CONTROL], double Info[UMFPACK_INFO])
182 {
183  return umfpack_dl_solve(sys,Ap,Ai,Ax,X,B,Numeric,Control,Info);
184 }
185 
186 inline SuiteSparse_long umfpack_solve(int sys, const SuiteSparse_long Ap[], const SuiteSparse_long Ai[], const std::complex<double> Ax[],
187  std::complex<double> X[], const std::complex<double> B[], void *Numeric,
188  const double Control[UMFPACK_CONTROL], double Info[UMFPACK_INFO])
189 {
190  return umfpack_zl_solve(sys,Ap,Ai,&numext::real_ref(Ax[0]),0,&numext::real_ref(X[0]),0,&numext::real_ref(B[0]),0,Numeric,Control,Info);
191 }
192 
193 // Get Lunz
194 inline int umfpack_get_lunz(int *lnz, int *unz, int *n_row, int *n_col, int *nz_udiag, void *Numeric, double)
195 {
196  return umfpack_di_get_lunz(lnz,unz,n_row,n_col,nz_udiag,Numeric);
197 }
198 
199 inline int umfpack_get_lunz(int *lnz, int *unz, int *n_row, int *n_col, int *nz_udiag, void *Numeric, std::complex<double>)
200 {
201  return umfpack_zi_get_lunz(lnz,unz,n_row,n_col,nz_udiag,Numeric);
202 }
203 
205  SuiteSparse_long *nz_udiag, void *Numeric, double)
206 {
207  return umfpack_dl_get_lunz(lnz,unz,n_row,n_col,nz_udiag,Numeric);
208 }
209 
211  SuiteSparse_long *nz_udiag, void *Numeric, std::complex<double>)
212 {
213  return umfpack_zl_get_lunz(lnz,unz,n_row,n_col,nz_udiag,Numeric);
214 }
215 
216 // Get Numeric
217 inline int umfpack_get_numeric(int Lp[], int Lj[], double Lx[], int Up[], int Ui[], double Ux[],
218  int P[], int Q[], double Dx[], int *do_recip, double Rs[], void *Numeric)
219 {
220  return umfpack_di_get_numeric(Lp,Lj,Lx,Up,Ui,Ux,P,Q,Dx,do_recip,Rs,Numeric);
221 }
222 
223 inline int umfpack_get_numeric(int Lp[], int Lj[], std::complex<double> Lx[], int Up[], int Ui[], std::complex<double> Ux[],
224  int P[], int Q[], std::complex<double> Dx[], int *do_recip, double Rs[], void *Numeric)
225 {
226  double& lx0_real = numext::real_ref(Lx[0]);
227  double& ux0_real = numext::real_ref(Ux[0]);
228  double& dx0_real = numext::real_ref(Dx[0]);
229  return umfpack_zi_get_numeric(Lp,Lj,Lx?&lx0_real:0,0,Up,Ui,Ux?&ux0_real:0,0,P,Q,
230  Dx?&dx0_real:0,0,do_recip,Rs,Numeric);
231 }
233  SuiteSparse_long P[], SuiteSparse_long Q[], double Dx[], SuiteSparse_long *do_recip, double Rs[], void *Numeric)
234 {
235  return umfpack_dl_get_numeric(Lp,Lj,Lx,Up,Ui,Ux,P,Q,Dx,do_recip,Rs,Numeric);
236 }
237 
238 inline SuiteSparse_long umfpack_get_numeric(SuiteSparse_long Lp[], SuiteSparse_long Lj[], std::complex<double> Lx[], SuiteSparse_long Up[], SuiteSparse_long Ui[], std::complex<double> Ux[],
239  SuiteSparse_long P[], SuiteSparse_long Q[], std::complex<double> Dx[], SuiteSparse_long *do_recip, double Rs[], void *Numeric)
240 {
241  double& lx0_real = numext::real_ref(Lx[0]);
242  double& ux0_real = numext::real_ref(Ux[0]);
243  double& dx0_real = numext::real_ref(Dx[0]);
244  return umfpack_zl_get_numeric(Lp,Lj,Lx?&lx0_real:0,0,Up,Ui,Ux?&ux0_real:0,0,P,Q,
245  Dx?&dx0_real:0,0,do_recip,Rs,Numeric);
246 }
247 
248 // Get Determinant
249 inline int umfpack_get_determinant(double *Mx, double *Ex, void *NumericHandle, double User_Info [UMFPACK_INFO], int)
250 {
251  return umfpack_di_get_determinant(Mx,Ex,NumericHandle,User_Info);
252 }
253 
254 inline int umfpack_get_determinant(std::complex<double> *Mx, double *Ex, void *NumericHandle, double User_Info [UMFPACK_INFO], int)
255 {
256  double& mx_real = numext::real_ref(*Mx);
257  return umfpack_zi_get_determinant(&mx_real,0,Ex,NumericHandle,User_Info);
258 }
259 
260 inline SuiteSparse_long umfpack_get_determinant(double *Mx, double *Ex, void *NumericHandle, double User_Info [UMFPACK_INFO], SuiteSparse_long)
261 {
262  return umfpack_dl_get_determinant(Mx,Ex,NumericHandle,User_Info);
263 }
264 
265 inline SuiteSparse_long umfpack_get_determinant(std::complex<double> *Mx, double *Ex, void *NumericHandle, double User_Info [UMFPACK_INFO], SuiteSparse_long)
266 {
267  double& mx_real = numext::real_ref(*Mx);
268  return umfpack_zl_get_determinant(&mx_real,0,Ex,NumericHandle,User_Info);
269 }
270 
271 
287 template<typename _MatrixType>
288 class UmfPackLU : public SparseSolverBase<UmfPackLU<_MatrixType> >
289 {
290  protected:
292  using Base::m_isInitialized;
293  public:
294  using Base::_solve_impl;
295  typedef _MatrixType MatrixType;
296  typedef typename MatrixType::Scalar Scalar;
298  typedef typename MatrixType::StorageIndex StorageIndex;
305  enum {
306  ColsAtCompileTime = MatrixType::ColsAtCompileTime,
307  MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime
308  };
309 
310  public:
311 
314 
316  : m_dummy(0,0), mp_matrix(m_dummy)
317  {
318  init();
319  }
320 
321  template<typename InputMatrixType>
322  explicit UmfPackLU(const InputMatrixType& matrix)
323  : mp_matrix(matrix)
324  {
325  init();
326  compute(matrix);
327  }
328 
330  {
333  }
334 
335  inline Index rows() const { return mp_matrix.rows(); }
336  inline Index cols() const { return mp_matrix.cols(); }
337 
344  {
345  eigen_assert(m_isInitialized && "Decomposition is not initialized.");
346  return m_info;
347  }
348 
349  inline const LUMatrixType& matrixL() const
350  {
352  return m_l;
353  }
354 
355  inline const LUMatrixType& matrixU() const
356  {
358  return m_u;
359  }
360 
361  inline const IntColVectorType& permutationP() const
362  {
364  return m_p;
365  }
366 
367  inline const IntRowVectorType& permutationQ() const
368  {
370  return m_q;
371  }
372 
377  template<typename InputMatrixType>
378  void compute(const InputMatrixType& matrix)
379  {
382  grab(matrix.derived());
384  factorize_impl();
385  }
386 
393  template<typename InputMatrixType>
394  void analyzePattern(const InputMatrixType& matrix)
395  {
398 
399  grab(matrix.derived());
400 
402  }
403 
409  inline int umfpackFactorizeReturncode() const
410  {
411  eigen_assert(m_numeric && "UmfPackLU: you must first call factorize()");
412  return m_fact_errorCode;
413  }
414 
421  inline const UmfpackControl& umfpackControl() const
422  {
423  return m_control;
424  }
425 
432  inline UmfpackControl& umfpackControl()
433  {
434  return m_control;
435  }
436 
443  template<typename InputMatrixType>
444  void factorize(const InputMatrixType& matrix)
445  {
446  eigen_assert(m_analysisIsOk && "UmfPackLU: you must first call analyzePattern()");
447  if(m_numeric)
449 
450  grab(matrix.derived());
451 
452  factorize_impl();
453  }
454 
460  {
462  }
463 
469  {
470  eigen_assert(m_analysisIsOk && "UmfPackLU: you must first call analyzePattern()");
472  }
473 
479  eigen_assert(m_analysisIsOk && "UmfPackLU: you must first call analyzePattern()");
481  }
482 
484  template<typename BDerived,typename XDerived>
486 
487  Scalar determinant() const;
488 
489  void extractData() const;
490 
491  protected:
492 
493  void init()
494  {
496  m_isInitialized = false;
497  m_numeric = 0;
498  m_symbolic = 0;
500 
502  }
503 
505  {
506  m_fact_errorCode = umfpack_symbolic(internal::convert_index<StorageIndex>(mp_matrix.rows()),
507  internal::convert_index<StorageIndex>(mp_matrix.cols()),
508  mp_matrix.outerIndexPtr(), mp_matrix.innerIndexPtr(), mp_matrix.valuePtr(),
510 
511  m_isInitialized = true;
513  m_analysisIsOk = true;
514  m_factorizationIsOk = false;
516  }
517 
519  {
520 
521  m_fact_errorCode = umfpack_numeric(mp_matrix.outerIndexPtr(), mp_matrix.innerIndexPtr(), mp_matrix.valuePtr(),
523 
524  m_info = m_fact_errorCode == UMFPACK_OK ? Success : NumericalIssue;
525  m_factorizationIsOk = true;
527  }
528 
529  template<typename MatrixDerived>
531  {
532  mp_matrix.~UmfpackMatrixRef();
533  ::new (&mp_matrix) UmfpackMatrixRef(A.derived());
534  }
535 
536  void grab(const UmfpackMatrixRef &A)
537  {
538  if(&(A.derived()) != &mp_matrix)
539  {
540  mp_matrix.~UmfpackMatrixRef();
541  ::new (&mp_matrix) UmfpackMatrixRef(A);
542  }
543  }
544 
545  // cached data to reduce reallocation, etc.
546  mutable LUMatrixType m_l;
547  StorageIndex m_fact_errorCode;
548  UmfpackControl m_control;
549  mutable UmfpackInfo m_umfpackInfo;
550 
551  mutable LUMatrixType m_u;
552  mutable IntColVectorType m_p;
553  mutable IntRowVectorType m_q;
554 
555  UmfpackMatrixType m_dummy;
556  UmfpackMatrixRef mp_matrix;
557 
558  void* m_numeric;
559  void* m_symbolic;
560 
565 
566  private:
567  UmfPackLU(const UmfPackLU& ) { }
568 };
569 
570 
571 template<typename MatrixType>
573 {
575  {
576  // get size of the data
577  StorageIndex lnz, unz, rows, cols, nz_udiag;
578  umfpack_get_lunz(&lnz, &unz, &rows, &cols, &nz_udiag, m_numeric, Scalar());
579 
580  // allocate data
581  m_l.resize(rows,(std::min)(rows,cols));
582  m_l.resizeNonZeros(lnz);
583 
584  m_u.resize((std::min)(rows,cols),cols);
585  m_u.resizeNonZeros(unz);
586 
587  m_p.resize(rows);
588  m_q.resize(cols);
589 
590  // extract
593  m_p.data(), m_q.data(), 0, 0, 0, m_numeric);
594 
595  m_extractedDataAreDirty = false;
596  }
597 }
598 
599 template<typename MatrixType>
601 {
602  Scalar det;
604  return det;
605 }
606 
607 template<typename MatrixType>
608 template<typename BDerived,typename XDerived>
610 {
611  Index rhsCols = b.cols();
612  eigen_assert((BDerived::Flags&RowMajorBit)==0 && "UmfPackLU backend does not support non col-major rhs yet");
613  eigen_assert((XDerived::Flags&RowMajorBit)==0 && "UmfPackLU backend does not support non col-major result yet");
614  eigen_assert(b.derived().data() != x.derived().data() && " Umfpack does not support inplace solve");
615 
616  Scalar* x_ptr = 0;
618  if(x.innerStride()!=1)
619  {
620  x_tmp.resize(x.rows());
621  x_ptr = x_tmp.data();
622  }
623  for (int j=0; j<rhsCols; ++j)
624  {
625  if(x.innerStride()==1)
626  x_ptr = &x.col(j).coeffRef(0);
627  StorageIndex errorCode = umfpack_solve(UMFPACK_A,
628  mp_matrix.outerIndexPtr(), mp_matrix.innerIndexPtr(), mp_matrix.valuePtr(),
629  x_ptr, &b.const_cast_derived().col(j).coeffRef(0),
631  if(x.innerStride()!=1)
632  x.col(j) = x_tmp;
633  if (errorCode!=0)
634  return false;
635  }
636 
637  return true;
638 }
639 
640 } // end namespace Eigen
641 
642 #endif // EIGEN_UMFPACKSUPPORT_H
void umfpack_defaults(double control[UMFPACK_CONTROL], double, int)
EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE ColXpr col(Index i)
Definition: DenseBase.h:1098
Matrix< int, MatrixType::RowsAtCompileTime, 1 > IntColVectorType
SCALAR Scalar
Definition: bench_gemm.cpp:46
Index rows() const
void umfpack_report_control(double control[UMFPACK_CONTROL], double, int)
int umfpack_get_lunz(int *lnz, int *unz, int *n_row, int *n_col, int *nz_udiag, void *Numeric, double)
void printUmfpackControl()
A sparse LU factorization and solver based on UmfPack.
EIGEN_DEVICE_FUNC internal::add_const_on_value_type< EIGEN_MATHFUNC_RETVAL(real_ref, Scalar) >::type real_ref(const Scalar &x)
void umfpack_report_info(double control[UMFPACK_CONTROL], double info[UMFPACK_INFO], double, int)
void factorize(const InputMatrixType &matrix)
Scalar * b
Definition: benchVecAdd.cpp:17
const StorageIndex * outerIndexPtr() const
Definition: SparseMatrix.h:168
const IntColVectorType & permutationP() const
EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const Scalar * data() const
#define min(a, b)
Definition: datatypes.h:19
void compute(const InputMatrixType &matrix)
void analyzePattern(const InputMatrixType &matrix)
IntRowVectorType m_q
void resize(Index rows, Index cols)
Definition: SparseMatrix.h:626
void umfpack_report_status(double control[UMFPACK_CONTROL], int status, double, int)
A base class for sparse solvers.
StorageIndex m_fact_errorCode
Matrix< Scalar, Dynamic, 1 > Vector
Namespace containing all symbols from the Eigen library.
Definition: jet.h:637
LUMatrixType m_u
ComputationInfo info() const
Reports whether previous computation was successful.
const unsigned int RowMajorBit
Definition: Constants.h:66
void resizeNonZeros(Index size)
Definition: SparseMatrix.h:649
Index cols() const
void umfpack_free_numeric(void **Numeric, double, int)
MatrixType::RealScalar RealScalar
UmfpackMatrixType m_dummy
SparseMatrix< Scalar > LUMatrixType
else if n * info
IntColVectorType m_p
void _solve_impl(const SparseMatrixBase< Rhs > &b, SparseMatrixBase< Dest > &dest) const
void grab(const UmfpackMatrixRef &A)
MatrixType::StorageIndex StorageIndex
void umfpack_free_symbolic(void **Symbolic, double, int)
int umfpack_get_determinant(double *Mx, double *Ex, void *NumericHandle, double User_Info [UMFPACK_INFO], int)
EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void resize(Index rows, Index cols)
ComputationInfo m_info
UmfpackControl m_control
_MatrixType MatrixType
UmfpackInfo m_umfpackInfo
#define SuiteSparse_long
const LUMatrixType & matrixU() const
int umfpack_numeric(const int Ap[], const int Ai[], const double Ax[], void *Symbolic, void **Numeric, const double Control[UMFPACK_CONTROL], double Info [UMFPACK_INFO])
EIGEN_DEFAULT_DENSE_INDEX_TYPE Index
The Index type as used for the API.
Definition: Meta.h:74
#define eigen_assert(x)
Definition: Macros.h:1037
UmfpackMatrixRef mp_matrix
const UmfpackControl & umfpackControl() const
NumTraits< Scalar >::Real RealScalar
Definition: bench_gemm.cpp:47
int umfpack_symbolic(int n_row, int n_col, const int Ap[], const int Ai[], const double Ax[], void **Symbolic, const double Control [UMFPACK_CONTROL], double Info [UMFPACK_INFO])
void grab(const EigenBase< MatrixDerived > &A)
const StorageIndex * innerIndexPtr() const
Definition: SparseMatrix.h:159
Ref< const UmfpackMatrixType, StandardCompressedFormat > UmfpackMatrixRef
Array< double, UMFPACK_CONTROL, 1 > UmfpackControl
void analyzePattern_impl()
int umfpackFactorizeReturncode() const
const LUMatrixType & matrixL() const
Matrix< int, 1, MatrixType::ColsAtCompileTime > IntRowVectorType
The quaternion class used to represent 3D orientations and rotations.
LUMatrixType m_l
int umfpack_get_numeric(int Lp[], int Lj[], double Lx[], int Up[], int Ui[], double Ux[], int P[], int Q[], double Dx[], int *do_recip, double Rs[], void *Numeric)
MatrixType::Scalar Scalar
UmfpackControl & umfpackControl()
SparseMatrix< Scalar, ColMajor, StorageIndex > UmfpackMatrixType
const IntRowVectorType & permutationQ() const
bool _solve_impl(const MatrixBase< BDerived > &b, MatrixBase< XDerived > &x) const
void extractData() const
int umfpack_solve(int sys, const int Ap[], const int Ai[], const double Ax[], double X[], const double B[], void *Numeric, const double Control[UMFPACK_CONTROL], double Info[UMFPACK_INFO])
Map< Matrix< T, Dynamic, Dynamic, ColMajor >, 0, OuterStride<> > matrix(T *data, int rows, int cols, int stride)
UmfPackLU(const UmfPackLU &)
#define X
Definition: icosphere.cpp:20
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
Array< double, UMFPACK_INFO, 1 > UmfpackInfo
ComputationInfo
Definition: Constants.h:440
UmfPackLU(const InputMatrixType &matrix)
EIGEN_DEVICE_FUNC Derived & derived()
Definition: EigenBase.h:46
Scalar determinant() const
SparseSolverBase< UmfPackLU< _MatrixType > > Base
Base class for all dense matrices, vectors, and expressions.
Definition: MatrixBase.h:48
const Scalar * valuePtr() const
Definition: SparseMatrix.h:150
std::ptrdiff_t j


gtsam
Author(s):
autogenerated on Tue Jul 4 2023 02:40:41