Branch data Line data Source code
1 : : // This file is part of Eigen, a lightweight C++ template library
2 : : // for linear algebra.
3 : : //
4 : : // Copyright (C) 2006-2008 Benoit Jacob <jacob.benoit.1@gmail.com>
5 : : // Copyright (C) 2008-2010 Gael Guennebaud <gael.guennebaud@inria.fr>
6 : : // Copyright (C) 2011 Jitse Niesen <jitse@maths.leeds.ac.uk>
7 : : //
8 : : // This Source Code Form is subject to the terms of the Mozilla
9 : : // Public License v. 2.0. If a copy of the MPL was not distributed
10 : : // with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
11 : :
12 : :
13 : : #ifndef EIGEN_PRODUCTEVALUATORS_H
14 : : #define EIGEN_PRODUCTEVALUATORS_H
15 : :
16 : : namespace Eigen {
17 : :
18 : : namespace internal {
19 : :
20 : : /** \internal
21 : : * Evaluator of a product expression.
22 : : * Since products require special treatments to handle all possible cases,
23 : : * we simply defer the evaluation logic to a product_evaluator class
24 : : * which offers more partial specialization possibilities.
25 : : *
26 : : * \sa class product_evaluator
27 : : */
28 : : template<typename Lhs, typename Rhs, int Options>
29 : : struct evaluator<Product<Lhs, Rhs, Options> >
30 : : : public product_evaluator<Product<Lhs, Rhs, Options> >
31 : : {
32 : : typedef Product<Lhs, Rhs, Options> XprType;
33 : : typedef product_evaluator<XprType> Base;
34 : :
35 : 2068 : EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE explicit evaluator(const XprType& xpr) : Base(xpr) {}
36 : : };
37 : :
38 : : // Catch "scalar * ( A * B )" and transform it to "(A*scalar) * B"
39 : : // TODO we should apply that rule only if that's really helpful
40 : : template<typename Lhs, typename Rhs, typename Scalar1, typename Scalar2, typename Plain1>
41 : : struct evaluator_assume_aliasing<CwiseBinaryOp<internal::scalar_product_op<Scalar1,Scalar2>,
42 : : const CwiseNullaryOp<internal::scalar_constant_op<Scalar1>, Plain1>,
43 : : const Product<Lhs, Rhs, DefaultProduct> > >
44 : : {
45 : : static const bool value = true;
46 : : };
47 : : template<typename Lhs, typename Rhs, typename Scalar1, typename Scalar2, typename Plain1>
48 : : struct evaluator<CwiseBinaryOp<internal::scalar_product_op<Scalar1,Scalar2>,
49 : : const CwiseNullaryOp<internal::scalar_constant_op<Scalar1>, Plain1>,
50 : : const Product<Lhs, Rhs, DefaultProduct> > >
51 : : : public evaluator<Product<EIGEN_SCALAR_BINARYOP_EXPR_RETURN_TYPE(Scalar1,Lhs,product), Rhs, DefaultProduct> >
52 : : {
53 : : typedef CwiseBinaryOp<internal::scalar_product_op<Scalar1,Scalar2>,
54 : : const CwiseNullaryOp<internal::scalar_constant_op<Scalar1>, Plain1>,
55 : : const Product<Lhs, Rhs, DefaultProduct> > XprType;
56 : : typedef evaluator<Product<EIGEN_SCALAR_BINARYOP_EXPR_RETURN_TYPE(Scalar1,Lhs,product), Rhs, DefaultProduct> > Base;
57 : :
58 : 304 : EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE explicit evaluator(const XprType& xpr)
59 : 304 : : Base(xpr.lhs().functor().m_other * xpr.rhs().lhs() * xpr.rhs().rhs())
60 : 304 : {}
61 : : };
62 : :
63 : :
64 : : template<typename Lhs, typename Rhs, int DiagIndex>
65 : : struct evaluator<Diagonal<const Product<Lhs, Rhs, DefaultProduct>, DiagIndex> >
66 : : : public evaluator<Diagonal<const Product<Lhs, Rhs, LazyProduct>, DiagIndex> >
67 : : {
68 : : typedef Diagonal<const Product<Lhs, Rhs, DefaultProduct>, DiagIndex> XprType;
69 : : typedef evaluator<Diagonal<const Product<Lhs, Rhs, LazyProduct>, DiagIndex> > Base;
70 : :
71 : : EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE explicit evaluator(const XprType& xpr)
72 : : : Base(Diagonal<const Product<Lhs, Rhs, LazyProduct>, DiagIndex>(
73 : : Product<Lhs, Rhs, LazyProduct>(xpr.nestedExpression().lhs(), xpr.nestedExpression().rhs()),
74 : : xpr.index() ))
75 : : {}
76 : : };
77 : :
78 : :
79 : : // Helper class to perform a matrix product with the destination at hand.
80 : : // Depending on the sizes of the factors, there are different evaluation strategies
81 : : // as controlled by internal::product_type.
82 : : template< typename Lhs, typename Rhs,
83 : : typename LhsShape = typename evaluator_traits<Lhs>::Shape,
84 : : typename RhsShape = typename evaluator_traits<Rhs>::Shape,
85 : : int ProductType = internal::product_type<Lhs,Rhs>::value>
86 : : struct generic_product_impl;
87 : :
88 : : template<typename Lhs, typename Rhs>
89 : : struct evaluator_assume_aliasing<Product<Lhs, Rhs, DefaultProduct> > {
90 : : static const bool value = true;
91 : : };
92 : :
93 : : // This is the default evaluator implementation for products:
94 : : // It creates a temporary and call generic_product_impl
95 : : template<typename Lhs, typename Rhs, int Options, int ProductTag, typename LhsShape, typename RhsShape>
96 : : struct product_evaluator<Product<Lhs, Rhs, Options>, ProductTag, LhsShape, RhsShape>
97 : : : public evaluator<typename Product<Lhs, Rhs, Options>::PlainObject>
98 : : {
99 : : typedef Product<Lhs, Rhs, Options> XprType;
100 : : typedef typename XprType::PlainObject PlainObject;
101 : : typedef evaluator<PlainObject> Base;
102 : : enum {
103 : : Flags = Base::Flags | EvalBeforeNestingBit
104 : : };
105 : :
106 : : EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE
107 : 152 : explicit product_evaluator(const XprType& xpr)
108 : 152 : : m_result(xpr.rows(), xpr.cols())
109 : : {
110 : 152 : ::new (static_cast<Base*>(this)) Base(m_result);
111 : :
112 : : // FIXME shall we handle nested_eval here?,
113 : : // if so, then we must take care at removing the call to nested_eval in the specializations (e.g., in permutation_matrix_product, transposition_matrix_product, etc.)
114 : : // typedef typename internal::nested_eval<Lhs,Rhs::ColsAtCompileTime>::type LhsNested;
115 : : // typedef typename internal::nested_eval<Rhs,Lhs::RowsAtCompileTime>::type RhsNested;
116 : : // typedef typename internal::remove_all<LhsNested>::type LhsNestedCleaned;
117 : : // typedef typename internal::remove_all<RhsNested>::type RhsNestedCleaned;
118 : : //
119 : : // const LhsNested lhs(xpr.lhs());
120 : : // const RhsNested rhs(xpr.rhs());
121 : : //
122 : : // generic_product_impl<LhsNestedCleaned, RhsNestedCleaned>::evalTo(m_result, lhs, rhs);
123 : :
124 : 152 : generic_product_impl<Lhs, Rhs, LhsShape, RhsShape, ProductTag>::evalTo(m_result, xpr.lhs(), xpr.rhs());
125 : 152 : }
126 : :
127 : : protected:
128 : : PlainObject m_result;
129 : : };
130 : :
131 : : // The following three shortcuts are enabled only if the scalar types match exactly.
132 : : // TODO: we could enable them for different scalar types when the product is not vectorized.
133 : :
134 : : // Dense = Product
135 : : template< typename DstXprType, typename Lhs, typename Rhs, int Options, typename Scalar>
136 : : struct Assignment<DstXprType, Product<Lhs,Rhs,Options>, internal::assign_op<Scalar,Scalar>, Dense2Dense,
137 : : typename enable_if<(Options==DefaultProduct || Options==AliasFreeProduct)>::type>
138 : : {
139 : : typedef Product<Lhs,Rhs,Options> SrcXprType;
140 : : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE
141 : 1612 : void run(DstXprType &dst, const SrcXprType &src, const internal::assign_op<Scalar,Scalar> &)
142 : : {
143 : 1612 : Index dstRows = src.rows();
144 : 1612 : Index dstCols = src.cols();
145 : 1612 : if((dst.rows()!=dstRows) || (dst.cols()!=dstCols))
146 : 0 : dst.resize(dstRows, dstCols);
147 : : // FIXME shall we handle nested_eval here?
148 : 1612 : generic_product_impl<Lhs, Rhs>::evalTo(dst, src.lhs(), src.rhs());
149 : 1612 : }
150 : : };
151 : :
152 : : // Dense += Product
153 : : template< typename DstXprType, typename Lhs, typename Rhs, int Options, typename Scalar>
154 : : struct Assignment<DstXprType, Product<Lhs,Rhs,Options>, internal::add_assign_op<Scalar,Scalar>, Dense2Dense,
155 : : typename enable_if<(Options==DefaultProduct || Options==AliasFreeProduct)>::type>
156 : : {
157 : : typedef Product<Lhs,Rhs,Options> SrcXprType;
158 : : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE
159 : : void run(DstXprType &dst, const SrcXprType &src, const internal::add_assign_op<Scalar,Scalar> &)
160 : : {
161 : : eigen_assert(dst.rows() == src.rows() && dst.cols() == src.cols());
162 : : // FIXME shall we handle nested_eval here?
163 : : generic_product_impl<Lhs, Rhs>::addTo(dst, src.lhs(), src.rhs());
164 : : }
165 : : };
166 : :
167 : : // Dense -= Product
168 : : template< typename DstXprType, typename Lhs, typename Rhs, int Options, typename Scalar>
169 : : struct Assignment<DstXprType, Product<Lhs,Rhs,Options>, internal::sub_assign_op<Scalar,Scalar>, Dense2Dense,
170 : : typename enable_if<(Options==DefaultProduct || Options==AliasFreeProduct)>::type>
171 : : {
172 : : typedef Product<Lhs,Rhs,Options> SrcXprType;
173 : : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE
174 : : void run(DstXprType &dst, const SrcXprType &src, const internal::sub_assign_op<Scalar,Scalar> &)
175 : : {
176 : : eigen_assert(dst.rows() == src.rows() && dst.cols() == src.cols());
177 : : // FIXME shall we handle nested_eval here?
178 : : generic_product_impl<Lhs, Rhs>::subTo(dst, src.lhs(), src.rhs());
179 : : }
180 : : };
181 : :
182 : :
183 : : // Dense ?= scalar * Product
184 : : // TODO we should apply that rule if that's really helpful
185 : : // for instance, this is not good for inner products
186 : : template< typename DstXprType, typename Lhs, typename Rhs, typename AssignFunc, typename Scalar, typename ScalarBis, typename Plain>
187 : : struct Assignment<DstXprType, CwiseBinaryOp<internal::scalar_product_op<ScalarBis,Scalar>, const CwiseNullaryOp<internal::scalar_constant_op<ScalarBis>,Plain>,
188 : : const Product<Lhs,Rhs,DefaultProduct> >, AssignFunc, Dense2Dense>
189 : : {
190 : : typedef CwiseBinaryOp<internal::scalar_product_op<ScalarBis,Scalar>,
191 : : const CwiseNullaryOp<internal::scalar_constant_op<ScalarBis>,Plain>,
192 : : const Product<Lhs,Rhs,DefaultProduct> > SrcXprType;
193 : : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE
194 : : void run(DstXprType &dst, const SrcXprType &src, const AssignFunc& func)
195 : : {
196 : : call_assignment_no_alias(dst, (src.lhs().functor().m_other * src.rhs().lhs())*src.rhs().rhs(), func);
197 : : }
198 : : };
199 : :
200 : : //----------------------------------------
201 : : // Catch "Dense ?= xpr + Product<>" expression to save one temporary
202 : : // FIXME we could probably enable these rules for any product, i.e., not only Dense and DefaultProduct
203 : :
204 : : template<typename OtherXpr, typename Lhs, typename Rhs>
205 : : struct evaluator_assume_aliasing<CwiseBinaryOp<internal::scalar_sum_op<typename OtherXpr::Scalar,typename Product<Lhs,Rhs,DefaultProduct>::Scalar>, const OtherXpr,
206 : : const Product<Lhs,Rhs,DefaultProduct> >, DenseShape > {
207 : : static const bool value = true;
208 : : };
209 : :
210 : : template<typename OtherXpr, typename Lhs, typename Rhs>
211 : : struct evaluator_assume_aliasing<CwiseBinaryOp<internal::scalar_difference_op<typename OtherXpr::Scalar,typename Product<Lhs,Rhs,DefaultProduct>::Scalar>, const OtherXpr,
212 : : const Product<Lhs,Rhs,DefaultProduct> >, DenseShape > {
213 : : static const bool value = true;
214 : : };
215 : :
216 : : template<typename DstXprType, typename OtherXpr, typename ProductType, typename Func1, typename Func2>
217 : : struct assignment_from_xpr_op_product
218 : : {
219 : : template<typename SrcXprType, typename InitialFunc>
220 : : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE
221 : : void run(DstXprType &dst, const SrcXprType &src, const InitialFunc& /*func*/)
222 : : {
223 : : call_assignment_no_alias(dst, src.lhs(), Func1());
224 : : call_assignment_no_alias(dst, src.rhs(), Func2());
225 : : }
226 : : };
227 : :
228 : : #define EIGEN_CATCH_ASSIGN_XPR_OP_PRODUCT(ASSIGN_OP,BINOP,ASSIGN_OP2) \
229 : : template< typename DstXprType, typename OtherXpr, typename Lhs, typename Rhs, typename DstScalar, typename SrcScalar, typename OtherScalar,typename ProdScalar> \
230 : : struct Assignment<DstXprType, CwiseBinaryOp<internal::BINOP<OtherScalar,ProdScalar>, const OtherXpr, \
231 : : const Product<Lhs,Rhs,DefaultProduct> >, internal::ASSIGN_OP<DstScalar,SrcScalar>, Dense2Dense> \
232 : : : assignment_from_xpr_op_product<DstXprType, OtherXpr, Product<Lhs,Rhs,DefaultProduct>, internal::ASSIGN_OP<DstScalar,OtherScalar>, internal::ASSIGN_OP2<DstScalar,ProdScalar> > \
233 : : {}
234 : :
235 : : EIGEN_CATCH_ASSIGN_XPR_OP_PRODUCT(assign_op, scalar_sum_op,add_assign_op);
236 : : EIGEN_CATCH_ASSIGN_XPR_OP_PRODUCT(add_assign_op,scalar_sum_op,add_assign_op);
237 : : EIGEN_CATCH_ASSIGN_XPR_OP_PRODUCT(sub_assign_op,scalar_sum_op,sub_assign_op);
238 : :
239 : : EIGEN_CATCH_ASSIGN_XPR_OP_PRODUCT(assign_op, scalar_difference_op,sub_assign_op);
240 : : EIGEN_CATCH_ASSIGN_XPR_OP_PRODUCT(add_assign_op,scalar_difference_op,sub_assign_op);
241 : : EIGEN_CATCH_ASSIGN_XPR_OP_PRODUCT(sub_assign_op,scalar_difference_op,add_assign_op);
242 : :
243 : : //----------------------------------------
244 : :
245 : : template<typename Lhs, typename Rhs>
246 : : struct generic_product_impl<Lhs,Rhs,DenseShape,DenseShape,InnerProduct>
247 : : {
248 : : template<typename Dst>
249 : : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void evalTo(Dst& dst, const Lhs& lhs, const Rhs& rhs)
250 : : {
251 : : dst.coeffRef(0,0) = (lhs.transpose().cwiseProduct(rhs)).sum();
252 : : }
253 : :
254 : : template<typename Dst>
255 : : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void addTo(Dst& dst, const Lhs& lhs, const Rhs& rhs)
256 : : {
257 : : dst.coeffRef(0,0) += (lhs.transpose().cwiseProduct(rhs)).sum();
258 : : }
259 : :
260 : : template<typename Dst>
261 : : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void subTo(Dst& dst, const Lhs& lhs, const Rhs& rhs)
262 : : { dst.coeffRef(0,0) -= (lhs.transpose().cwiseProduct(rhs)).sum(); }
263 : : };
264 : :
265 : :
266 : : /***********************************************************************
267 : : * Implementation of outer dense * dense vector product
268 : : ***********************************************************************/
269 : :
270 : : // Column major result
271 : : template<typename Dst, typename Lhs, typename Rhs, typename Func>
272 : : void EIGEN_DEVICE_FUNC outer_product_selector_run(Dst& dst, const Lhs &lhs, const Rhs &rhs, const Func& func, const false_type&)
273 : : {
274 : : evaluator<Rhs> rhsEval(rhs);
275 : : ei_declare_local_nested_eval(Lhs,lhs,Rhs::SizeAtCompileTime,actual_lhs);
276 : : // FIXME if cols is large enough, then it might be useful to make sure that lhs is sequentially stored
277 : : // FIXME not very good if rhs is real and lhs complex while alpha is real too
278 : : const Index cols = dst.cols();
279 : : for (Index j=0; j<cols; ++j)
280 : : func(dst.col(j), rhsEval.coeff(Index(0),j) * actual_lhs);
281 : : }
282 : :
283 : : // Row major result
284 : : template<typename Dst, typename Lhs, typename Rhs, typename Func>
285 : : void EIGEN_DEVICE_FUNC outer_product_selector_run(Dst& dst, const Lhs &lhs, const Rhs &rhs, const Func& func, const true_type&)
286 : : {
287 : : evaluator<Lhs> lhsEval(lhs);
288 : : ei_declare_local_nested_eval(Rhs,rhs,Lhs::SizeAtCompileTime,actual_rhs);
289 : : // FIXME if rows is large enough, then it might be useful to make sure that rhs is sequentially stored
290 : : // FIXME not very good if lhs is real and rhs complex while alpha is real too
291 : : const Index rows = dst.rows();
292 : : for (Index i=0; i<rows; ++i)
293 : : func(dst.row(i), lhsEval.coeff(i,Index(0)) * actual_rhs);
294 : : }
295 : :
296 : : template<typename Lhs, typename Rhs>
297 : : struct generic_product_impl<Lhs,Rhs,DenseShape,DenseShape,OuterProduct>
298 : : {
299 : : template<typename T> struct is_row_major : internal::conditional<(int(T::Flags)&RowMajorBit), internal::true_type, internal::false_type>::type {};
300 : : typedef typename Product<Lhs,Rhs>::Scalar Scalar;
301 : :
302 : : // TODO it would be nice to be able to exploit our *_assign_op functors for that purpose
303 : : struct set { template<typename Dst, typename Src> EIGEN_DEVICE_FUNC void operator()(const Dst& dst, const Src& src) const { dst.const_cast_derived() = src; } };
304 : : struct add { template<typename Dst, typename Src> EIGEN_DEVICE_FUNC void operator()(const Dst& dst, const Src& src) const { dst.const_cast_derived() += src; } };
305 : : struct sub { template<typename Dst, typename Src> EIGEN_DEVICE_FUNC void operator()(const Dst& dst, const Src& src) const { dst.const_cast_derived() -= src; } };
306 : : struct adds {
307 : : Scalar m_scale;
308 : : explicit adds(const Scalar& s) : m_scale(s) {}
309 : : template<typename Dst, typename Src> void EIGEN_DEVICE_FUNC operator()(const Dst& dst, const Src& src) const {
310 : : dst.const_cast_derived() += m_scale * src;
311 : : }
312 : : };
313 : :
314 : : template<typename Dst>
315 : : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void evalTo(Dst& dst, const Lhs& lhs, const Rhs& rhs)
316 : : {
317 : : internal::outer_product_selector_run(dst, lhs, rhs, set(), is_row_major<Dst>());
318 : : }
319 : :
320 : : template<typename Dst>
321 : : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void addTo(Dst& dst, const Lhs& lhs, const Rhs& rhs)
322 : : {
323 : : internal::outer_product_selector_run(dst, lhs, rhs, add(), is_row_major<Dst>());
324 : : }
325 : :
326 : : template<typename Dst>
327 : : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void subTo(Dst& dst, const Lhs& lhs, const Rhs& rhs)
328 : : {
329 : : internal::outer_product_selector_run(dst, lhs, rhs, sub(), is_row_major<Dst>());
330 : : }
331 : :
332 : : template<typename Dst>
333 : : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void scaleAndAddTo(Dst& dst, const Lhs& lhs, const Rhs& rhs, const Scalar& alpha)
334 : : {
335 : : internal::outer_product_selector_run(dst, lhs, rhs, adds(alpha), is_row_major<Dst>());
336 : : }
337 : :
338 : : };
339 : :
340 : :
341 : : // This base class provides default implementations for evalTo, addTo, subTo, in terms of scaleAndAddTo
342 : : template<typename Lhs, typename Rhs, typename Derived>
343 : : struct generic_product_impl_base
344 : : {
345 : : typedef typename Product<Lhs,Rhs>::Scalar Scalar;
346 : :
347 : : template<typename Dst>
348 : : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void evalTo(Dst& dst, const Lhs& lhs, const Rhs& rhs)
349 : : { dst.setZero(); scaleAndAddTo(dst, lhs, rhs, Scalar(1)); }
350 : :
351 : : template<typename Dst>
352 : : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void addTo(Dst& dst, const Lhs& lhs, const Rhs& rhs)
353 : : { scaleAndAddTo(dst,lhs, rhs, Scalar(1)); }
354 : :
355 : : template<typename Dst>
356 : : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void subTo(Dst& dst, const Lhs& lhs, const Rhs& rhs)
357 : : { scaleAndAddTo(dst, lhs, rhs, Scalar(-1)); }
358 : :
359 : : template<typename Dst>
360 : : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void scaleAndAddTo(Dst& dst, const Lhs& lhs, const Rhs& rhs, const Scalar& alpha)
361 : : { Derived::scaleAndAddTo(dst,lhs,rhs,alpha); }
362 : :
363 : : };
364 : :
365 : : template<typename Lhs, typename Rhs>
366 : : struct generic_product_impl<Lhs,Rhs,DenseShape,DenseShape,GemvProduct>
367 : : : generic_product_impl_base<Lhs,Rhs,generic_product_impl<Lhs,Rhs,DenseShape,DenseShape,GemvProduct> >
368 : : {
369 : : typedef typename nested_eval<Lhs,1>::type LhsNested;
370 : : typedef typename nested_eval<Rhs,1>::type RhsNested;
371 : : typedef typename Product<Lhs,Rhs>::Scalar Scalar;
372 : : enum { Side = Lhs::IsVectorAtCompileTime ? OnTheLeft : OnTheRight };
373 : : typedef typename internal::remove_all<typename internal::conditional<int(Side)==OnTheRight,LhsNested,RhsNested>::type>::type MatrixType;
374 : :
375 : : template<typename Dest>
376 : : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void scaleAndAddTo(Dest& dst, const Lhs& lhs, const Rhs& rhs, const Scalar& alpha)
377 : : {
378 : : // Fallback to inner product if both the lhs and rhs is a runtime vector.
379 : : if (lhs.rows() == 1 && rhs.cols() == 1) {
380 : : dst.coeffRef(0,0) += alpha * lhs.row(0).conjugate().dot(rhs.col(0));
381 : : return;
382 : : }
383 : : LhsNested actual_lhs(lhs);
384 : : RhsNested actual_rhs(rhs);
385 : : internal::gemv_dense_selector<Side,
386 : : (int(MatrixType::Flags)&RowMajorBit) ? RowMajor : ColMajor,
387 : : bool(internal::blas_traits<MatrixType>::HasUsableDirectAccess)
388 : : >::run(actual_lhs, actual_rhs, dst, alpha);
389 : : }
390 : : };
391 : :
392 : : template<typename Lhs, typename Rhs>
393 : : struct generic_product_impl<Lhs,Rhs,DenseShape,DenseShape,CoeffBasedProductMode>
394 : : {
395 : : typedef typename Product<Lhs,Rhs>::Scalar Scalar;
396 : :
397 : : template<typename Dst>
398 : 1764 : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void evalTo(Dst& dst, const Lhs& lhs, const Rhs& rhs)
399 : : {
400 : : // Same as: dst.noalias() = lhs.lazyProduct(rhs);
401 : : // but easier on the compiler side
402 : 1764 : call_assignment_no_alias(dst, lhs.lazyProduct(rhs), internal::assign_op<typename Dst::Scalar,Scalar>());
403 : 1764 : }
404 : :
405 : : template<typename Dst>
406 : : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void addTo(Dst& dst, const Lhs& lhs, const Rhs& rhs)
407 : : {
408 : : // dst.noalias() += lhs.lazyProduct(rhs);
409 : : call_assignment_no_alias(dst, lhs.lazyProduct(rhs), internal::add_assign_op<typename Dst::Scalar,Scalar>());
410 : : }
411 : :
412 : : template<typename Dst>
413 : : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void subTo(Dst& dst, const Lhs& lhs, const Rhs& rhs)
414 : : {
415 : : // dst.noalias() -= lhs.lazyProduct(rhs);
416 : : call_assignment_no_alias(dst, lhs.lazyProduct(rhs), internal::sub_assign_op<typename Dst::Scalar,Scalar>());
417 : : }
418 : :
419 : : // This is a special evaluation path called from generic_product_impl<...,GemmProduct> in file GeneralMatrixMatrix.h
420 : : // This variant tries to extract scalar multiples from both the LHS and RHS and factor them out. For instance:
421 : : // dst {,+,-}= (s1*A)*(B*s2)
422 : : // will be rewritten as:
423 : : // dst {,+,-}= (s1*s2) * (A.lazyProduct(B))
424 : : // There are at least four benefits of doing so:
425 : : // 1 - huge performance gain for heap-allocated matrix types as it save costly allocations.
426 : : // 2 - it is faster than simply by-passing the heap allocation through stack allocation.
427 : : // 3 - it makes this fallback consistent with the heavy GEMM routine.
428 : : // 4 - it fully by-passes huge stack allocation attempts when multiplying huge fixed-size matrices.
429 : : // (see https://stackoverflow.com/questions/54738495)
430 : : // For small fixed sizes matrices, howver, the gains are less obvious, it is sometimes x2 faster, but sometimes x3 slower,
431 : : // and the behavior depends also a lot on the compiler... This is why this re-writting strategy is currently
432 : : // enabled only when falling back from the main GEMM.
433 : : template<typename Dst, typename Func>
434 : : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE
435 : : void eval_dynamic(Dst& dst, const Lhs& lhs, const Rhs& rhs, const Func &func)
436 : : {
437 : : enum {
438 : : HasScalarFactor = blas_traits<Lhs>::HasScalarFactor || blas_traits<Rhs>::HasScalarFactor,
439 : : ConjLhs = blas_traits<Lhs>::NeedToConjugate,
440 : : ConjRhs = blas_traits<Rhs>::NeedToConjugate
441 : : };
442 : : // FIXME: in c++11 this should be auto, and extractScalarFactor should also return auto
443 : : // this is important for real*complex_mat
444 : : Scalar actualAlpha = combine_scalar_factors<Scalar>(lhs, rhs);
445 : :
446 : : eval_dynamic_impl(dst,
447 : : blas_traits<Lhs>::extract(lhs).template conjugateIf<ConjLhs>(),
448 : : blas_traits<Rhs>::extract(rhs).template conjugateIf<ConjRhs>(),
449 : : func,
450 : : actualAlpha,
451 : : typename conditional<HasScalarFactor,true_type,false_type>::type());
452 : : }
453 : :
454 : : protected:
455 : :
456 : : template<typename Dst, typename LhsT, typename RhsT, typename Func, typename Scalar>
457 : : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE
458 : : void eval_dynamic_impl(Dst& dst, const LhsT& lhs, const RhsT& rhs, const Func &func, const Scalar& s /* == 1 */, false_type)
459 : : {
460 : : EIGEN_UNUSED_VARIABLE(s);
461 : : eigen_internal_assert(s==Scalar(1));
462 : : call_restricted_packet_assignment_no_alias(dst, lhs.lazyProduct(rhs), func);
463 : : }
464 : :
465 : : template<typename Dst, typename LhsT, typename RhsT, typename Func, typename Scalar>
466 : : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE
467 : : void eval_dynamic_impl(Dst& dst, const LhsT& lhs, const RhsT& rhs, const Func &func, const Scalar& s, true_type)
468 : : {
469 : : call_restricted_packet_assignment_no_alias(dst, s * lhs.lazyProduct(rhs), func);
470 : : }
471 : : };
472 : :
473 : : // This specialization enforces the use of a coefficient-based evaluation strategy
474 : : template<typename Lhs, typename Rhs>
475 : : struct generic_product_impl<Lhs,Rhs,DenseShape,DenseShape,LazyCoeffBasedProductMode>
476 : : : generic_product_impl<Lhs,Rhs,DenseShape,DenseShape,CoeffBasedProductMode> {};
477 : :
478 : : // Case 2: Evaluate coeff by coeff
479 : : //
480 : : // This is mostly taken from CoeffBasedProduct.h
481 : : // The main difference is that we add an extra argument to the etor_product_*_impl::run() function
482 : : // for the inner dimension of the product, because evaluator object do not know their size.
483 : :
484 : : template<int Traversal, int UnrollingIndex, typename Lhs, typename Rhs, typename RetScalar>
485 : : struct etor_product_coeff_impl;
486 : :
487 : : template<int StorageOrder, int UnrollingIndex, typename Lhs, typename Rhs, typename Packet, int LoadMode>
488 : : struct etor_product_packet_impl;
489 : :
490 : : template<typename Lhs, typename Rhs, int ProductTag>
491 : : struct product_evaluator<Product<Lhs, Rhs, LazyProduct>, ProductTag, DenseShape, DenseShape>
492 : : : evaluator_base<Product<Lhs, Rhs, LazyProduct> >
493 : : {
494 : : typedef Product<Lhs, Rhs, LazyProduct> XprType;
495 : : typedef typename XprType::Scalar Scalar;
496 : : typedef typename XprType::CoeffReturnType CoeffReturnType;
497 : :
498 : : EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE
499 : 1916 : explicit product_evaluator(const XprType& xpr)
500 : 1916 : : m_lhs(xpr.lhs()),
501 : 1916 : m_rhs(xpr.rhs()),
502 : 1916 : m_lhsImpl(m_lhs), // FIXME the creation of the evaluator objects should result in a no-op, but check that!
503 : 1916 : m_rhsImpl(m_rhs), // Moreover, they are only useful for the packet path, so we could completely disable them when not needed,
504 : : // or perhaps declare them on the fly on the packet method... We have experiment to check what's best.
505 : 3832 : m_innerDim(xpr.lhs().cols())
506 : : {
507 : : EIGEN_INTERNAL_CHECK_COST_VALUE(NumTraits<Scalar>::MulCost);
508 : : EIGEN_INTERNAL_CHECK_COST_VALUE(NumTraits<Scalar>::AddCost);
509 : : EIGEN_INTERNAL_CHECK_COST_VALUE(CoeffReadCost);
510 : : #if 0
511 : : std::cerr << "LhsOuterStrideBytes= " << LhsOuterStrideBytes << "\n";
512 : : std::cerr << "RhsOuterStrideBytes= " << RhsOuterStrideBytes << "\n";
513 : : std::cerr << "LhsAlignment= " << LhsAlignment << "\n";
514 : : std::cerr << "RhsAlignment= " << RhsAlignment << "\n";
515 : : std::cerr << "CanVectorizeLhs= " << CanVectorizeLhs << "\n";
516 : : std::cerr << "CanVectorizeRhs= " << CanVectorizeRhs << "\n";
517 : : std::cerr << "CanVectorizeInner= " << CanVectorizeInner << "\n";
518 : : std::cerr << "EvalToRowMajor= " << EvalToRowMajor << "\n";
519 : : std::cerr << "Alignment= " << Alignment << "\n";
520 : : std::cerr << "Flags= " << Flags << "\n";
521 : : #endif
522 : 1916 : }
523 : :
524 : : // Everything below here is taken from CoeffBasedProduct.h
525 : :
526 : : typedef typename internal::nested_eval<Lhs,Rhs::ColsAtCompileTime>::type LhsNested;
527 : : typedef typename internal::nested_eval<Rhs,Lhs::RowsAtCompileTime>::type RhsNested;
528 : :
529 : : typedef typename internal::remove_all<LhsNested>::type LhsNestedCleaned;
530 : : typedef typename internal::remove_all<RhsNested>::type RhsNestedCleaned;
531 : :
532 : : typedef evaluator<LhsNestedCleaned> LhsEtorType;
533 : : typedef evaluator<RhsNestedCleaned> RhsEtorType;
534 : :
535 : : enum {
536 : : RowsAtCompileTime = LhsNestedCleaned::RowsAtCompileTime,
537 : : ColsAtCompileTime = RhsNestedCleaned::ColsAtCompileTime,
538 : : InnerSize = EIGEN_SIZE_MIN_PREFER_FIXED(LhsNestedCleaned::ColsAtCompileTime, RhsNestedCleaned::RowsAtCompileTime),
539 : : MaxRowsAtCompileTime = LhsNestedCleaned::MaxRowsAtCompileTime,
540 : : MaxColsAtCompileTime = RhsNestedCleaned::MaxColsAtCompileTime
541 : : };
542 : :
543 : : typedef typename find_best_packet<Scalar,RowsAtCompileTime>::type LhsVecPacketType;
544 : : typedef typename find_best_packet<Scalar,ColsAtCompileTime>::type RhsVecPacketType;
545 : :
546 : : enum {
547 : :
548 : : LhsCoeffReadCost = LhsEtorType::CoeffReadCost,
549 : : RhsCoeffReadCost = RhsEtorType::CoeffReadCost,
550 : : CoeffReadCost = InnerSize==0 ? NumTraits<Scalar>::ReadCost
551 : : : InnerSize == Dynamic ? HugeCost
552 : : : InnerSize * (NumTraits<Scalar>::MulCost + int(LhsCoeffReadCost) + int(RhsCoeffReadCost))
553 : : + (InnerSize - 1) * NumTraits<Scalar>::AddCost,
554 : :
555 : : Unroll = CoeffReadCost <= EIGEN_UNROLLING_LIMIT,
556 : :
557 : : LhsFlags = LhsEtorType::Flags,
558 : : RhsFlags = RhsEtorType::Flags,
559 : :
560 : : LhsRowMajor = LhsFlags & RowMajorBit,
561 : : RhsRowMajor = RhsFlags & RowMajorBit,
562 : :
563 : : LhsVecPacketSize = unpacket_traits<LhsVecPacketType>::size,
564 : : RhsVecPacketSize = unpacket_traits<RhsVecPacketType>::size,
565 : :
566 : : // Here, we don't care about alignment larger than the usable packet size.
567 : : LhsAlignment = EIGEN_PLAIN_ENUM_MIN(LhsEtorType::Alignment,LhsVecPacketSize*int(sizeof(typename LhsNestedCleaned::Scalar))),
568 : : RhsAlignment = EIGEN_PLAIN_ENUM_MIN(RhsEtorType::Alignment,RhsVecPacketSize*int(sizeof(typename RhsNestedCleaned::Scalar))),
569 : :
570 : : SameType = is_same<typename LhsNestedCleaned::Scalar,typename RhsNestedCleaned::Scalar>::value,
571 : :
572 : : CanVectorizeRhs = bool(RhsRowMajor) && (RhsFlags & PacketAccessBit) && (ColsAtCompileTime!=1),
573 : : CanVectorizeLhs = (!LhsRowMajor) && (LhsFlags & PacketAccessBit) && (RowsAtCompileTime!=1),
574 : :
575 : : EvalToRowMajor = (MaxRowsAtCompileTime==1&&MaxColsAtCompileTime!=1) ? 1
576 : : : (MaxColsAtCompileTime==1&&MaxRowsAtCompileTime!=1) ? 0
577 : : : (bool(RhsRowMajor) && !CanVectorizeLhs),
578 : :
579 : : Flags = ((int(LhsFlags) | int(RhsFlags)) & HereditaryBits & ~RowMajorBit)
580 : : | (EvalToRowMajor ? RowMajorBit : 0)
581 : : // TODO enable vectorization for mixed types
582 : : | (SameType && (CanVectorizeLhs || CanVectorizeRhs) ? PacketAccessBit : 0)
583 : : | (XprType::IsVectorAtCompileTime ? LinearAccessBit : 0),
584 : :
585 : : LhsOuterStrideBytes = int(LhsNestedCleaned::OuterStrideAtCompileTime) * int(sizeof(typename LhsNestedCleaned::Scalar)),
586 : : RhsOuterStrideBytes = int(RhsNestedCleaned::OuterStrideAtCompileTime) * int(sizeof(typename RhsNestedCleaned::Scalar)),
587 : :
588 : : Alignment = bool(CanVectorizeLhs) ? (LhsOuterStrideBytes<=0 || (int(LhsOuterStrideBytes) % EIGEN_PLAIN_ENUM_MAX(1,LhsAlignment))!=0 ? 0 : LhsAlignment)
589 : : : bool(CanVectorizeRhs) ? (RhsOuterStrideBytes<=0 || (int(RhsOuterStrideBytes) % EIGEN_PLAIN_ENUM_MAX(1,RhsAlignment))!=0 ? 0 : RhsAlignment)
590 : : : 0,
591 : :
592 : : /* CanVectorizeInner deserves special explanation. It does not affect the product flags. It is not used outside
593 : : * of Product. If the Product itself is not a packet-access expression, there is still a chance that the inner
594 : : * loop of the product might be vectorized. This is the meaning of CanVectorizeInner. Since it doesn't affect
595 : : * the Flags, it is safe to make this value depend on ActualPacketAccessBit, that doesn't affect the ABI.
596 : : */
597 : : CanVectorizeInner = SameType
598 : : && LhsRowMajor
599 : : && (!RhsRowMajor)
600 : : && (int(LhsFlags) & int(RhsFlags) & ActualPacketAccessBit)
601 : : && (int(InnerSize) % packet_traits<Scalar>::size == 0)
602 : : };
603 : :
604 : 3130 : EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const CoeffReturnType coeff(Index row, Index col) const
605 : : {
606 : 3130 : return (m_lhs.row(row).transpose().cwiseProduct( m_rhs.col(col) )).sum();
607 : : }
608 : :
609 : : /* Allow index-based non-packet access. It is impossible though to allow index-based packed access,
610 : : * which is why we don't set the LinearAccessBit.
611 : : * TODO: this seems possible when the result is a vector
612 : : */
613 : : EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE
614 : : const CoeffReturnType coeff(Index index) const
615 : : {
616 : : const Index row = (RowsAtCompileTime == 1 || MaxRowsAtCompileTime==1) ? 0 : index;
617 : : const Index col = (RowsAtCompileTime == 1 || MaxRowsAtCompileTime==1) ? index : 0;
618 : : return (m_lhs.row(row).transpose().cwiseProduct( m_rhs.col(col) )).sum();
619 : : }
620 : :
621 : : template<int LoadMode, typename PacketType>
622 : : EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE
623 : 2068 : const PacketType packet(Index row, Index col) const
624 : : {
625 : : PacketType res;
626 : : typedef etor_product_packet_impl<bool(int(Flags)&RowMajorBit) ? RowMajor : ColMajor,
627 : : Unroll ? int(InnerSize) : Dynamic,
628 : : LhsEtorType, RhsEtorType, PacketType, LoadMode> PacketImpl;
629 : 2068 : PacketImpl::run(row, col, m_lhsImpl, m_rhsImpl, m_innerDim, res);
630 : 2068 : return res;
631 : : }
632 : :
633 : : template<int LoadMode, typename PacketType>
634 : : EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE
635 : : const PacketType packet(Index index) const
636 : : {
637 : : const Index row = (RowsAtCompileTime == 1 || MaxRowsAtCompileTime==1) ? 0 : index;
638 : : const Index col = (RowsAtCompileTime == 1 || MaxRowsAtCompileTime==1) ? index : 0;
639 : : return packet<LoadMode,PacketType>(row,col);
640 : : }
641 : :
642 : : protected:
643 : : typename internal::add_const_on_value_type<LhsNested>::type m_lhs;
644 : : typename internal::add_const_on_value_type<RhsNested>::type m_rhs;
645 : :
646 : : LhsEtorType m_lhsImpl;
647 : : RhsEtorType m_rhsImpl;
648 : :
649 : : // TODO: Get rid of m_innerDim if known at compile time
650 : : Index m_innerDim;
651 : : };
652 : :
653 : : template<typename Lhs, typename Rhs>
654 : : struct product_evaluator<Product<Lhs, Rhs, DefaultProduct>, LazyCoeffBasedProductMode, DenseShape, DenseShape>
655 : : : product_evaluator<Product<Lhs, Rhs, LazyProduct>, CoeffBasedProductMode, DenseShape, DenseShape>
656 : : {
657 : : typedef Product<Lhs, Rhs, DefaultProduct> XprType;
658 : : typedef Product<Lhs, Rhs, LazyProduct> BaseProduct;
659 : : typedef product_evaluator<BaseProduct, CoeffBasedProductMode, DenseShape, DenseShape> Base;
660 : : enum {
661 : : Flags = Base::Flags | EvalBeforeNestingBit
662 : : };
663 : : EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE
664 : 152 : explicit product_evaluator(const XprType& xpr)
665 : 152 : : Base(BaseProduct(xpr.lhs(),xpr.rhs()))
666 : 152 : {}
667 : : };
668 : :
669 : : /****************************************
670 : : *** Coeff based product, Packet path ***
671 : : ****************************************/
672 : :
673 : : template<int UnrollingIndex, typename Lhs, typename Rhs, typename Packet, int LoadMode>
674 : : struct etor_product_packet_impl<RowMajor, UnrollingIndex, Lhs, Rhs, Packet, LoadMode>
675 : : {
676 : : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void run(Index row, Index col, const Lhs& lhs, const Rhs& rhs, Index innerDim, Packet &res)
677 : : {
678 : : etor_product_packet_impl<RowMajor, UnrollingIndex-1, Lhs, Rhs, Packet, LoadMode>::run(row, col, lhs, rhs, innerDim, res);
679 : : res = pmadd(pset1<Packet>(lhs.coeff(row, Index(UnrollingIndex-1))), rhs.template packet<LoadMode,Packet>(Index(UnrollingIndex-1), col), res);
680 : : }
681 : : };
682 : :
683 : : template<int UnrollingIndex, typename Lhs, typename Rhs, typename Packet, int LoadMode>
684 : : struct etor_product_packet_impl<ColMajor, UnrollingIndex, Lhs, Rhs, Packet, LoadMode>
685 : : {
686 : 3830 : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void run(Index row, Index col, const Lhs& lhs, const Rhs& rhs, Index innerDim, Packet &res)
687 : : {
688 : 3830 : etor_product_packet_impl<ColMajor, UnrollingIndex-1, Lhs, Rhs, Packet, LoadMode>::run(row, col, lhs, rhs, innerDim, res);
689 : 3830 : res = pmadd(lhs.template packet<LoadMode,Packet>(row, Index(UnrollingIndex-1)), pset1<Packet>(rhs.coeff(Index(UnrollingIndex-1), col)), res);
690 : 3830 : }
691 : : };
692 : :
693 : : template<typename Lhs, typename Rhs, typename Packet, int LoadMode>
694 : : struct etor_product_packet_impl<RowMajor, 1, Lhs, Rhs, Packet, LoadMode>
695 : : {
696 : : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void run(Index row, Index col, const Lhs& lhs, const Rhs& rhs, Index /*innerDim*/, Packet &res)
697 : : {
698 : : res = pmul(pset1<Packet>(lhs.coeff(row, Index(0))),rhs.template packet<LoadMode,Packet>(Index(0), col));
699 : : }
700 : : };
701 : :
702 : : template<typename Lhs, typename Rhs, typename Packet, int LoadMode>
703 : : struct etor_product_packet_impl<ColMajor, 1, Lhs, Rhs, Packet, LoadMode>
704 : : {
705 : 2068 : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void run(Index row, Index col, const Lhs& lhs, const Rhs& rhs, Index /*innerDim*/, Packet &res)
706 : : {
707 : 2068 : res = pmul(lhs.template packet<LoadMode,Packet>(row, Index(0)), pset1<Packet>(rhs.coeff(Index(0), col)));
708 : 2068 : }
709 : : };
710 : :
711 : : template<typename Lhs, typename Rhs, typename Packet, int LoadMode>
712 : : struct etor_product_packet_impl<RowMajor, 0, Lhs, Rhs, Packet, LoadMode>
713 : : {
714 : : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void run(Index /*row*/, Index /*col*/, const Lhs& /*lhs*/, const Rhs& /*rhs*/, Index /*innerDim*/, Packet &res)
715 : : {
716 : : res = pset1<Packet>(typename unpacket_traits<Packet>::type(0));
717 : : }
718 : : };
719 : :
720 : : template<typename Lhs, typename Rhs, typename Packet, int LoadMode>
721 : : struct etor_product_packet_impl<ColMajor, 0, Lhs, Rhs, Packet, LoadMode>
722 : : {
723 : : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void run(Index /*row*/, Index /*col*/, const Lhs& /*lhs*/, const Rhs& /*rhs*/, Index /*innerDim*/, Packet &res)
724 : : {
725 : : res = pset1<Packet>(typename unpacket_traits<Packet>::type(0));
726 : : }
727 : : };
728 : :
729 : : template<typename Lhs, typename Rhs, typename Packet, int LoadMode>
730 : : struct etor_product_packet_impl<RowMajor, Dynamic, Lhs, Rhs, Packet, LoadMode>
731 : : {
732 : : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void run(Index row, Index col, const Lhs& lhs, const Rhs& rhs, Index innerDim, Packet& res)
733 : : {
734 : : res = pset1<Packet>(typename unpacket_traits<Packet>::type(0));
735 : : for(Index i = 0; i < innerDim; ++i)
736 : : res = pmadd(pset1<Packet>(lhs.coeff(row, i)), rhs.template packet<LoadMode,Packet>(i, col), res);
737 : : }
738 : : };
739 : :
740 : : template<typename Lhs, typename Rhs, typename Packet, int LoadMode>
741 : : struct etor_product_packet_impl<ColMajor, Dynamic, Lhs, Rhs, Packet, LoadMode>
742 : : {
743 : : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void run(Index row, Index col, const Lhs& lhs, const Rhs& rhs, Index innerDim, Packet& res)
744 : : {
745 : : res = pset1<Packet>(typename unpacket_traits<Packet>::type(0));
746 : : for(Index i = 0; i < innerDim; ++i)
747 : : res = pmadd(lhs.template packet<LoadMode,Packet>(row, i), pset1<Packet>(rhs.coeff(i, col)), res);
748 : : }
749 : : };
750 : :
751 : :
752 : : /***************************************************************************
753 : : * Triangular products
754 : : ***************************************************************************/
755 : : template<int Mode, bool LhsIsTriangular,
756 : : typename Lhs, bool LhsIsVector,
757 : : typename Rhs, bool RhsIsVector>
758 : : struct triangular_product_impl;
759 : :
760 : : template<typename Lhs, typename Rhs, int ProductTag>
761 : : struct generic_product_impl<Lhs,Rhs,TriangularShape,DenseShape,ProductTag>
762 : : : generic_product_impl_base<Lhs,Rhs,generic_product_impl<Lhs,Rhs,TriangularShape,DenseShape,ProductTag> >
763 : : {
764 : : typedef typename Product<Lhs,Rhs>::Scalar Scalar;
765 : :
766 : : template<typename Dest>
767 : : static void scaleAndAddTo(Dest& dst, const Lhs& lhs, const Rhs& rhs, const Scalar& alpha)
768 : : {
769 : : triangular_product_impl<Lhs::Mode,true,typename Lhs::MatrixType,false,Rhs, Rhs::ColsAtCompileTime==1>
770 : : ::run(dst, lhs.nestedExpression(), rhs, alpha);
771 : : }
772 : : };
773 : :
774 : : template<typename Lhs, typename Rhs, int ProductTag>
775 : : struct generic_product_impl<Lhs,Rhs,DenseShape,TriangularShape,ProductTag>
776 : : : generic_product_impl_base<Lhs,Rhs,generic_product_impl<Lhs,Rhs,DenseShape,TriangularShape,ProductTag> >
777 : : {
778 : : typedef typename Product<Lhs,Rhs>::Scalar Scalar;
779 : :
780 : : template<typename Dest>
781 : : static void scaleAndAddTo(Dest& dst, const Lhs& lhs, const Rhs& rhs, const Scalar& alpha)
782 : : {
783 : : triangular_product_impl<Rhs::Mode,false,Lhs,Lhs::RowsAtCompileTime==1, typename Rhs::MatrixType, false>::run(dst, lhs, rhs.nestedExpression(), alpha);
784 : : }
785 : : };
786 : :
787 : :
788 : : /***************************************************************************
789 : : * SelfAdjoint products
790 : : ***************************************************************************/
791 : : template <typename Lhs, int LhsMode, bool LhsIsVector,
792 : : typename Rhs, int RhsMode, bool RhsIsVector>
793 : : struct selfadjoint_product_impl;
794 : :
795 : : template<typename Lhs, typename Rhs, int ProductTag>
796 : : struct generic_product_impl<Lhs,Rhs,SelfAdjointShape,DenseShape,ProductTag>
797 : : : generic_product_impl_base<Lhs,Rhs,generic_product_impl<Lhs,Rhs,SelfAdjointShape,DenseShape,ProductTag> >
798 : : {
799 : : typedef typename Product<Lhs,Rhs>::Scalar Scalar;
800 : :
801 : : template<typename Dest>
802 : : static EIGEN_DEVICE_FUNC
803 : : void scaleAndAddTo(Dest& dst, const Lhs& lhs, const Rhs& rhs, const Scalar& alpha)
804 : : {
805 : : selfadjoint_product_impl<typename Lhs::MatrixType,Lhs::Mode,false,Rhs,0,Rhs::IsVectorAtCompileTime>::run(dst, lhs.nestedExpression(), rhs, alpha);
806 : : }
807 : : };
808 : :
809 : : template<typename Lhs, typename Rhs, int ProductTag>
810 : : struct generic_product_impl<Lhs,Rhs,DenseShape,SelfAdjointShape,ProductTag>
811 : : : generic_product_impl_base<Lhs,Rhs,generic_product_impl<Lhs,Rhs,DenseShape,SelfAdjointShape,ProductTag> >
812 : : {
813 : : typedef typename Product<Lhs,Rhs>::Scalar Scalar;
814 : :
815 : : template<typename Dest>
816 : : static void scaleAndAddTo(Dest& dst, const Lhs& lhs, const Rhs& rhs, const Scalar& alpha)
817 : : {
818 : : selfadjoint_product_impl<Lhs,0,Lhs::IsVectorAtCompileTime,typename Rhs::MatrixType,Rhs::Mode,false>::run(dst, lhs, rhs.nestedExpression(), alpha);
819 : : }
820 : : };
821 : :
822 : :
823 : : /***************************************************************************
824 : : * Diagonal products
825 : : ***************************************************************************/
826 : :
827 : : template<typename MatrixType, typename DiagonalType, typename Derived, int ProductOrder>
828 : : struct diagonal_product_evaluator_base
829 : : : evaluator_base<Derived>
830 : : {
831 : : typedef typename ScalarBinaryOpTraits<typename MatrixType::Scalar, typename DiagonalType::Scalar>::ReturnType Scalar;
832 : : public:
833 : : enum {
834 : : CoeffReadCost = int(NumTraits<Scalar>::MulCost) + int(evaluator<MatrixType>::CoeffReadCost) + int(evaluator<DiagonalType>::CoeffReadCost),
835 : :
836 : : MatrixFlags = evaluator<MatrixType>::Flags,
837 : : DiagFlags = evaluator<DiagonalType>::Flags,
838 : :
839 : : _StorageOrder = (Derived::MaxRowsAtCompileTime==1 && Derived::MaxColsAtCompileTime!=1) ? RowMajor
840 : : : (Derived::MaxColsAtCompileTime==1 && Derived::MaxRowsAtCompileTime!=1) ? ColMajor
841 : : : MatrixFlags & RowMajorBit ? RowMajor : ColMajor,
842 : : _SameStorageOrder = _StorageOrder == (MatrixFlags & RowMajorBit ? RowMajor : ColMajor),
843 : :
844 : : _ScalarAccessOnDiag = !((int(_StorageOrder) == ColMajor && int(ProductOrder) == OnTheLeft)
845 : : ||(int(_StorageOrder) == RowMajor && int(ProductOrder) == OnTheRight)),
846 : : _SameTypes = is_same<typename MatrixType::Scalar, typename DiagonalType::Scalar>::value,
847 : : // FIXME currently we need same types, but in the future the next rule should be the one
848 : : //_Vectorizable = bool(int(MatrixFlags)&PacketAccessBit) && ((!_PacketOnDiag) || (_SameTypes && bool(int(DiagFlags)&PacketAccessBit))),
849 : : _Vectorizable = bool(int(MatrixFlags)&PacketAccessBit)
850 : : && _SameTypes
851 : : && (_SameStorageOrder || (MatrixFlags&LinearAccessBit)==LinearAccessBit)
852 : : && (_ScalarAccessOnDiag || (bool(int(DiagFlags)&PacketAccessBit))),
853 : : _LinearAccessMask = (MatrixType::RowsAtCompileTime==1 || MatrixType::ColsAtCompileTime==1) ? LinearAccessBit : 0,
854 : : Flags = ((HereditaryBits|_LinearAccessMask) & (unsigned int)(MatrixFlags)) | (_Vectorizable ? PacketAccessBit : 0),
855 : : Alignment = evaluator<MatrixType>::Alignment,
856 : :
857 : : AsScalarProduct = (DiagonalType::SizeAtCompileTime==1)
858 : : || (DiagonalType::SizeAtCompileTime==Dynamic && MatrixType::RowsAtCompileTime==1 && ProductOrder==OnTheLeft)
859 : : || (DiagonalType::SizeAtCompileTime==Dynamic && MatrixType::ColsAtCompileTime==1 && ProductOrder==OnTheRight)
860 : : };
861 : :
862 : : EIGEN_DEVICE_FUNC diagonal_product_evaluator_base(const MatrixType &mat, const DiagonalType &diag)
863 : : : m_diagImpl(diag), m_matImpl(mat)
864 : : {
865 : : EIGEN_INTERNAL_CHECK_COST_VALUE(NumTraits<Scalar>::MulCost);
866 : : EIGEN_INTERNAL_CHECK_COST_VALUE(CoeffReadCost);
867 : : }
868 : :
869 : : EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const Scalar coeff(Index idx) const
870 : : {
871 : : if(AsScalarProduct)
872 : : return m_diagImpl.coeff(0) * m_matImpl.coeff(idx);
873 : : else
874 : : return m_diagImpl.coeff(idx) * m_matImpl.coeff(idx);
875 : : }
876 : :
877 : : protected:
878 : : template<int LoadMode,typename PacketType>
879 : : EIGEN_STRONG_INLINE PacketType packet_impl(Index row, Index col, Index id, internal::true_type) const
880 : : {
881 : : return internal::pmul(m_matImpl.template packet<LoadMode,PacketType>(row, col),
882 : : internal::pset1<PacketType>(m_diagImpl.coeff(id)));
883 : : }
884 : :
885 : : template<int LoadMode,typename PacketType>
886 : : EIGEN_STRONG_INLINE PacketType packet_impl(Index row, Index col, Index id, internal::false_type) const
887 : : {
888 : : enum {
889 : : InnerSize = (MatrixType::Flags & RowMajorBit) ? MatrixType::ColsAtCompileTime : MatrixType::RowsAtCompileTime,
890 : : DiagonalPacketLoadMode = EIGEN_PLAIN_ENUM_MIN(LoadMode,((InnerSize%16) == 0) ? int(Aligned16) : int(evaluator<DiagonalType>::Alignment)) // FIXME hardcoded 16!!
891 : : };
892 : : return internal::pmul(m_matImpl.template packet<LoadMode,PacketType>(row, col),
893 : : m_diagImpl.template packet<DiagonalPacketLoadMode,PacketType>(id));
894 : : }
895 : :
896 : : evaluator<DiagonalType> m_diagImpl;
897 : : evaluator<MatrixType> m_matImpl;
898 : : };
899 : :
900 : : // diagonal * dense
901 : : template<typename Lhs, typename Rhs, int ProductKind, int ProductTag>
902 : : struct product_evaluator<Product<Lhs, Rhs, ProductKind>, ProductTag, DiagonalShape, DenseShape>
903 : : : diagonal_product_evaluator_base<Rhs, typename Lhs::DiagonalVectorType, Product<Lhs, Rhs, LazyProduct>, OnTheLeft>
904 : : {
905 : : typedef diagonal_product_evaluator_base<Rhs, typename Lhs::DiagonalVectorType, Product<Lhs, Rhs, LazyProduct>, OnTheLeft> Base;
906 : : using Base::m_diagImpl;
907 : : using Base::m_matImpl;
908 : : using Base::coeff;
909 : : typedef typename Base::Scalar Scalar;
910 : :
911 : : typedef Product<Lhs, Rhs, ProductKind> XprType;
912 : : typedef typename XprType::PlainObject PlainObject;
913 : : typedef typename Lhs::DiagonalVectorType DiagonalType;
914 : :
915 : :
916 : : enum { StorageOrder = Base::_StorageOrder };
917 : :
918 : : EIGEN_DEVICE_FUNC explicit product_evaluator(const XprType& xpr)
919 : : : Base(xpr.rhs(), xpr.lhs().diagonal())
920 : : {
921 : : }
922 : :
923 : : EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const Scalar coeff(Index row, Index col) const
924 : : {
925 : : return m_diagImpl.coeff(row) * m_matImpl.coeff(row, col);
926 : : }
927 : :
928 : : #ifndef EIGEN_GPUCC
929 : : template<int LoadMode,typename PacketType>
930 : : EIGEN_STRONG_INLINE PacketType packet(Index row, Index col) const
931 : : {
932 : : // FIXME: NVCC used to complain about the template keyword, but we have to check whether this is still the case.
933 : : // See also similar calls below.
934 : : return this->template packet_impl<LoadMode,PacketType>(row,col, row,
935 : : typename internal::conditional<int(StorageOrder)==RowMajor, internal::true_type, internal::false_type>::type());
936 : : }
937 : :
938 : : template<int LoadMode,typename PacketType>
939 : : EIGEN_STRONG_INLINE PacketType packet(Index idx) const
940 : : {
941 : : return packet<LoadMode,PacketType>(int(StorageOrder)==ColMajor?idx:0,int(StorageOrder)==ColMajor?0:idx);
942 : : }
943 : : #endif
944 : : };
945 : :
946 : : // dense * diagonal
947 : : template<typename Lhs, typename Rhs, int ProductKind, int ProductTag>
948 : : struct product_evaluator<Product<Lhs, Rhs, ProductKind>, ProductTag, DenseShape, DiagonalShape>
949 : : : diagonal_product_evaluator_base<Lhs, typename Rhs::DiagonalVectorType, Product<Lhs, Rhs, LazyProduct>, OnTheRight>
950 : : {
951 : : typedef diagonal_product_evaluator_base<Lhs, typename Rhs::DiagonalVectorType, Product<Lhs, Rhs, LazyProduct>, OnTheRight> Base;
952 : : using Base::m_diagImpl;
953 : : using Base::m_matImpl;
954 : : using Base::coeff;
955 : : typedef typename Base::Scalar Scalar;
956 : :
957 : : typedef Product<Lhs, Rhs, ProductKind> XprType;
958 : : typedef typename XprType::PlainObject PlainObject;
959 : :
960 : : enum { StorageOrder = Base::_StorageOrder };
961 : :
962 : : EIGEN_DEVICE_FUNC explicit product_evaluator(const XprType& xpr)
963 : : : Base(xpr.lhs(), xpr.rhs().diagonal())
964 : : {
965 : : }
966 : :
967 : : EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const Scalar coeff(Index row, Index col) const
968 : : {
969 : : return m_matImpl.coeff(row, col) * m_diagImpl.coeff(col);
970 : : }
971 : :
972 : : #ifndef EIGEN_GPUCC
973 : : template<int LoadMode,typename PacketType>
974 : : EIGEN_STRONG_INLINE PacketType packet(Index row, Index col) const
975 : : {
976 : : return this->template packet_impl<LoadMode,PacketType>(row,col, col,
977 : : typename internal::conditional<int(StorageOrder)==ColMajor, internal::true_type, internal::false_type>::type());
978 : : }
979 : :
980 : : template<int LoadMode,typename PacketType>
981 : : EIGEN_STRONG_INLINE PacketType packet(Index idx) const
982 : : {
983 : : return packet<LoadMode,PacketType>(int(StorageOrder)==ColMajor?idx:0,int(StorageOrder)==ColMajor?0:idx);
984 : : }
985 : : #endif
986 : : };
987 : :
988 : : /***************************************************************************
989 : : * Products with permutation matrices
990 : : ***************************************************************************/
991 : :
992 : : /** \internal
993 : : * \class permutation_matrix_product
994 : : * Internal helper class implementing the product between a permutation matrix and a matrix.
995 : : * This class is specialized for DenseShape below and for SparseShape in SparseCore/SparsePermutation.h
996 : : */
997 : : template<typename ExpressionType, int Side, bool Transposed, typename ExpressionShape>
998 : : struct permutation_matrix_product;
999 : :
1000 : : template<typename ExpressionType, int Side, bool Transposed>
1001 : : struct permutation_matrix_product<ExpressionType, Side, Transposed, DenseShape>
1002 : : {
1003 : : typedef typename nested_eval<ExpressionType, 1>::type MatrixType;
1004 : : typedef typename remove_all<MatrixType>::type MatrixTypeCleaned;
1005 : :
1006 : : template<typename Dest, typename PermutationType>
1007 : : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void run(Dest& dst, const PermutationType& perm, const ExpressionType& xpr)
1008 : : {
1009 : : MatrixType mat(xpr);
1010 : : const Index n = Side==OnTheLeft ? mat.rows() : mat.cols();
1011 : : // FIXME we need an is_same for expression that is not sensitive to constness. For instance
1012 : : // is_same_xpr<Block<const Matrix>, Block<Matrix> >::value should be true.
1013 : : //if(is_same<MatrixTypeCleaned,Dest>::value && extract_data(dst) == extract_data(mat))
1014 : : if(is_same_dense(dst, mat))
1015 : : {
1016 : : // apply the permutation inplace
1017 : : Matrix<bool,PermutationType::RowsAtCompileTime,1,0,PermutationType::MaxRowsAtCompileTime> mask(perm.size());
1018 : : mask.fill(false);
1019 : : Index r = 0;
1020 : : while(r < perm.size())
1021 : : {
1022 : : // search for the next seed
1023 : : while(r<perm.size() && mask[r]) r++;
1024 : : if(r>=perm.size())
1025 : : break;
1026 : : // we got one, let's follow it until we are back to the seed
1027 : : Index k0 = r++;
1028 : : Index kPrev = k0;
1029 : : mask.coeffRef(k0) = true;
1030 : : for(Index k=perm.indices().coeff(k0); k!=k0; k=perm.indices().coeff(k))
1031 : : {
1032 : : Block<Dest, Side==OnTheLeft ? 1 : Dest::RowsAtCompileTime, Side==OnTheRight ? 1 : Dest::ColsAtCompileTime>(dst, k)
1033 : : .swap(Block<Dest, Side==OnTheLeft ? 1 : Dest::RowsAtCompileTime, Side==OnTheRight ? 1 : Dest::ColsAtCompileTime>
1034 : : (dst,((Side==OnTheLeft) ^ Transposed) ? k0 : kPrev));
1035 : :
1036 : : mask.coeffRef(k) = true;
1037 : : kPrev = k;
1038 : : }
1039 : : }
1040 : : }
1041 : : else
1042 : : {
1043 : : for(Index i = 0; i < n; ++i)
1044 : : {
1045 : : Block<Dest, Side==OnTheLeft ? 1 : Dest::RowsAtCompileTime, Side==OnTheRight ? 1 : Dest::ColsAtCompileTime>
1046 : : (dst, ((Side==OnTheLeft) ^ Transposed) ? perm.indices().coeff(i) : i)
1047 : :
1048 : : =
1049 : :
1050 : : Block<const MatrixTypeCleaned,Side==OnTheLeft ? 1 : MatrixTypeCleaned::RowsAtCompileTime,Side==OnTheRight ? 1 : MatrixTypeCleaned::ColsAtCompileTime>
1051 : : (mat, ((Side==OnTheRight) ^ Transposed) ? perm.indices().coeff(i) : i);
1052 : : }
1053 : : }
1054 : : }
1055 : : };
1056 : :
1057 : : template<typename Lhs, typename Rhs, int ProductTag, typename MatrixShape>
1058 : : struct generic_product_impl<Lhs, Rhs, PermutationShape, MatrixShape, ProductTag>
1059 : : {
1060 : : template<typename Dest>
1061 : : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void evalTo(Dest& dst, const Lhs& lhs, const Rhs& rhs)
1062 : : {
1063 : : permutation_matrix_product<Rhs, OnTheLeft, false, MatrixShape>::run(dst, lhs, rhs);
1064 : : }
1065 : : };
1066 : :
1067 : : template<typename Lhs, typename Rhs, int ProductTag, typename MatrixShape>
1068 : : struct generic_product_impl<Lhs, Rhs, MatrixShape, PermutationShape, ProductTag>
1069 : : {
1070 : : template<typename Dest>
1071 : : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void evalTo(Dest& dst, const Lhs& lhs, const Rhs& rhs)
1072 : : {
1073 : : permutation_matrix_product<Lhs, OnTheRight, false, MatrixShape>::run(dst, rhs, lhs);
1074 : : }
1075 : : };
1076 : :
1077 : : template<typename Lhs, typename Rhs, int ProductTag, typename MatrixShape>
1078 : : struct generic_product_impl<Inverse<Lhs>, Rhs, PermutationShape, MatrixShape, ProductTag>
1079 : : {
1080 : : template<typename Dest>
1081 : : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void evalTo(Dest& dst, const Inverse<Lhs>& lhs, const Rhs& rhs)
1082 : : {
1083 : : permutation_matrix_product<Rhs, OnTheLeft, true, MatrixShape>::run(dst, lhs.nestedExpression(), rhs);
1084 : : }
1085 : : };
1086 : :
1087 : : template<typename Lhs, typename Rhs, int ProductTag, typename MatrixShape>
1088 : : struct generic_product_impl<Lhs, Inverse<Rhs>, MatrixShape, PermutationShape, ProductTag>
1089 : : {
1090 : : template<typename Dest>
1091 : : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void evalTo(Dest& dst, const Lhs& lhs, const Inverse<Rhs>& rhs)
1092 : : {
1093 : : permutation_matrix_product<Lhs, OnTheRight, true, MatrixShape>::run(dst, rhs.nestedExpression(), lhs);
1094 : : }
1095 : : };
1096 : :
1097 : :
1098 : : /***************************************************************************
1099 : : * Products with transpositions matrices
1100 : : ***************************************************************************/
1101 : :
1102 : : // FIXME could we unify Transpositions and Permutation into a single "shape"??
1103 : :
1104 : : /** \internal
1105 : : * \class transposition_matrix_product
1106 : : * Internal helper class implementing the product between a permutation matrix and a matrix.
1107 : : */
1108 : : template<typename ExpressionType, int Side, bool Transposed, typename ExpressionShape>
1109 : : struct transposition_matrix_product
1110 : : {
1111 : : typedef typename nested_eval<ExpressionType, 1>::type MatrixType;
1112 : : typedef typename remove_all<MatrixType>::type MatrixTypeCleaned;
1113 : :
1114 : : template<typename Dest, typename TranspositionType>
1115 : : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void run(Dest& dst, const TranspositionType& tr, const ExpressionType& xpr)
1116 : : {
1117 : : MatrixType mat(xpr);
1118 : : typedef typename TranspositionType::StorageIndex StorageIndex;
1119 : : const Index size = tr.size();
1120 : : StorageIndex j = 0;
1121 : :
1122 : : if(!is_same_dense(dst,mat))
1123 : : dst = mat;
1124 : :
1125 : : for(Index k=(Transposed?size-1:0) ; Transposed?k>=0:k<size ; Transposed?--k:++k)
1126 : : if(Index(j=tr.coeff(k))!=k)
1127 : : {
1128 : : if(Side==OnTheLeft) dst.row(k).swap(dst.row(j));
1129 : : else if(Side==OnTheRight) dst.col(k).swap(dst.col(j));
1130 : : }
1131 : : }
1132 : : };
1133 : :
1134 : : template<typename Lhs, typename Rhs, int ProductTag, typename MatrixShape>
1135 : : struct generic_product_impl<Lhs, Rhs, TranspositionsShape, MatrixShape, ProductTag>
1136 : : {
1137 : : template<typename Dest>
1138 : : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void evalTo(Dest& dst, const Lhs& lhs, const Rhs& rhs)
1139 : : {
1140 : : transposition_matrix_product<Rhs, OnTheLeft, false, MatrixShape>::run(dst, lhs, rhs);
1141 : : }
1142 : : };
1143 : :
1144 : : template<typename Lhs, typename Rhs, int ProductTag, typename MatrixShape>
1145 : : struct generic_product_impl<Lhs, Rhs, MatrixShape, TranspositionsShape, ProductTag>
1146 : : {
1147 : : template<typename Dest>
1148 : : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void evalTo(Dest& dst, const Lhs& lhs, const Rhs& rhs)
1149 : : {
1150 : : transposition_matrix_product<Lhs, OnTheRight, false, MatrixShape>::run(dst, rhs, lhs);
1151 : : }
1152 : : };
1153 : :
1154 : :
1155 : : template<typename Lhs, typename Rhs, int ProductTag, typename MatrixShape>
1156 : : struct generic_product_impl<Transpose<Lhs>, Rhs, TranspositionsShape, MatrixShape, ProductTag>
1157 : : {
1158 : : template<typename Dest>
1159 : : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void evalTo(Dest& dst, const Transpose<Lhs>& lhs, const Rhs& rhs)
1160 : : {
1161 : : transposition_matrix_product<Rhs, OnTheLeft, true, MatrixShape>::run(dst, lhs.nestedExpression(), rhs);
1162 : : }
1163 : : };
1164 : :
1165 : : template<typename Lhs, typename Rhs, int ProductTag, typename MatrixShape>
1166 : : struct generic_product_impl<Lhs, Transpose<Rhs>, MatrixShape, TranspositionsShape, ProductTag>
1167 : : {
1168 : : template<typename Dest>
1169 : : static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void evalTo(Dest& dst, const Lhs& lhs, const Transpose<Rhs>& rhs)
1170 : : {
1171 : : transposition_matrix_product<Lhs, OnTheRight, true, MatrixShape>::run(dst, rhs.nestedExpression(), lhs);
1172 : : }
1173 : : };
1174 : :
1175 : : } // end namespace internal
1176 : :
1177 : : } // end namespace Eigen
1178 : :
1179 : : #endif // EIGEN_PRODUCT_EVALUATORS_H
|