PaStiXSupport.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 #ifndef EIGEN_PASTIXSUPPORT_H
11 #define EIGEN_PASTIXSUPPORT_H
12 
13 namespace Eigen {
14 
23 template<typename _MatrixType, bool IsStrSym = false> class PastixLU;
24 template<typename _MatrixType, int Options> class PastixLLT;
25 template<typename _MatrixType, int Options> class PastixLDLT;
26 
27 namespace internal
28 {
29 
30  template<class Pastix> struct pastix_traits;
31 
32  template<typename _MatrixType>
33  struct pastix_traits< PastixLU<_MatrixType> >
34  {
35  typedef _MatrixType MatrixType;
36  typedef typename _MatrixType::Scalar Scalar;
37  typedef typename _MatrixType::RealScalar RealScalar;
38  typedef typename _MatrixType::Index Index;
39  };
40 
41  template<typename _MatrixType, int Options>
42  struct pastix_traits< PastixLLT<_MatrixType,Options> >
43  {
44  typedef _MatrixType MatrixType;
45  typedef typename _MatrixType::Scalar Scalar;
46  typedef typename _MatrixType::RealScalar RealScalar;
47  typedef typename _MatrixType::Index Index;
48  };
49 
50  template<typename _MatrixType, int Options>
51  struct pastix_traits< PastixLDLT<_MatrixType,Options> >
52  {
53  typedef _MatrixType MatrixType;
54  typedef typename _MatrixType::Scalar Scalar;
55  typedef typename _MatrixType::RealScalar RealScalar;
56  typedef typename _MatrixType::Index Index;
57  };
58 
59  void eigen_pastix(pastix_data_t **pastix_data, int pastix_comm, int n, int *ptr, int *idx, float *vals, int *perm, int * invp, float *x, int nbrhs, int *iparm, double *dparm)
60  {
61  if (n == 0) { ptr = NULL; idx = NULL; vals = NULL; }
62  if (nbrhs == 0) {x = NULL; nbrhs=1;}
63  s_pastix(pastix_data, pastix_comm, n, ptr, idx, vals, perm, invp, x, nbrhs, iparm, dparm);
64  }
65 
66  void eigen_pastix(pastix_data_t **pastix_data, int pastix_comm, int n, int *ptr, int *idx, double *vals, int *perm, int * invp, double *x, int nbrhs, int *iparm, double *dparm)
67  {
68  if (n == 0) { ptr = NULL; idx = NULL; vals = NULL; }
69  if (nbrhs == 0) {x = NULL; nbrhs=1;}
70  d_pastix(pastix_data, pastix_comm, n, ptr, idx, vals, perm, invp, x, nbrhs, iparm, dparm);
71  }
72 
73  void eigen_pastix(pastix_data_t **pastix_data, int pastix_comm, int n, int *ptr, int *idx, std::complex<float> *vals, int *perm, int * invp, std::complex<float> *x, int nbrhs, int *iparm, double *dparm)
74  {
75  if (n == 0) { ptr = NULL; idx = NULL; vals = NULL; }
76  if (nbrhs == 0) {x = NULL; nbrhs=1;}
77  c_pastix(pastix_data, pastix_comm, n, ptr, idx, reinterpret_cast<COMPLEX*>(vals), perm, invp, reinterpret_cast<COMPLEX*>(x), nbrhs, iparm, dparm);
78  }
79 
80  void eigen_pastix(pastix_data_t **pastix_data, int pastix_comm, int n, int *ptr, int *idx, std::complex<double> *vals, int *perm, int * invp, std::complex<double> *x, int nbrhs, int *iparm, double *dparm)
81  {
82  if (n == 0) { ptr = NULL; idx = NULL; vals = NULL; }
83  if (nbrhs == 0) {x = NULL; nbrhs=1;}
84  z_pastix(pastix_data, pastix_comm, n, ptr, idx, reinterpret_cast<DCOMPLEX*>(vals), perm, invp, reinterpret_cast<DCOMPLEX*>(x), nbrhs, iparm, dparm);
85  }
86 
87  // Convert the matrix to Fortran-style Numbering
88  template <typename MatrixType>
89  void c_to_fortran_numbering (MatrixType& mat)
90  {
91  if ( !(mat.outerIndexPtr()[0]) )
92  {
93  int i;
94  for(i = 0; i <= mat.rows(); ++i)
95  ++mat.outerIndexPtr()[i];
96  for(i = 0; i < mat.nonZeros(); ++i)
97  ++mat.innerIndexPtr()[i];
98  }
99  }
100 
101  // Convert to C-style Numbering
102  template <typename MatrixType>
103  void fortran_to_c_numbering (MatrixType& mat)
104  {
105  // Check the Numbering
106  if ( mat.outerIndexPtr()[0] == 1 )
107  { // Convert to C-style numbering
108  int i;
109  for(i = 0; i <= mat.rows(); ++i)
110  --mat.outerIndexPtr()[i];
111  for(i = 0; i < mat.nonZeros(); ++i)
112  --mat.innerIndexPtr()[i];
113  }
114  }
115 }
116 
117 // This is the base class to interface with PaStiX functions.
118 // Users should not used this class directly.
119 template <class Derived>
121 {
122  public:
124  typedef _MatrixType MatrixType;
125  typedef typename MatrixType::Scalar Scalar;
126  typedef typename MatrixType::RealScalar RealScalar;
127  typedef typename MatrixType::Index Index;
130 
131  public:
132 
133  PastixBase() : m_initisOk(false), m_analysisIsOk(false), m_factorizationIsOk(false), m_isInitialized(false), m_pastixdata(0), m_size(0)
134  {
135  init();
136  }
137 
139  {
140  clean();
141  }
142 
147  template<typename Rhs>
149  solve(const MatrixBase<Rhs>& b) const
150  {
151  eigen_assert(m_isInitialized && "Pastix solver is not initialized.");
152  eigen_assert(rows()==b.rows()
153  && "PastixBase::solve(): invalid number of rows of the right hand side matrix b");
154  return internal::solve_retval<PastixBase, Rhs>(*this, b.derived());
155  }
156 
157  template<typename Rhs,typename Dest>
158  bool _solve (const MatrixBase<Rhs> &b, MatrixBase<Dest> &x) const;
159 
160  Derived& derived()
161  {
162  return *static_cast<Derived*>(this);
163  }
164  const Derived& derived() const
165  {
166  return *static_cast<const Derived*>(this);
167  }
168 
175  {
176  return m_iparm;
177  }
178 
183  int& iparm(int idxparam)
184  {
185  return m_iparm(idxparam);
186  }
187 
193  {
194  return m_dparm;
195  }
196 
197 
201  double& dparm(int idxparam)
202  {
203  return m_dparm(idxparam);
204  }
205 
206  inline Index cols() const { return m_size; }
207  inline Index rows() const { return m_size; }
208 
218  {
219  eigen_assert(m_isInitialized && "Decomposition is not initialized.");
220  return m_info;
221  }
222 
227  template<typename Rhs>
229  solve(const SparseMatrixBase<Rhs>& b) const
230  {
231  eigen_assert(m_isInitialized && "Pastix LU, LLT or LDLT is not initialized.");
232  eigen_assert(rows()==b.rows()
233  && "PastixBase::solve(): invalid number of rows of the right hand side matrix b");
235  }
236 
237  protected:
238 
239  // Initialize the Pastix data structure, check the matrix
240  void init();
241 
242  // Compute the ordering and the symbolic factorization
243  void analyzePattern(ColSpMatrix& mat);
244 
245  // Compute the numerical factorization
246  void factorize(ColSpMatrix& mat);
247 
248  // Free all the data allocated by Pastix
249  void clean()
250  {
251  eigen_assert(m_initisOk && "The Pastix structure should be allocated first");
252  m_iparm(IPARM_START_TASK) = API_TASK_CLEAN;
253  m_iparm(IPARM_END_TASK) = API_TASK_CLEAN;
254  internal::eigen_pastix(&m_pastixdata, MPI_COMM_WORLD, 0, 0, 0, (Scalar*)0,
255  m_perm.data(), m_invp.data(), 0, 0, m_iparm.data(), m_dparm.data());
256  }
257 
258  void compute(ColSpMatrix& mat);
259 
265  mutable pastix_data_t *m_pastixdata; // Data structure for pastix
266  mutable int m_comm; // The MPI communicator identifier
267  mutable Matrix<int,IPARM_SIZE,1> m_iparm; // integer vector for the input parameters
268  mutable Matrix<double,DPARM_SIZE,1> m_dparm; // Scalar vector for the input parameters
269  mutable Matrix<Index,Dynamic,1> m_perm; // Permutation vector
270  mutable Matrix<Index,Dynamic,1> m_invp; // Inverse permutation vector
271  mutable int m_size; // Size of the matrix
272 };
273 
278 template <class Derived>
280 {
281  m_size = 0;
282  m_iparm.setZero(IPARM_SIZE);
283  m_dparm.setZero(DPARM_SIZE);
284 
285  m_iparm(IPARM_MODIFY_PARAMETER) = API_NO;
286  pastix(&m_pastixdata, MPI_COMM_WORLD,
287  0, 0, 0, 0,
288  0, 0, 0, 1, m_iparm.data(), m_dparm.data());
289 
290  m_iparm[IPARM_MATRIX_VERIFICATION] = API_NO;
291  m_iparm[IPARM_VERBOSE] = 2;
292  m_iparm[IPARM_ORDERING] = API_ORDER_SCOTCH;
293  m_iparm[IPARM_INCOMPLETE] = API_NO;
294  m_iparm[IPARM_OOC_LIMIT] = 2000;
295  m_iparm[IPARM_RHS_MAKING] = API_RHS_B;
296  m_iparm(IPARM_MATRIX_VERIFICATION) = API_NO;
297 
298  m_iparm(IPARM_START_TASK) = API_TASK_INIT;
299  m_iparm(IPARM_END_TASK) = API_TASK_INIT;
300  internal::eigen_pastix(&m_pastixdata, MPI_COMM_WORLD, 0, 0, 0, (Scalar*)0,
301  0, 0, 0, 0, m_iparm.data(), m_dparm.data());
302 
303  // Check the returned error
304  if(m_iparm(IPARM_ERROR_NUMBER)) {
305  m_info = InvalidInput;
306  m_initisOk = false;
307  }
308  else {
309  m_info = Success;
310  m_initisOk = true;
311  }
312 }
313 
314 template <class Derived>
316 {
317  eigen_assert(mat.rows() == mat.cols() && "The input matrix should be squared");
318 
319  analyzePattern(mat);
320  factorize(mat);
321 
322  m_iparm(IPARM_MATRIX_VERIFICATION) = API_NO;
323  m_isInitialized = m_factorizationIsOk;
324 }
325 
326 
327 template <class Derived>
329 {
330  eigen_assert(m_initisOk && "The initialization of PaSTiX failed");
331 
332  // clean previous calls
333  if(m_size>0)
334  clean();
335 
336  m_size = mat.rows();
337  m_perm.resize(m_size);
338  m_invp.resize(m_size);
339 
340  m_iparm(IPARM_START_TASK) = API_TASK_ORDERING;
341  m_iparm(IPARM_END_TASK) = API_TASK_ANALYSE;
342  internal::eigen_pastix(&m_pastixdata, MPI_COMM_WORLD, m_size, mat.outerIndexPtr(), mat.innerIndexPtr(),
343  mat.valuePtr(), m_perm.data(), m_invp.data(), 0, 0, m_iparm.data(), m_dparm.data());
344 
345  // Check the returned error
346  if(m_iparm(IPARM_ERROR_NUMBER))
347  {
348  m_info = NumericalIssue;
349  m_analysisIsOk = false;
350  }
351  else
352  {
353  m_info = Success;
354  m_analysisIsOk = true;
355  }
356 }
357 
358 template <class Derived>
360 {
361 // if(&m_cpyMat != &mat) m_cpyMat = mat;
362  eigen_assert(m_analysisIsOk && "The analysis phase should be called before the factorization phase");
363  m_iparm(IPARM_START_TASK) = API_TASK_NUMFACT;
364  m_iparm(IPARM_END_TASK) = API_TASK_NUMFACT;
365  m_size = mat.rows();
366 
367  internal::eigen_pastix(&m_pastixdata, MPI_COMM_WORLD, m_size, mat.outerIndexPtr(), mat.innerIndexPtr(),
368  mat.valuePtr(), m_perm.data(), m_invp.data(), 0, 0, m_iparm.data(), m_dparm.data());
369 
370  // Check the returned error
371  if(m_iparm(IPARM_ERROR_NUMBER))
372  {
373  m_info = NumericalIssue;
374  m_factorizationIsOk = false;
375  m_isInitialized = false;
376  }
377  else
378  {
379  m_info = Success;
380  m_factorizationIsOk = true;
381  m_isInitialized = true;
382  }
383 }
384 
385 /* Solve the system */
386 template<typename Base>
387 template<typename Rhs,typename Dest>
389 {
390  eigen_assert(m_isInitialized && "The matrix should be factorized first");
391  EIGEN_STATIC_ASSERT((Dest::Flags&RowMajorBit)==0,
392  THIS_METHOD_IS_ONLY_FOR_COLUMN_MAJOR_MATRICES);
393  int rhs = 1;
394 
395  x = b; /* on return, x is overwritten by the computed solution */
396 
397  for (int i = 0; i < b.cols(); i++){
398  m_iparm[IPARM_START_TASK] = API_TASK_SOLVE;
399  m_iparm[IPARM_END_TASK] = API_TASK_REFINE;
400 
401  internal::eigen_pastix(&m_pastixdata, MPI_COMM_WORLD, x.rows(), 0, 0, 0,
402  m_perm.data(), m_invp.data(), &x(0, i), rhs, m_iparm.data(), m_dparm.data());
403  }
404 
405  // Check the returned error
406  m_info = m_iparm(IPARM_ERROR_NUMBER)==0 ? Success : NumericalIssue;
407 
408  return m_iparm(IPARM_ERROR_NUMBER)==0;
409 }
410 
430 template<typename _MatrixType, bool IsStrSym>
431 class PastixLU : public PastixBase< PastixLU<_MatrixType> >
432 {
433  public:
436  typedef typename Base::ColSpMatrix ColSpMatrix;
437  typedef typename MatrixType::Index Index;
438 
439  public:
440  PastixLU() : Base()
441  {
442  init();
443  }
444 
445  PastixLU(const MatrixType& matrix):Base()
446  {
447  init();
448  compute(matrix);
449  }
455  void compute (const MatrixType& matrix)
456  {
457  m_structureIsUptodate = false;
458  ColSpMatrix temp;
459  grabMatrix(matrix, temp);
460  Base::compute(temp);
461  }
467  void analyzePattern(const MatrixType& matrix)
468  {
469  m_structureIsUptodate = false;
470  ColSpMatrix temp;
471  grabMatrix(matrix, temp);
472  Base::analyzePattern(temp);
473  }
474 
480  void factorize(const MatrixType& matrix)
481  {
482  ColSpMatrix temp;
483  grabMatrix(matrix, temp);
484  Base::factorize(temp);
485  }
486  protected:
487 
488  void init()
489  {
490  m_structureIsUptodate = false;
491  m_iparm(IPARM_SYM) = API_SYM_NO;
492  m_iparm(IPARM_FACTORIZATION) = API_FACT_LU;
493  }
494 
495  void grabMatrix(const MatrixType& matrix, ColSpMatrix& out)
496  {
497  if(IsStrSym)
498  out = matrix;
499  else
500  {
501  if(!m_structureIsUptodate)
502  {
503  // update the transposed structure
504  m_transposedStructure = matrix.transpose();
505 
506  // Set the elements of the matrix to zero
507  for (Index j=0; j<m_transposedStructure.outerSize(); ++j)
508  for(typename ColSpMatrix::InnerIterator it(m_transposedStructure, j); it; ++it)
509  it.valueRef() = 0.0;
510 
511  m_structureIsUptodate = true;
512  }
513 
514  out = m_transposedStructure + matrix;
515  }
517  }
518 
519  using Base::m_iparm;
520  using Base::m_dparm;
521 
524 };
525 
540 template<typename _MatrixType, int _UpLo>
541 class PastixLLT : public PastixBase< PastixLLT<_MatrixType, _UpLo> >
542 {
543  public:
546  typedef typename Base::ColSpMatrix ColSpMatrix;
547 
548  public:
549  enum { UpLo = _UpLo };
550  PastixLLT() : Base()
551  {
552  init();
553  }
554 
555  PastixLLT(const MatrixType& matrix):Base()
556  {
557  init();
558  compute(matrix);
559  }
560 
564  void compute (const MatrixType& matrix)
565  {
566  ColSpMatrix temp;
567  grabMatrix(matrix, temp);
568  Base::compute(temp);
569  }
570 
575  void analyzePattern(const MatrixType& matrix)
576  {
577  ColSpMatrix temp;
578  grabMatrix(matrix, temp);
579  Base::analyzePattern(temp);
580  }
584  void factorize(const MatrixType& matrix)
585  {
586  ColSpMatrix temp;
587  grabMatrix(matrix, temp);
588  Base::factorize(temp);
589  }
590  protected:
591  using Base::m_iparm;
592 
593  void init()
594  {
595  m_iparm(IPARM_SYM) = API_SYM_YES;
596  m_iparm(IPARM_FACTORIZATION) = API_FACT_LLT;
597  }
598 
599  void grabMatrix(const MatrixType& matrix, ColSpMatrix& out)
600  {
601  // Pastix supports only lower, column-major matrices
602  out.template selfadjointView<Lower>() = matrix.template selfadjointView<UpLo>();
604  }
605 };
606 
621 template<typename _MatrixType, int _UpLo>
622 class PastixLDLT : public PastixBase< PastixLDLT<_MatrixType, _UpLo> >
623 {
624  public:
627  typedef typename Base::ColSpMatrix ColSpMatrix;
628 
629  public:
630  enum { UpLo = _UpLo };
631  PastixLDLT():Base()
632  {
633  init();
634  }
635 
636  PastixLDLT(const MatrixType& matrix):Base()
637  {
638  init();
639  compute(matrix);
640  }
641 
645  void compute (const MatrixType& matrix)
646  {
647  ColSpMatrix temp;
648  grabMatrix(matrix, temp);
649  Base::compute(temp);
650  }
651 
656  void analyzePattern(const MatrixType& matrix)
657  {
658  ColSpMatrix temp;
659  grabMatrix(matrix, temp);
660  Base::analyzePattern(temp);
661  }
665  void factorize(const MatrixType& matrix)
666  {
667  ColSpMatrix temp;
668  grabMatrix(matrix, temp);
669  Base::factorize(temp);
670  }
671 
672  protected:
673  using Base::m_iparm;
674 
675  void init()
676  {
677  m_iparm(IPARM_SYM) = API_SYM_YES;
678  m_iparm(IPARM_FACTORIZATION) = API_FACT_LDLT;
679  }
680 
681  void grabMatrix(const MatrixType& matrix, ColSpMatrix& out)
682  {
683  // Pastix supports only lower, column-major matrices
684  out.template selfadjointView<Lower>() = matrix.template selfadjointView<UpLo>();
686  }
687 };
688 
689 namespace internal {
690 
691 template<typename _MatrixType, typename Rhs>
692 struct solve_retval<PastixBase<_MatrixType>, Rhs>
693  : solve_retval_base<PastixBase<_MatrixType>, Rhs>
694 {
696  EIGEN_MAKE_SOLVE_HELPERS(Dec,Rhs)
697 
698  template<typename Dest> void evalTo(Dest& dst) const
699  {
700  dec()._solve(rhs(),dst);
701  }
702 };
703 
704 template<typename _MatrixType, typename Rhs>
705 struct sparse_solve_retval<PastixBase<_MatrixType>, Rhs>
706  : sparse_solve_retval_base<PastixBase<_MatrixType>, Rhs>
707 {
710 
711  template<typename Dest> void evalTo(Dest& dst) const
712  {
713  this->defaultEvalTo(dst);
714  }
715 };
716 
717 } // end namespace internal
718 
719 } // end namespace Eigen
720 
721 #endif
SparseMatrix< Scalar, ColMajor > ColSpMatrix
void grabMatrix(const MatrixType &matrix, ColSpMatrix &out)
Matrix< Index, Dynamic, 1 > m_perm
ComputationInfo m_info
Index cols() const
Definition: SparseMatrix.h:121
void analyzePattern(ColSpMatrix &mat)
PastixLU(const MatrixType &matrix)
void compute(const MatrixType &matrix)
Base::ColSpMatrix ColSpMatrix
PastixLDLT(const MatrixType &matrix)
MatrixType::Index Index
Matrix< int, IPARM_SIZE, 1 > m_iparm
void factorize(const MatrixType &matrix)
iterative scaling algorithm to equilibrate rows and column norms in matrices
Definition: matrix.hpp:471
PastixLLT(const MatrixType &matrix)
void factorize(const MatrixType &matrix)
MatrixType::RealScalar RealScalar
#define EIGEN_STATIC_ASSERT(CONDITION, MSG)
Definition: StaticAssert.h:111
void compute(const MatrixType &matrix)
pastix_data_t * m_pastixdata
const unsigned int RowMajorBit
void compute(const MatrixType &matrix)
const Index * outerIndexPtr() const
Definition: SparseMatrix.h:149
Array< RealScalar, IPARM_SIZE, 1 > & dparm()
const Derived & derived() const
double & dparm(int idxparam)
void fortran_to_c_numbering(MatrixType &mat)
bool _solve(const MatrixBase< Rhs > &b, MatrixBase< Dest > &x) const
void compute(ColSpMatrix &mat)
MatrixType::Index Index
void analyzePattern(const MatrixType &matrix)
_MatrixType MatrixType
Base class of any sparse matrices or sparse expressions.
Matrix< Scalar, Dynamic, 1 > Vector
A sparse direct supernodal Cholesky (LLT) factorization and solver based on the PaStiX library...
Definition: PaStiXSupport.h:25
PastixBase< PastixLLT< MatrixType, _UpLo > > Base
_MatrixType MatrixType
bool m_structureIsUptodate
Provides a generic way to set and pass user-specified options.
Definition: options.hpp:65
PastixBase< PastixLU< MatrixType > > Base
int & iparm(int idxparam)
const internal::solve_retval< PastixBase, Rhs > solve(const MatrixBase< Rhs > &b) const
Transpose< Derived > transpose()
_MatrixType MatrixType
#define EIGEN_MAKE_SPARSE_SOLVE_HELPERS(DecompositionType, Rhs)
Definition: SparseSolve.h:71
void c_to_fortran_numbering(MatrixType &mat)
Definition: PaStiXSupport.h:89
void analyzePattern(const MatrixType &matrix)
void out(const real_t *x, real_t *f)
void rhs(const real_t *x, real_t *f)
void eigen_pastix(pastix_data_t **pastix_data, int pastix_comm, int n, int *ptr, int *idx, float *vals, int *perm, int *invp, float *x, int nbrhs, int *iparm, double *dparm)
Definition: PaStiXSupport.h:59
void grabMatrix(const MatrixType &matrix, ColSpMatrix &out)
Base::ColSpMatrix ColSpMatrix
Array< Index, IPARM_SIZE, 1 > & iparm()
const Derived & derived() const
Base::ColSpMatrix ColSpMatrix
void factorize(ColSpMatrix &mat)
void factorize(const MatrixType &matrix)
const internal::sparse_solve_retval< PastixBase, Rhs > solve(const SparseMatrixBase< Rhs > &b) const
General-purpose arrays with easy API for coefficient-wise operations.
Definition: Array.h:42
Index rows() const
const Index * innerIndexPtr() const
Definition: SparseMatrix.h:140
Interface to the PaStix solver.
Definition: PaStiXSupport.h:23
Matrix< Index, Dynamic, 1 > m_invp
#define EIGEN_MAKE_SOLVE_HELPERS(DecompositionType, Rhs)
Definition: Solve.h:61
A sparse direct supernodal Cholesky (LLT) factorization and solver based on the PaStiX library...
Definition: PaStiXSupport.h:24
Index cols() const
ComputationInfo info() const
Reports whether previous computation was successful.
ColSpMatrix m_transposedStructure
void grabMatrix(const MatrixType &matrix, ColSpMatrix &out)
#define eigen_assert(x)
internal::pastix_traits< Derived >::MatrixType _MatrixType
const Scalar * valuePtr() const
Definition: SparseMatrix.h:131
_MatrixType MatrixType
void analyzePattern(const MatrixType &matrix)
Matrix< double, DPARM_SIZE, 1 > m_dparm
Base class for all dense matrices, vectors, and expressions.
Definition: MatrixBase.h:48
void init(int nV, int nC, SymmetricMatrix *H, real_t *g, Matrix *A, const real_t *const lb, const real_t *const ub, const real_t *const lbA, const real_t *const ubA, int nWSR, const real_t *const x0, Options *options, int nOutputs, mxArray *plhs[])
Index rows() const
Definition: SparseMatrix.h:119
PastixBase< PastixLDLT< MatrixType, _UpLo > > Base
Derived & derived()
MatrixType::Scalar Scalar


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