10 #ifndef EIGEN_SPARSEDENSEPRODUCT_H 11 #define EIGEN_SPARSEDENSEPRODUCT_H 17 template <>
struct product_promote_storage_type<Sparse,Dense, OuterProduct> {
typedef Sparse ret; };
18 template <>
struct product_promote_storage_type<Dense,Sparse, OuterProduct> {
typedef Sparse ret; };
20 template<
typename SparseLhsType,
typename DenseRhsType,
typename DenseResType,
23 bool ColPerCol = ((DenseRhsType::Flags&
RowMajorBit)==0) || DenseRhsType::ColsAtCompileTime==1>
24 struct sparse_time_dense_product_impl;
26 template<
typename SparseLhsType,
typename DenseRhsType,
typename DenseResType>
27 struct sparse_time_dense_product_impl<SparseLhsType,DenseRhsType,DenseResType, typename DenseResType::Scalar,
RowMajor, true>
29 typedef typename internal::remove_all<SparseLhsType>::type Lhs;
30 typedef typename internal::remove_all<DenseRhsType>::type Rhs;
31 typedef typename internal::remove_all<DenseResType>::type Res;
32 typedef typename evaluator<Lhs>::InnerIterator LhsInnerIterator;
33 typedef evaluator<Lhs> LhsEval;
34 static void run(
const SparseLhsType& lhs,
const DenseRhsType& rhs, DenseResType& res,
const typename Res::Scalar& alpha)
38 Index n = lhs.outerSize();
39 #ifdef EIGEN_HAS_OPENMP 40 Eigen::initParallel();
41 Index threads = Eigen::nbThreads();
44 for(
Index c=0; c<rhs.cols(); ++c)
46 #ifdef EIGEN_HAS_OPENMP 49 if(threads>1 && lhsEval.nonZerosEstimate() > 20000)
51 #pragma omp parallel for schedule(dynamic,(n+threads*4-1)/(threads*4)) num_threads(threads) 52 for(
Index i=0; i<n; ++i)
53 processRow(lhsEval,rhs,res,alpha,i,c);
58 for(
Index i=0; i<n; ++i)
59 processRow(lhsEval,rhs,res,alpha,i,c);
64 static void processRow(
const LhsEval& lhsEval,
const DenseRhsType& rhs, DenseResType& res,
const typename Res::Scalar& alpha,
Index i,
Index col)
66 typename Res::Scalar tmp(0);
67 for(LhsInnerIterator it(lhsEval,i); it ;++it)
68 tmp += it.value() * rhs.coeff(it.index(),col);
69 res.coeffRef(i,col) += alpha * tmp;
85 template<
typename SparseLhsType,
typename DenseRhsType,
typename DenseResType,
typename AlphaType>
86 struct sparse_time_dense_product_impl<SparseLhsType,DenseRhsType,DenseResType, AlphaType,
ColMajor, true>
88 typedef typename internal::remove_all<SparseLhsType>::type Lhs;
89 typedef typename internal::remove_all<DenseRhsType>::type Rhs;
90 typedef typename internal::remove_all<DenseResType>::type Res;
91 typedef evaluator<Lhs> LhsEval;
92 typedef typename LhsEval::InnerIterator LhsInnerIterator;
93 static void run(
const SparseLhsType& lhs,
const DenseRhsType& rhs, DenseResType& res,
const AlphaType& alpha)
96 for(
Index c=0; c<rhs.cols(); ++c)
98 for(
Index j=0; j<lhs.outerSize(); ++j)
101 typename ScalarBinaryOpTraits<AlphaType, typename Rhs::Scalar>::ReturnType rhs_j(alpha * rhs.coeff(j,c));
102 for(LhsInnerIterator it(lhsEval,j); it ;++it)
103 res.coeffRef(it.index(),c) += it.value() * rhs_j;
109 template<
typename SparseLhsType,
typename DenseRhsType,
typename DenseResType>
110 struct sparse_time_dense_product_impl<SparseLhsType,DenseRhsType,DenseResType, typename DenseResType::Scalar,
RowMajor, false>
112 typedef typename internal::remove_all<SparseLhsType>::type Lhs;
113 typedef typename internal::remove_all<DenseRhsType>::type Rhs;
114 typedef typename internal::remove_all<DenseResType>::type Res;
115 typedef evaluator<Lhs> LhsEval;
116 typedef typename LhsEval::InnerIterator LhsInnerIterator;
117 static void run(
const SparseLhsType& lhs,
const DenseRhsType& rhs, DenseResType& res,
const typename Res::Scalar& alpha)
119 Index n = lhs.rows();
120 LhsEval lhsEval(lhs);
122 #ifdef EIGEN_HAS_OPENMP 123 Eigen::initParallel();
124 Index threads = Eigen::nbThreads();
127 if(threads>1 && lhsEval.nonZerosEstimate()*rhs.cols() > 20000)
129 #pragma omp parallel for schedule(dynamic,(n+threads*4-1)/(threads*4)) num_threads(threads) 130 for(
Index i=0; i<n; ++i)
131 processRow(lhsEval,rhs,res,alpha,i);
136 for(
Index i=0; i<n; ++i)
137 processRow(lhsEval, rhs, res, alpha, i);
141 static void processRow(
const LhsEval& lhsEval,
const DenseRhsType& rhs, Res& res,
const typename Res::Scalar& alpha,
Index i)
143 typename Res::RowXpr res_i(res.row(i));
144 for(LhsInnerIterator it(lhsEval,i); it ;++it)
145 res_i += (alpha*it.value()) * rhs.row(it.index());
149 template<
typename SparseLhsType,
typename DenseRhsType,
typename DenseResType>
150 struct sparse_time_dense_product_impl<SparseLhsType,DenseRhsType,DenseResType, typename DenseResType::Scalar,
ColMajor, false>
152 typedef typename internal::remove_all<SparseLhsType>::type Lhs;
153 typedef typename internal::remove_all<DenseRhsType>::type Rhs;
154 typedef typename internal::remove_all<DenseResType>::type Res;
155 typedef typename evaluator<Lhs>::InnerIterator LhsInnerIterator;
156 static void run(
const SparseLhsType& lhs,
const DenseRhsType& rhs, DenseResType& res,
const typename Res::Scalar& alpha)
158 evaluator<Lhs> lhsEval(lhs);
159 for(
Index j=0; j<lhs.outerSize(); ++j)
161 typename Rhs::ConstRowXpr rhs_j(rhs.row(j));
162 for(LhsInnerIterator it(lhsEval,j); it ;++it)
163 res.row(it.index()) += (alpha*it.value()) * rhs_j;
168 template<
typename SparseLhsType,
typename DenseRhsType,
typename DenseResType,
typename AlphaType>
169 inline void sparse_time_dense_product(
const SparseLhsType& lhs,
const DenseRhsType& rhs, DenseResType& res,
const AlphaType& alpha)
171 sparse_time_dense_product_impl<SparseLhsType,DenseRhsType,DenseResType, AlphaType>::run(lhs, rhs, res, alpha);
178 template<
typename Lhs,
typename Rhs,
int ProductType>
179 struct generic_product_impl<Lhs, Rhs, SparseShape, DenseShape, ProductType>
180 : generic_product_impl_base<Lhs,Rhs,generic_product_impl<Lhs,Rhs,SparseShape,DenseShape,ProductType> >
182 typedef typename Product<Lhs,Rhs>::Scalar Scalar;
184 template<
typename Dest>
185 static void scaleAndAddTo(Dest& dst,
const Lhs& lhs,
const Rhs& rhs,
const Scalar& alpha)
187 typedef typename nested_eval<Lhs,((Rhs::Flags&RowMajorBit)==0) ? 1 : Rhs::ColsAtCompileTime>::type LhsNested;
188 typedef typename nested_eval<Rhs,((Lhs::Flags&RowMajorBit)==0) ? 1 : Dynamic>::type RhsNested;
189 LhsNested lhsNested(lhs);
190 RhsNested rhsNested(rhs);
191 internal::sparse_time_dense_product(lhsNested, rhsNested, dst, alpha);
195 template<
typename Lhs,
typename Rhs,
int ProductType>
196 struct generic_product_impl<Lhs, Rhs, SparseTriangularShape, DenseShape, ProductType>
197 : generic_product_impl<Lhs, Rhs, SparseShape, DenseShape, ProductType>
200 template<
typename Lhs,
typename Rhs,
int ProductType>
201 struct generic_product_impl<Lhs, Rhs, DenseShape, SparseShape, ProductType>
202 : generic_product_impl_base<Lhs,Rhs,generic_product_impl<Lhs,Rhs,DenseShape,SparseShape,ProductType> >
204 typedef typename Product<Lhs,Rhs>::Scalar Scalar;
206 template<
typename Dst>
207 static void scaleAndAddTo(Dst& dst,
const Lhs& lhs,
const Rhs& rhs,
const Scalar& alpha)
209 typedef typename nested_eval<Lhs,((Rhs::Flags&RowMajorBit)==0) ? Dynamic : 1>::type LhsNested;
210 typedef typename nested_eval<Rhs,((Lhs::Flags&RowMajorBit)==RowMajorBit) ? 1 : Lhs::RowsAtCompileTime>::type RhsNested;
211 LhsNested lhsNested(lhs);
212 RhsNested rhsNested(rhs);
215 Transpose<Dst> dstT(dst);
216 internal::sparse_time_dense_product(rhsNested.transpose(), lhsNested.transpose(), dstT, alpha);
220 template<
typename Lhs,
typename Rhs,
int ProductType>
221 struct generic_product_impl<Lhs, Rhs, DenseShape, SparseTriangularShape, ProductType>
222 : generic_product_impl<Lhs, Rhs, DenseShape, SparseShape, ProductType>
225 template<
typename LhsT,
typename RhsT,
bool NeedToTranspose>
226 struct sparse_dense_outer_product_evaluator
229 typedef typename conditional<NeedToTranspose,RhsT,LhsT>::type Lhs1;
230 typedef typename conditional<NeedToTranspose,LhsT,RhsT>::type ActualRhs;
231 typedef Product<LhsT,RhsT,DefaultProduct> ProdXprType;
235 typedef typename conditional<is_same<typename internal::traits<Lhs1>::StorageKind,Sparse>::value,
236 Lhs1, SparseView<Lhs1> >::type ActualLhs;
237 typedef typename conditional<is_same<typename internal::traits<Lhs1>::StorageKind,Sparse>::value,
238 Lhs1
const&, SparseView<Lhs1> >::type LhsArg;
240 typedef evaluator<ActualLhs> LhsEval;
241 typedef evaluator<ActualRhs> RhsEval;
242 typedef typename evaluator<ActualLhs>::InnerIterator LhsIterator;
243 typedef typename ProdXprType::Scalar Scalar;
251 class InnerIterator :
public LhsIterator
254 InnerIterator(
const sparse_dense_outer_product_evaluator &xprEval,
Index outer)
255 : LhsIterator(xprEval.m_lhsXprImpl, 0),
258 m_factor(get(xprEval.m_rhsXprImpl, outer, typename
internal::traits<ActualRhs>::StorageKind() ))
261 EIGEN_STRONG_INLINE
Index outer()
const {
return m_outer; }
262 EIGEN_STRONG_INLINE
Index row()
const {
return NeedToTranspose ? m_outer : LhsIterator::index(); }
263 EIGEN_STRONG_INLINE
Index col()
const {
return NeedToTranspose ? LhsIterator::index() : m_outer; }
265 EIGEN_STRONG_INLINE Scalar value()
const {
return LhsIterator::value() * m_factor; }
266 EIGEN_STRONG_INLINE
operator bool()
const {
return LhsIterator::operator bool() && (!m_empty); }
269 Scalar
get(
const RhsEval &rhs,
Index outer, Dense = Dense())
const 271 return rhs.coeff(outer);
274 Scalar
get(
const RhsEval &rhs,
Index outer, Sparse = Sparse())
276 typename RhsEval::InnerIterator it(rhs, outer);
277 if (it && it.index()==0 && it.value()!=Scalar(0))
288 sparse_dense_outer_product_evaluator(
const Lhs1 &lhs,
const ActualRhs &rhs)
289 : m_lhs(lhs), m_lhsXprImpl(m_lhs), m_rhsXprImpl(rhs)
291 EIGEN_INTERNAL_CHECK_COST_VALUE(CoeffReadCost);
295 sparse_dense_outer_product_evaluator(
const ActualRhs &rhs,
const Lhs1 &lhs)
296 : m_lhs(lhs), m_lhsXprImpl(m_lhs), m_rhsXprImpl(rhs)
298 EIGEN_INTERNAL_CHECK_COST_VALUE(CoeffReadCost);
303 evaluator<ActualLhs> m_lhsXprImpl;
304 evaluator<ActualRhs> m_rhsXprImpl;
308 template<
typename Lhs,
typename Rhs>
309 struct product_evaluator<Product<Lhs, Rhs, DefaultProduct>, OuterProduct, SparseShape, DenseShape>
310 : sparse_dense_outer_product_evaluator<Lhs,Rhs, Lhs::IsRowMajor>
312 typedef sparse_dense_outer_product_evaluator<Lhs,Rhs, Lhs::IsRowMajor> Base;
314 typedef Product<Lhs, Rhs> XprType;
315 typedef typename XprType::PlainObject PlainObject;
317 explicit product_evaluator(
const XprType& xpr)
318 : Base(xpr.lhs(), xpr.rhs())
323 template<
typename Lhs,
typename Rhs>
324 struct product_evaluator<Product<Lhs, Rhs, DefaultProduct>, OuterProduct, DenseShape, SparseShape>
325 : sparse_dense_outer_product_evaluator<Lhs,Rhs, Rhs::IsRowMajor>
327 typedef sparse_dense_outer_product_evaluator<Lhs,Rhs, Rhs::IsRowMajor> Base;
329 typedef Product<Lhs, Rhs> XprType;
330 typedef typename XprType::PlainObject PlainObject;
332 explicit product_evaluator(
const XprType& xpr)
333 : Base(xpr.lhs(), xpr.rhs())
342 #endif // EIGEN_SPARSEDENSEPRODUCT_H Definition: Constants.h:319
const int HugeCost
Definition: Constants.h:44
Namespace containing all symbols from the Eigen library.
Definition: Core:141
const unsigned int RowMajorBit
Definition: Constants.h:66
EIGEN_DEFAULT_DENSE_INDEX_TYPE Index
The Index type as used for the API.
Definition: Meta.h:74
Definition: Eigen_Colamd.h:50
Definition: Constants.h:321