10 #ifndef EIGEN_UMFPACKSUPPORT_H 11 #define EIGEN_UMFPACKSUPPORT_H 15 #ifndef SuiteSparse_long 17 #define SuiteSparse_long UF_long 19 #error neither SuiteSparse_long nor UF_long are defined 32 { umfpack_di_defaults(control); }
34 inline void umfpack_defaults(
double control[UMFPACK_CONTROL], std::complex<double>,
int)
35 { umfpack_zi_defaults(control); }
38 { umfpack_dl_defaults(control); }
41 { umfpack_zl_defaults(control); }
45 { umfpack_di_report_info(control,
info);}
48 { umfpack_zi_report_info(control,
info);}
51 { umfpack_dl_report_info(control,
info);}
54 { umfpack_zl_report_info(control,
info);}
58 { umfpack_di_report_status(control, status);}
61 { umfpack_zi_report_status(control, status);}
64 { umfpack_dl_report_status(control, status);}
67 { umfpack_zl_report_status(control, status);}
71 { umfpack_di_report_control(control);}
74 { umfpack_zi_report_control(control);}
77 { umfpack_dl_report_control(control);}
80 { umfpack_zl_report_control(control);}
84 { umfpack_di_free_numeric(Numeric); *Numeric = 0; }
87 { umfpack_zi_free_numeric(Numeric); *Numeric = 0; }
90 { umfpack_dl_free_numeric(Numeric); *Numeric = 0; }
93 { umfpack_zl_free_numeric(Numeric); *Numeric = 0; }
97 { umfpack_di_free_symbolic(Symbolic); *Symbolic = 0; }
100 { umfpack_zi_free_symbolic(Symbolic); *Symbolic = 0; }
103 { umfpack_dl_free_symbolic(Symbolic); *Symbolic = 0; }
106 { umfpack_zl_free_symbolic(Symbolic); *Symbolic = 0; }
110 const int Ap[],
const int Ai[],
const double Ax[],
void **Symbolic,
111 const double Control [UMFPACK_CONTROL],
double Info [UMFPACK_INFO])
113 return umfpack_di_symbolic(n_row,n_col,Ap,Ai,Ax,Symbolic,Control,Info);
117 const int Ap[],
const int Ai[],
const std::complex<double> Ax[],
void **Symbolic,
118 const double Control [UMFPACK_CONTROL],
double Info [UMFPACK_INFO])
120 return umfpack_zi_symbolic(n_row,n_col,Ap,Ai,&
numext::real_ref(Ax[0]),0,Symbolic,Control,Info);
124 const double Control [UMFPACK_CONTROL],
double Info [UMFPACK_INFO])
126 return umfpack_dl_symbolic(n_row,n_col,Ap,Ai,Ax,Symbolic,Control,Info);
131 const double Control [UMFPACK_CONTROL],
double Info [UMFPACK_INFO])
133 return umfpack_zl_symbolic(n_row,n_col,Ap,Ai,&
numext::real_ref(Ax[0]),0,Symbolic,Control,Info);
138 void *Symbolic,
void **Numeric,
139 const double Control[UMFPACK_CONTROL],
double Info [UMFPACK_INFO])
141 return umfpack_di_numeric(Ap,Ai,Ax,Symbolic,Numeric,Control,Info);
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])
148 return umfpack_zi_numeric(Ap,Ai,&
numext::real_ref(Ax[0]),0,Symbolic,Numeric,Control,Info);
151 void *Symbolic,
void **Numeric,
152 const double Control[UMFPACK_CONTROL],
double Info [UMFPACK_INFO])
154 return umfpack_dl_numeric(Ap,Ai,Ax,Symbolic,Numeric,Control,Info);
158 void *Symbolic,
void **Numeric,
159 const double Control[UMFPACK_CONTROL],
double Info [UMFPACK_INFO])
161 return umfpack_zl_numeric(Ap,Ai,&
numext::real_ref(Ax[0]),0,Symbolic,Numeric,Control,Info);
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])
169 return umfpack_di_solve(sys,Ap,Ai,Ax,X,B,Numeric,Control,Info);
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])
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);
180 double X[],
const double B[],
void *Numeric,
181 const double Control[UMFPACK_CONTROL],
double Info[UMFPACK_INFO])
183 return umfpack_dl_solve(sys,Ap,Ai,Ax,X,B,Numeric,Control,Info);
187 std::complex<double>
X[],
const std::complex<double>
B[],
void *Numeric,
188 const double Control[UMFPACK_CONTROL],
double Info[UMFPACK_INFO])
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);
194 inline int umfpack_get_lunz(
int *lnz,
int *unz,
int *n_row,
int *n_col,
int *nz_udiag,
void *Numeric,
double)
196 return umfpack_di_get_lunz(lnz,unz,n_row,n_col,nz_udiag,Numeric);
199 inline int umfpack_get_lunz(
int *lnz,
int *unz,
int *n_row,
int *n_col,
int *nz_udiag,
void *Numeric, std::complex<double>)
201 return umfpack_zi_get_lunz(lnz,unz,n_row,n_col,nz_udiag,Numeric);
207 return umfpack_dl_get_lunz(lnz,unz,n_row,n_col,nz_udiag,Numeric);
213 return umfpack_zl_get_lunz(lnz,unz,n_row,n_col,nz_udiag,Numeric);
218 int P[],
int Q[],
double Dx[],
int *do_recip,
double Rs[],
void *Numeric)
220 return umfpack_di_get_numeric(Lp,Lj,Lx,Up,Ui,Ux,P,Q,Dx,do_recip,Rs,Numeric);
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)
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);
235 return umfpack_dl_get_numeric(Lp,Lj,Lx,Up,Ui,Ux,P,Q,Dx,do_recip,Rs,Numeric);
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);
251 return umfpack_di_get_determinant(Mx,Ex,NumericHandle,User_Info);
254 inline int umfpack_get_determinant(std::complex<double> *Mx,
double *Ex,
void *NumericHandle,
double User_Info [UMFPACK_INFO],
int)
257 return umfpack_zi_get_determinant(&mx_real,0,Ex,NumericHandle,User_Info);
262 return umfpack_dl_get_determinant(Mx,Ex,NumericHandle,User_Info);
268 return umfpack_zl_get_determinant(&mx_real,0,Ex,NumericHandle,User_Info);
287 template<
typename _MatrixType>
321 template<
typename InputMatrixType>
377 template<
typename InputMatrixType>
382 grab(matrix.derived());
393 template<
typename InputMatrixType>
399 grab(matrix.derived());
443 template<
typename InputMatrixType>
450 grab(matrix.derived());
484 template<
typename BDerived,
typename XDerived>
507 internal::convert_index<StorageIndex>(
mp_matrix.cols()),
529 template<
typename MatrixDerived>
536 void grab(
const UmfpackMatrixRef &
A)
552 mutable IntColVectorType
m_p;
553 mutable IntRowVectorType
m_q;
571 template<
typename MatrixType>
599 template<
typename MatrixType>
607 template<
typename MatrixType>
608 template<
typename BDerived,
typename XDerived>
611 Index rhsCols = b.cols();
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");
618 if(x.innerStride()!=1)
621 x_ptr = x_tmp.
data();
623 for (
int j=0;
j<rhsCols; ++
j)
625 if(x.innerStride()==1)
626 x_ptr = &x.
col(
j).coeffRef(0);
629 x_ptr, &b.const_cast_derived().
col(
j).coeffRef(0),
631 if(x.innerStride()!=1)
642 #endif // EIGEN_UMFPACKSUPPORT_H void umfpack_defaults(double control[UMFPACK_CONTROL], double, int)
EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE ColXpr col(Index i)
Matrix< int, MatrixType::RowsAtCompileTime, 1 > IntColVectorType
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)
const StorageIndex * outerIndexPtr() const
const IntColVectorType & permutationP() const
EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const Scalar * data() const
void compute(const InputMatrixType &matrix)
void analyzePattern(const InputMatrixType &matrix)
void resize(Index rows, Index cols)
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.
void printUmfpackStatus()
ComputationInfo info() const
Reports whether previous computation was successful.
const unsigned int RowMajorBit
void resizeNonZeros(Index size)
void umfpack_free_numeric(void **Numeric, double, int)
MatrixType::RealScalar RealScalar
UmfpackMatrixType m_dummy
SparseMatrix< Scalar > LUMatrixType
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)
UmfpackInfo m_umfpackInfo
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.
UmfpackMatrixRef mp_matrix
const UmfpackControl & umfpackControl() const
NumTraits< Scalar >::Real RealScalar
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
Ref< const UmfpackMatrixType, StandardCompressedFormat > UmfpackMatrixRef
Array< double, UMFPACK_CONTROL, 1 > UmfpackControl
void analyzePattern_impl()
int umfpackFactorizeReturncode() const
bool m_extractedDataAreDirty
const LUMatrixType & matrixL() const
Matrix< int, 1, MatrixType::ColsAtCompileTime > IntRowVectorType
The quaternion class used to represent 3D orientations and rotations.
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
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 &)
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
UmfPackLU(const InputMatrixType &matrix)
EIGEN_DEVICE_FUNC Derived & derived()
Scalar determinant() const
SparseSolverBase< UmfPackLU< _MatrixType > > Base
Base class for all dense matrices, vectors, and expressions.
const Scalar * valuePtr() const