diff options
Diffstat (limited to 'Eigen/src')
-rw-r--r-- | Eigen/src/Sparse/SparseLU.h | 20 | ||||
-rw-r--r-- | Eigen/src/Sparse/SuperLUSupport.h | 16 |
2 files changed, 31 insertions, 5 deletions
diff --git a/Eigen/src/Sparse/SparseLU.h b/Eigen/src/Sparse/SparseLU.h index 142592050..9cab67a51 100644 --- a/Eigen/src/Sparse/SparseLU.h +++ b/Eigen/src/Sparse/SparseLU.h @@ -25,6 +25,12 @@ #ifndef EIGEN_SPARSELU_H #define EIGEN_SPARSELU_H +enum { + SvNoTrans = 0, + SvTranspose = 1, + SvAdjoint = 2 +}; + /** \ingroup Sparse_Module * * \class SparseLU @@ -115,7 +121,8 @@ class SparseLU //inline const MatrixType& matrixU() const { return m_matrixU; } template<typename BDerived, typename XDerived> - bool solve(const MatrixBase<BDerived> &b, MatrixBase<XDerived>* x) const; + bool solve(const MatrixBase<BDerived> &b, MatrixBase<XDerived>* x, + const int transposed = SvNoTrans) const; /** \returns true if the factorization succeeded */ inline bool succeeded(void) const { return m_succeeded; } @@ -136,10 +143,17 @@ void SparseLU<MatrixType,Backend>::compute(const MatrixType& a) ei_assert(false && "not implemented yet"); } -/** Computes *x = U^-1 L^-1 b */ +/** Computes *x = U^-1 L^-1 b + * + * If \a transpose is set to SvTranspose or SvAdjoint, the solution + * of the transposed/adjoint system is computed instead. + * + * Not all backends implement the solution of the transposed or + * adjoint system. + */ template<typename MatrixType, int Backend> template<typename BDerived, typename XDerived> -bool SparseLU<MatrixType,Backend>::solve(const MatrixBase<BDerived> &b, MatrixBase<XDerived>* x) const +bool SparseLU<MatrixType,Backend>::solve(const MatrixBase<BDerived> &b, MatrixBase<XDerived>* x, const int transposed) const { ei_assert(false && "not implemented yet"); return false; diff --git a/Eigen/src/Sparse/SuperLUSupport.h b/Eigen/src/Sparse/SuperLUSupport.h index 6df94e35b..2aeccc534 100644 --- a/Eigen/src/Sparse/SuperLUSupport.h +++ b/Eigen/src/Sparse/SuperLUSupport.h @@ -318,7 +318,7 @@ class SparseLU<MatrixType,SuperLU> : public SparseLU<MatrixType> Scalar determinant() const; template<typename BDerived, typename XDerived> - bool solve(const MatrixBase<BDerived> &b, MatrixBase<XDerived>* x) const; + bool solve(const MatrixBase<BDerived> &b, MatrixBase<XDerived>* x, const int transposed = SvNoTrans) const; void compute(const MatrixType& matrix); @@ -413,12 +413,22 @@ void SparseLU<MatrixType,SuperLU>::compute(const MatrixType& a) template<typename MatrixType> template<typename BDerived,typename XDerived> -bool SparseLU<MatrixType,SuperLU>::solve(const MatrixBase<BDerived> &b, MatrixBase<XDerived> *x) const +bool SparseLU<MatrixType,SuperLU>::solve(const MatrixBase<BDerived> &b, + MatrixBase<XDerived> *x, const int transposed) const { const int size = m_matrix.rows(); const int rhsCols = b.cols(); ei_assert(size==b.rows()); + switch (transposed) { + case SvNoTrans : m_sluOptions.Trans = NOTRANS; break; + case SvTranspose : m_sluOptions.Trans = TRANS; break; + case SvAdjoint : m_sluOptions.Trans = CONJ; break; + default: + std::cerr << "Eigen: tranpsiotion option \"" << transposed << "\" not supported by the SuperLU backend\n"; + m_sluOptions.Trans = NOTRANS; + } + m_sluOptions.Fact = FACTORED; m_sluOptions.IterRefine = NOREFINE; @@ -443,6 +453,8 @@ bool SparseLU<MatrixType,SuperLU>::solve(const MatrixBase<BDerived> &b, MatrixBa &m_sluStat, &info, Scalar()); StatFree(&m_sluStat); + // reset to previous state + m_sluOptions.Trans = NOTRANS; return info==0; } |