Eigen  3.4.90 (git rev a4098ac676528a83cfb73d4d26ce1b42ec05f47c)
SparseSolverBase.h
1// This file is part of Eigen, a lightweight C++ template library
2// for linear algebra.
3//
4// Copyright (C) 2014 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_SPARSESOLVERBASE_H
11#define EIGEN_SPARSESOLVERBASE_H
12
13#include "./InternalHeaderCheck.h"
14
15namespace Eigen {
16
17namespace internal {
18
23template<typename Decomposition, typename Rhs, typename Dest>
24typename enable_if<Rhs::ColsAtCompileTime!=1 && Dest::ColsAtCompileTime!=1>::type
25solve_sparse_through_dense_panels(const Decomposition &dec, const Rhs& rhs, Dest &dest)
26{
27 EIGEN_STATIC_ASSERT((Dest::Flags&RowMajorBit)==0,THIS_METHOD_IS_ONLY_FOR_COLUMN_MAJOR_MATRICES);
28 typedef typename Dest::Scalar DestScalar;
29 // we process the sparse rhs per block of NbColsAtOnce columns temporarily stored into a dense matrix.
30 static const Index NbColsAtOnce = 4;
31 Index rhsCols = rhs.cols();
32 Index size = rhs.rows();
33 // the temporary matrices do not need more columns than NbColsAtOnce:
34 Index tmpCols = (std::min)(rhsCols, NbColsAtOnce);
37 for(Index k=0; k<rhsCols; k+=NbColsAtOnce)
38 {
39 Index actualCols = std::min<Index>(rhsCols-k, NbColsAtOnce);
40 tmp.leftCols(actualCols) = rhs.middleCols(k,actualCols);
41 tmpX.leftCols(actualCols) = dec.solve(tmp.leftCols(actualCols));
42 dest.middleCols(k,actualCols) = tmpX.leftCols(actualCols).sparseView();
43 }
44}
45
46// Overload for vector as rhs
47template<typename Decomposition, typename Rhs, typename Dest>
48typename enable_if<Rhs::ColsAtCompileTime==1 || Dest::ColsAtCompileTime==1>::type
49solve_sparse_through_dense_panels(const Decomposition &dec, const Rhs& rhs, Dest &dest)
50{
51 typedef typename Dest::Scalar DestScalar;
52 Index size = rhs.rows();
55 dest_dense = dec.solve(rhs_dense);
56 dest = dest_dense.sparseView();
57}
58
59} // end namespace internal
60
68template<typename Derived>
69class SparseSolverBase : internal::noncopyable
70{
71 public:
72
75 : m_isInitialized(false)
76 {}
77
79 {}
80
81 Derived& derived() { return *static_cast<Derived*>(this); }
82 const Derived& derived() const { return *static_cast<const Derived*>(this); }
83
88 template<typename Rhs>
89 inline const Solve<Derived, Rhs>
90 solve(const MatrixBase<Rhs>& b) const
91 {
92 eigen_assert(m_isInitialized && "Solver is not initialized.");
93 eigen_assert(derived().rows()==b.rows() && "solve(): invalid number of rows of the right hand side matrix b");
94 return Solve<Derived, Rhs>(derived(), b.derived());
95 }
96
101 template<typename Rhs>
102 inline const Solve<Derived, Rhs>
104 {
105 eigen_assert(m_isInitialized && "Solver is not initialized.");
106 eigen_assert(derived().rows()==b.rows() && "solve(): invalid number of rows of the right hand side matrix b");
107 return Solve<Derived, Rhs>(derived(), b.derived());
108 }
109
110 #ifndef EIGEN_PARSED_BY_DOXYGEN
112 template<typename Rhs,typename Dest>
113 void _solve_impl(const SparseMatrixBase<Rhs> &b, SparseMatrixBase<Dest> &dest) const
114 {
115 internal::solve_sparse_through_dense_panels(derived(), b.derived(), dest.derived());
116 }
117 #endif // EIGEN_PARSED_BY_DOXYGEN
118
119 protected:
120
121 mutable bool m_isInitialized;
122};
123
124} // end namespace Eigen
125
126#endif // EIGEN_SPARSESOLVERBASE_H
Derived & derived()
Definition: EigenBase.h:48
EIGEN_CONSTEXPR Index rows() const EIGEN_NOEXCEPT
Definition: EigenBase.h:62
Base class for all dense matrices, vectors, and expressions.
Definition: MatrixBase.h:52
The matrix class, also used for vectors and row-vectors.
Definition: Matrix.h:182
Pseudo expression representing a solving operation.
Definition: Solve.h:65
Base class of any sparse matrices or sparse expressions.
Definition: SparseMatrixBase.h:30
Index rows() const
Definition: SparseMatrixBase.h:178
A base class for sparse solvers.
Definition: SparseSolverBase.h:70
const Solve< Derived, Rhs > solve(const MatrixBase< Rhs > &b) const
Definition: SparseSolverBase.h:90
const Solve< Derived, Rhs > solve(const SparseMatrixBase< Rhs > &b) const
Definition: SparseSolverBase.h:103
SparseSolverBase()
Definition: SparseSolverBase.h:74
const unsigned int RowMajorBit
Definition: Constants.h:68
Namespace containing all symbols from the Eigen library.
Definition: B01_Experimental.dox:1
EIGEN_DEFAULT_DENSE_INDEX_TYPE Index
The Index type as used for the API.
Definition: Meta.h:59
Derived & derived()
Definition: EigenBase.h:48