10#ifndef EIGEN_SPARSESOLVERBASE_H
11#define EIGEN_SPARSESOLVERBASE_H
21template<
typename Decomposition,
typename Rhs,
typename Dest>
22typename enable_if<Rhs::ColsAtCompileTime!=1 && Dest::ColsAtCompileTime!=1>::type
23solve_sparse_through_dense_panels(
const Decomposition &dec,
const Rhs& rhs, Dest &dest)
25 EIGEN_STATIC_ASSERT((Dest::Flags&
RowMajorBit)==0,THIS_METHOD_IS_ONLY_FOR_COLUMN_MAJOR_MATRICES);
26 typedef typename Dest::Scalar DestScalar;
28 static const Index NbColsAtOnce = 4;
29 Index rhsCols = rhs.cols();
30 Index size = rhs.rows();
32 Index tmpCols = (std::min)(rhsCols, NbColsAtOnce);
35 for(
Index k=0; k<rhsCols; k+=NbColsAtOnce)
37 Index actualCols = std::min<Index>(rhsCols-k, NbColsAtOnce);
38 tmp.leftCols(actualCols) = rhs.middleCols(k,actualCols);
39 tmpX.leftCols(actualCols) = dec.solve(tmp.leftCols(actualCols));
40 dest.middleCols(k,actualCols) = tmpX.leftCols(actualCols).sparseView();
45template<
typename Decomposition,
typename Rhs,
typename Dest>
46typename enable_if<Rhs::ColsAtCompileTime==1 || Dest::ColsAtCompileTime==1>::type
47solve_sparse_through_dense_panels(
const Decomposition &dec,
const Rhs& rhs, Dest &dest)
49 typedef typename Dest::Scalar DestScalar;
50 Index size = rhs.rows();
53 dest_dense = dec.solve(rhs_dense);
54 dest = dest_dense.sparseView();
66template<
typename Derived>
73 : m_isInitialized(false)
79 Derived& derived() {
return *
static_cast<Derived*
>(
this); }
80 const Derived& derived()
const {
return *
static_cast<const Derived*
>(
this); }
86 template<
typename Rhs>
87 inline const Solve<Derived, Rhs>
90 eigen_assert(m_isInitialized &&
"Solver is not initialized.");
91 eigen_assert(derived().rows()==b.
rows() &&
"solve(): invalid number of rows of the right hand side matrix b");
99 template<
typename Rhs>
103 eigen_assert(m_isInitialized &&
"Solver is not initialized.");
104 eigen_assert(derived().rows()==b.
rows() &&
"solve(): invalid number of rows of the right hand side matrix b");
108 #ifndef EIGEN_PARSED_BY_DOXYGEN
110 template<
typename Rhs,
typename Dest>
113 internal::solve_sparse_through_dense_panels(derived(), b.
derived(), dest.
derived());
119 mutable bool m_isInitialized;
Derived & derived()
Definition: EigenBase.h:46
EIGEN_CONSTEXPR Index rows() const EIGEN_NOEXCEPT
Definition: EigenBase.h:60
Base class for all dense matrices, vectors, and expressions.
Definition: MatrixBase.h:50
The matrix class, also used for vectors and row-vectors.
Definition: Matrix.h:180
Pseudo expression representing a solving operation.
Definition: Solve.h:63
Base class of any sparse matrices or sparse expressions.
Definition: SparseMatrixBase.h:28
Index rows() const
Definition: SparseMatrixBase.h:176
A base class for sparse solvers.
Definition: SparseSolverBase.h:68
const Solve< Derived, Rhs > solve(const MatrixBase< Rhs > &b) const
Definition: SparseSolverBase.h:88
const Solve< Derived, Rhs > solve(const SparseMatrixBase< Rhs > &b) const
Definition: SparseSolverBase.h:101
SparseSolverBase()
Definition: SparseSolverBase.h:72
const unsigned int RowMajorBit
Definition: Constants.h:66
Namespace containing all symbols from the Eigen library.
Definition: Core:141
EIGEN_DEFAULT_DENSE_INDEX_TYPE Index
The Index type as used for the API.
Definition: Meta.h:74
Derived & derived()
Definition: EigenBase.h:46