00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025 #ifndef EIGEN_SPARSEDENSEPRODUCT_H
00026 #define EIGEN_SPARSEDENSEPRODUCT_H
00027
00028 template<typename Lhs, typename Rhs, int InnerSize> struct SparseDenseProductReturnType
00029 {
00030 typedef SparseTimeDenseProduct<Lhs,Rhs> Type;
00031 };
00032
00033 template<typename Lhs, typename Rhs> struct SparseDenseProductReturnType<Lhs,Rhs,1>
00034 {
00035 typedef SparseDenseOuterProduct<Lhs,Rhs,false> Type;
00036 };
00037
00038 template<typename Lhs, typename Rhs, int InnerSize> struct DenseSparseProductReturnType
00039 {
00040 typedef DenseTimeSparseProduct<Lhs,Rhs> Type;
00041 };
00042
00043 template<typename Lhs, typename Rhs> struct DenseSparseProductReturnType<Lhs,Rhs,1>
00044 {
00045 typedef SparseDenseOuterProduct<Rhs,Lhs,true> Type;
00046 };
00047
00048 namespace internal {
00049
00050 template<typename Lhs, typename Rhs, bool Tr>
00051 struct traits<SparseDenseOuterProduct<Lhs,Rhs,Tr> >
00052 {
00053 typedef Sparse StorageKind;
00054 typedef typename scalar_product_traits<typename traits<Lhs>::Scalar,
00055 typename traits<Rhs>::Scalar>::ReturnType Scalar;
00056 typedef typename Lhs::Index Index;
00057 typedef typename Lhs::Nested LhsNested;
00058 typedef typename Rhs::Nested RhsNested;
00059 typedef typename remove_all<LhsNested>::type _LhsNested;
00060 typedef typename remove_all<RhsNested>::type _RhsNested;
00061
00062 enum {
00063 LhsCoeffReadCost = traits<_LhsNested>::CoeffReadCost,
00064 RhsCoeffReadCost = traits<_RhsNested>::CoeffReadCost,
00065
00066 RowsAtCompileTime = Tr ? int(traits<Rhs>::RowsAtCompileTime) : int(traits<Lhs>::RowsAtCompileTime),
00067 ColsAtCompileTime = Tr ? int(traits<Lhs>::ColsAtCompileTime) : int(traits<Rhs>::ColsAtCompileTime),
00068 MaxRowsAtCompileTime = Tr ? int(traits<Rhs>::MaxRowsAtCompileTime) : int(traits<Lhs>::MaxRowsAtCompileTime),
00069 MaxColsAtCompileTime = Tr ? int(traits<Lhs>::MaxColsAtCompileTime) : int(traits<Rhs>::MaxColsAtCompileTime),
00070
00071 Flags = Tr ? RowMajorBit : 0,
00072
00073 CoeffReadCost = LhsCoeffReadCost + RhsCoeffReadCost + NumTraits<Scalar>::MulCost
00074 };
00075 };
00076
00077 }
00078
00079 template<typename Lhs, typename Rhs, bool Tr>
00080 class SparseDenseOuterProduct
00081 : public SparseMatrixBase<SparseDenseOuterProduct<Lhs,Rhs,Tr> >
00082 {
00083 public:
00084
00085 typedef SparseMatrixBase<SparseDenseOuterProduct> Base;
00086 EIGEN_DENSE_PUBLIC_INTERFACE(SparseDenseOuterProduct)
00087 typedef internal::traits<SparseDenseOuterProduct> Traits;
00088
00089 private:
00090
00091 typedef typename Traits::LhsNested LhsNested;
00092 typedef typename Traits::RhsNested RhsNested;
00093 typedef typename Traits::_LhsNested _LhsNested;
00094 typedef typename Traits::_RhsNested _RhsNested;
00095
00096 public:
00097
00098 class InnerIterator;
00099
00100 EIGEN_STRONG_INLINE SparseDenseOuterProduct(const Lhs& lhs, const Rhs& rhs)
00101 : m_lhs(lhs), m_rhs(rhs)
00102 {
00103 EIGEN_STATIC_ASSERT(!Tr,YOU_MADE_A_PROGRAMMING_MISTAKE);
00104 }
00105
00106 EIGEN_STRONG_INLINE SparseDenseOuterProduct(const Rhs& rhs, const Lhs& lhs)
00107 : m_lhs(lhs), m_rhs(rhs)
00108 {
00109 EIGEN_STATIC_ASSERT(Tr,YOU_MADE_A_PROGRAMMING_MISTAKE);
00110 }
00111
00112 EIGEN_STRONG_INLINE Index rows() const { return Tr ? m_rhs.rows() : m_lhs.rows(); }
00113 EIGEN_STRONG_INLINE Index cols() const { return Tr ? m_lhs.cols() : m_rhs.cols(); }
00114
00115 EIGEN_STRONG_INLINE const _LhsNested& lhs() const { return m_lhs; }
00116 EIGEN_STRONG_INLINE const _RhsNested& rhs() const { return m_rhs; }
00117
00118 protected:
00119 LhsNested m_lhs;
00120 RhsNested m_rhs;
00121 };
00122
00123 template<typename Lhs, typename Rhs, bool Transpose>
00124 class SparseDenseOuterProduct<Lhs,Rhs,Transpose>::InnerIterator : public _LhsNested::InnerIterator
00125 {
00126 typedef typename _LhsNested::InnerIterator Base;
00127 public:
00128 EIGEN_STRONG_INLINE InnerIterator(const SparseDenseOuterProduct& prod, Index outer)
00129 : Base(prod.lhs(), 0), m_outer(outer), m_factor(prod.rhs().coeff(outer))
00130 {
00131 }
00132
00133 inline Index outer() const { return m_outer; }
00134 inline Index row() const { return Transpose ? Base::row() : m_outer; }
00135 inline Index col() const { return Transpose ? m_outer : Base::row(); }
00136
00137 inline Scalar value() const { return Base::value() * m_factor; }
00138
00139 protected:
00140 int m_outer;
00141 Scalar m_factor;
00142 };
00143
00144 namespace internal {
00145 template<typename Lhs, typename Rhs>
00146 struct traits<SparseTimeDenseProduct<Lhs,Rhs> >
00147 : traits<ProductBase<SparseTimeDenseProduct<Lhs,Rhs>, Lhs, Rhs> >
00148 {
00149 typedef Dense StorageKind;
00150 typedef MatrixXpr XprKind;
00151 };
00152 }
00153
00154 template<typename Lhs, typename Rhs>
00155 class SparseTimeDenseProduct
00156 : public ProductBase<SparseTimeDenseProduct<Lhs,Rhs>, Lhs, Rhs>
00157 {
00158 public:
00159 EIGEN_PRODUCT_PUBLIC_INTERFACE(SparseTimeDenseProduct)
00160
00161 SparseTimeDenseProduct(const Lhs& lhs, const Rhs& rhs) : Base(lhs,rhs)
00162 {}
00163
00164 template<typename Dest> void scaleAndAddTo(Dest& dest, Scalar alpha) const
00165 {
00166 typedef typename internal::remove_all<Lhs>::type _Lhs;
00167 typedef typename internal::remove_all<Rhs>::type _Rhs;
00168 typedef typename _Lhs::InnerIterator LhsInnerIterator;
00169 enum { LhsIsRowMajor = (_Lhs::Flags&RowMajorBit)==RowMajorBit };
00170 for(Index j=0; j<m_lhs.outerSize(); ++j)
00171 {
00172 typename Rhs::Scalar rhs_j = alpha * m_rhs.coeff(LhsIsRowMajor ? 0 : j,0);
00173 typename Dest::RowXpr dest_j(dest.row(LhsIsRowMajor ? j : 0));
00174 for(LhsInnerIterator it(m_lhs,j); it ;++it)
00175 {
00176 if(LhsIsRowMajor) dest_j += (alpha*it.value()) * m_rhs.row(it.index());
00177 else if(Rhs::ColsAtCompileTime==1) dest.coeffRef(it.index()) += it.value() * rhs_j;
00178 else dest.row(it.index()) += (alpha*it.value()) * m_rhs.row(j);
00179 }
00180 }
00181 }
00182
00183 private:
00184 SparseTimeDenseProduct& operator=(const SparseTimeDenseProduct&);
00185 };
00186
00187
00188
00189 namespace internal {
00190 template<typename Lhs, typename Rhs>
00191 struct traits<DenseTimeSparseProduct<Lhs,Rhs> >
00192 : traits<ProductBase<DenseTimeSparseProduct<Lhs,Rhs>, Lhs, Rhs> >
00193 {
00194 typedef Dense StorageKind;
00195 };
00196 }
00197
00198 template<typename Lhs, typename Rhs>
00199 class DenseTimeSparseProduct
00200 : public ProductBase<DenseTimeSparseProduct<Lhs,Rhs>, Lhs, Rhs>
00201 {
00202 public:
00203 EIGEN_PRODUCT_PUBLIC_INTERFACE(DenseTimeSparseProduct)
00204
00205 DenseTimeSparseProduct(const Lhs& lhs, const Rhs& rhs) : Base(lhs,rhs)
00206 {}
00207
00208 template<typename Dest> void scaleAndAddTo(Dest& dest, Scalar alpha) const
00209 {
00210 typedef typename internal::remove_all<Rhs>::type _Rhs;
00211 typedef typename _Rhs::InnerIterator RhsInnerIterator;
00212 enum { RhsIsRowMajor = (_Rhs::Flags&RowMajorBit)==RowMajorBit };
00213 for(Index j=0; j<m_rhs.outerSize(); ++j)
00214 for(RhsInnerIterator i(m_rhs,j); i; ++i)
00215 dest.col(RhsIsRowMajor ? i.index() : j) += (alpha*i.value()) * m_lhs.col(RhsIsRowMajor ? j : i.index());
00216 }
00217
00218 private:
00219 DenseTimeSparseProduct& operator=(const DenseTimeSparseProduct&);
00220 };
00221
00222
00223 template<typename Derived>
00224 template<typename OtherDerived>
00225 inline const typename SparseDenseProductReturnType<Derived,OtherDerived>::Type
00226 SparseMatrixBase<Derived>::operator*(const MatrixBase<OtherDerived> &other) const
00227 {
00228 return typename SparseDenseProductReturnType<Derived,OtherDerived>::Type(derived(), other.derived());
00229 }
00230
00231 #endif // EIGEN_SPARSEDENSEPRODUCT_H