1 // This file is part of Eigen, a lightweight C++ template library 2 // for linear algebra. 3 // 4 // Copyright (C) 2008-2019 Gael Guennebaud <[email protected]> 5 // Copyright (C) 2006-2008 Benoit Jacob <[email protected]> 6 // 7 // This Source Code Form is subject to the terms of the Mozilla 8 // Public License v. 2.0. If a copy of the MPL was not distributed 9 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/. 10 11 #ifndef EIGEN_PARTIAL_REDUX_H 12 #define EIGEN_PARTIAL_REDUX_H 13 14 namespace Eigen { 15 16 /** \class PartialReduxExpr 17 * \ingroup Core_Module 18 * 19 * \brief Generic expression of a partially reduxed matrix 20 * 21 * \tparam MatrixType the type of the matrix we are applying the redux operation 22 * \tparam MemberOp type of the member functor 23 * \tparam Direction indicates the direction of the redux (#Vertical or #Horizontal) 24 * 25 * This class represents an expression of a partial redux operator of a matrix. 26 * It is the return type of some VectorwiseOp functions, 27 * and most of the time this is the only way it is used. 28 * 29 * \sa class VectorwiseOp 30 */ 31 32 template< typename MatrixType, typename MemberOp, int Direction> 33 class PartialReduxExpr; 34 35 namespace internal { 36 template<typename MatrixType, typename MemberOp, int Direction> 37 struct traits<PartialReduxExpr<MatrixType, MemberOp, Direction> > 38 : traits<MatrixType> 39 { 40 typedef typename MemberOp::result_type Scalar; 41 typedef typename traits<MatrixType>::StorageKind StorageKind; 42 typedef typename traits<MatrixType>::XprKind XprKind; 43 typedef typename MatrixType::Scalar InputScalar; 44 enum { 45 RowsAtCompileTime = Direction==Vertical ? 1 : MatrixType::RowsAtCompileTime, 46 ColsAtCompileTime = Direction==Horizontal ? 1 : MatrixType::ColsAtCompileTime, 47 MaxRowsAtCompileTime = Direction==Vertical ? 1 : MatrixType::MaxRowsAtCompileTime, 48 MaxColsAtCompileTime = Direction==Horizontal ? 1 : MatrixType::MaxColsAtCompileTime, 49 Flags = RowsAtCompileTime == 1 ? RowMajorBit : 0, 50 TraversalSize = Direction==Vertical ? MatrixType::RowsAtCompileTime : MatrixType::ColsAtCompileTime 51 }; 52 }; 53 } 54 55 template< typename MatrixType, typename MemberOp, int Direction> 56 class PartialReduxExpr : public internal::dense_xpr_base< PartialReduxExpr<MatrixType, MemberOp, Direction> >::type, 57 internal::no_assignment_operator 58 { 59 public: 60 61 typedef typename internal::dense_xpr_base<PartialReduxExpr>::type Base; 62 EIGEN_DENSE_PUBLIC_INTERFACE(PartialReduxExpr) 63 64 EIGEN_DEVICE_FUNC 65 explicit PartialReduxExpr(const MatrixType& mat, const MemberOp& func = MemberOp()) 66 : m_matrix(mat), m_functor(func) {} 67 68 EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR 69 Index rows() const EIGEN_NOEXCEPT { return (Direction==Vertical ? 1 : m_matrix.rows()); } 70 EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR 71 Index cols() const EIGEN_NOEXCEPT { return (Direction==Horizontal ? 1 : m_matrix.cols()); } 72 73 EIGEN_DEVICE_FUNC 74 typename MatrixType::Nested nestedExpression() const { return m_matrix; } 75 76 EIGEN_DEVICE_FUNC 77 const MemberOp& functor() const { return m_functor; } 78 79 protected: 80 typename MatrixType::Nested m_matrix; 81 const MemberOp m_functor; 82 }; 83 84 template<typename A,typename B> struct partial_redux_dummy_func; 85 86 #define EIGEN_MAKE_PARTIAL_REDUX_FUNCTOR(MEMBER,COST,VECTORIZABLE,BINARYOP) \ 87 template <typename ResultType,typename Scalar> \ 88 struct member_##MEMBER { \ 89 EIGEN_EMPTY_STRUCT_CTOR(member_##MEMBER) \ 90 typedef ResultType result_type; \ 91 typedef BINARYOP<Scalar,Scalar> BinaryOp; \ 92 template<int Size> struct Cost { enum { value = COST }; }; \ 93 enum { Vectorizable = VECTORIZABLE }; \ 94 template<typename XprType> \ 95 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE \ 96 ResultType operator()(const XprType& mat) const \ 97 { return mat.MEMBER(); } \ 98 BinaryOp binaryFunc() const { return BinaryOp(); } \ 99 } 100 101 #define EIGEN_MEMBER_FUNCTOR(MEMBER,COST) \ 102 EIGEN_MAKE_PARTIAL_REDUX_FUNCTOR(MEMBER,COST,0,partial_redux_dummy_func) 103 104 namespace internal { 105 106 EIGEN_MEMBER_FUNCTOR(norm, (Size+5) * NumTraits<Scalar>::MulCost + (Size-1)*NumTraits<Scalar>::AddCost); 107 EIGEN_MEMBER_FUNCTOR(stableNorm, (Size+5) * NumTraits<Scalar>::MulCost + (Size-1)*NumTraits<Scalar>::AddCost); 108 EIGEN_MEMBER_FUNCTOR(blueNorm, (Size+5) * NumTraits<Scalar>::MulCost + (Size-1)*NumTraits<Scalar>::AddCost); 109 EIGEN_MEMBER_FUNCTOR(hypotNorm, (Size-1) * functor_traits<scalar_hypot_op<Scalar> >::Cost ); 110 EIGEN_MEMBER_FUNCTOR(all, (Size-1)*NumTraits<Scalar>::AddCost); 111 EIGEN_MEMBER_FUNCTOR(any, (Size-1)*NumTraits<Scalar>::AddCost); 112 EIGEN_MEMBER_FUNCTOR(count, (Size-1)*NumTraits<Scalar>::AddCost); 113 114 EIGEN_MAKE_PARTIAL_REDUX_FUNCTOR(sum, (Size-1)*NumTraits<Scalar>::AddCost, 1, internal::scalar_sum_op); 115 EIGEN_MAKE_PARTIAL_REDUX_FUNCTOR(minCoeff, (Size-1)*NumTraits<Scalar>::AddCost, 1, internal::scalar_min_op); 116 EIGEN_MAKE_PARTIAL_REDUX_FUNCTOR(maxCoeff, (Size-1)*NumTraits<Scalar>::AddCost, 1, internal::scalar_max_op); 117 EIGEN_MAKE_PARTIAL_REDUX_FUNCTOR(prod, (Size-1)*NumTraits<Scalar>::MulCost, 1, internal::scalar_product_op); 118 119 template <int p, typename ResultType,typename Scalar> 120 struct member_lpnorm { 121 typedef ResultType result_type; 122 enum { Vectorizable = 0 }; 123 template<int Size> struct Cost 124 { enum { value = (Size+5) * NumTraits<Scalar>::MulCost + (Size-1)*NumTraits<Scalar>::AddCost }; }; 125 EIGEN_DEVICE_FUNC member_lpnorm() {} 126 template<typename XprType> 127 EIGEN_DEVICE_FUNC inline ResultType operator()(const XprType& mat) const 128 { return mat.template lpNorm<p>(); } 129 }; 130 131 template <typename BinaryOpT, typename Scalar> 132 struct member_redux { 133 typedef BinaryOpT BinaryOp; 134 typedef typename result_of< 135 BinaryOp(const Scalar&,const Scalar&) 136 >::type result_type; 137 138 enum { Vectorizable = functor_traits<BinaryOp>::PacketAccess }; 139 template<int Size> struct Cost { enum { value = (Size-1) * functor_traits<BinaryOp>::Cost }; }; 140 EIGEN_DEVICE_FUNC explicit member_redux(const BinaryOp func) : m_functor(func) {} 141 template<typename Derived> 142 EIGEN_DEVICE_FUNC inline result_type operator()(const DenseBase<Derived>& mat) const 143 { return mat.redux(m_functor); } 144 const BinaryOp& binaryFunc() const { return m_functor; } 145 const BinaryOp m_functor; 146 }; 147 } 148 149 /** \class VectorwiseOp 150 * \ingroup Core_Module 151 * 152 * \brief Pseudo expression providing broadcasting and partial reduction operations 153 * 154 * \tparam ExpressionType the type of the object on which to do partial reductions 155 * \tparam Direction indicates whether to operate on columns (#Vertical) or rows (#Horizontal) 156 * 157 * This class represents a pseudo expression with broadcasting and partial reduction features. 158 * It is the return type of DenseBase::colwise() and DenseBase::rowwise() 159 * and most of the time this is the only way it is explicitly used. 160 * 161 * To understand the logic of rowwise/colwise expression, let's consider a generic case `A.colwise().foo()` 162 * where `foo` is any method of `VectorwiseOp`. This expression is equivalent to applying `foo()` to each 163 * column of `A` and then re-assemble the outputs in a matrix expression: 164 * \code [A.col(0).foo(), A.col(1).foo(), ..., A.col(A.cols()-1).foo()] \endcode 165 * 166 * Example: \include MatrixBase_colwise.cpp 167 * Output: \verbinclude MatrixBase_colwise.out 168 * 169 * The begin() and end() methods are obviously exceptions to the previous rule as they 170 * return STL-compatible begin/end iterators to the rows or columns of the nested expression. 171 * Typical use cases include for-range-loop and calls to STL algorithms: 172 * 173 * Example: \include MatrixBase_colwise_iterator_cxx11.cpp 174 * Output: \verbinclude MatrixBase_colwise_iterator_cxx11.out 175 * 176 * For a partial reduction on an empty input, some rules apply. 177 * For the sake of clarity, let's consider a vertical reduction: 178 * - If the number of columns is zero, then a 1x0 row-major vector expression is returned. 179 * - Otherwise, if the number of rows is zero, then 180 * - a row vector of zeros is returned for sum-like reductions (sum, squaredNorm, norm, etc.) 181 * - a row vector of ones is returned for a product reduction (e.g., <code>MatrixXd(n,0).colwise().prod()</code>) 182 * - an assert is triggered for all other reductions (minCoeff,maxCoeff,redux(bin_op)) 183 * 184 * \sa DenseBase::colwise(), DenseBase::rowwise(), class PartialReduxExpr 185 */ 186 template<typename ExpressionType, int Direction> class VectorwiseOp 187 { 188 public: 189 190 typedef typename ExpressionType::Scalar Scalar; 191 typedef typename ExpressionType::RealScalar RealScalar; 192 typedef Eigen::Index Index; ///< \deprecated since Eigen 3.3 193 typedef typename internal::ref_selector<ExpressionType>::non_const_type ExpressionTypeNested; 194 typedef typename internal::remove_all<ExpressionTypeNested>::type ExpressionTypeNestedCleaned; 195 196 template<template<typename OutScalar,typename InputScalar> class Functor, 197 typename ReturnScalar=Scalar> struct ReturnType 198 { 199 typedef PartialReduxExpr<ExpressionType, 200 Functor<ReturnScalar,Scalar>, 201 Direction 202 > Type; 203 }; 204 205 template<typename BinaryOp> struct ReduxReturnType 206 { 207 typedef PartialReduxExpr<ExpressionType, 208 internal::member_redux<BinaryOp,Scalar>, 209 Direction 210 > Type; 211 }; 212 213 enum { 214 isVertical = (Direction==Vertical) ? 1 : 0, 215 isHorizontal = (Direction==Horizontal) ? 1 : 0 216 }; 217 218 protected: 219 220 template<typename OtherDerived> struct ExtendedType { 221 typedef Replicate<OtherDerived, 222 isVertical ? 1 : ExpressionType::RowsAtCompileTime, 223 isHorizontal ? 1 : ExpressionType::ColsAtCompileTime> Type; 224 }; 225 226 /** \internal 227 * Replicates a vector to match the size of \c *this */ 228 template<typename OtherDerived> 229 EIGEN_DEVICE_FUNC 230 typename ExtendedType<OtherDerived>::Type 231 extendedTo(const DenseBase<OtherDerived>& other) const 232 { 233 EIGEN_STATIC_ASSERT(EIGEN_IMPLIES(isVertical, OtherDerived::MaxColsAtCompileTime==1), 234 YOU_PASSED_A_ROW_VECTOR_BUT_A_COLUMN_VECTOR_WAS_EXPECTED) 235 EIGEN_STATIC_ASSERT(EIGEN_IMPLIES(isHorizontal, OtherDerived::MaxRowsAtCompileTime==1), 236 YOU_PASSED_A_COLUMN_VECTOR_BUT_A_ROW_VECTOR_WAS_EXPECTED) 237 return typename ExtendedType<OtherDerived>::Type 238 (other.derived(), 239 isVertical ? 1 : m_matrix.rows(), 240 isHorizontal ? 1 : m_matrix.cols()); 241 } 242 243 template<typename OtherDerived> struct OppositeExtendedType { 244 typedef Replicate<OtherDerived, 245 isHorizontal ? 1 : ExpressionType::RowsAtCompileTime, 246 isVertical ? 1 : ExpressionType::ColsAtCompileTime> Type; 247 }; 248 249 /** \internal 250 * Replicates a vector in the opposite direction to match the size of \c *this */ 251 template<typename OtherDerived> 252 EIGEN_DEVICE_FUNC 253 typename OppositeExtendedType<OtherDerived>::Type 254 extendedToOpposite(const DenseBase<OtherDerived>& other) const 255 { 256 EIGEN_STATIC_ASSERT(EIGEN_IMPLIES(isHorizontal, OtherDerived::MaxColsAtCompileTime==1), 257 YOU_PASSED_A_ROW_VECTOR_BUT_A_COLUMN_VECTOR_WAS_EXPECTED) 258 EIGEN_STATIC_ASSERT(EIGEN_IMPLIES(isVertical, OtherDerived::MaxRowsAtCompileTime==1), 259 YOU_PASSED_A_COLUMN_VECTOR_BUT_A_ROW_VECTOR_WAS_EXPECTED) 260 return typename OppositeExtendedType<OtherDerived>::Type 261 (other.derived(), 262 isHorizontal ? 1 : m_matrix.rows(), 263 isVertical ? 1 : m_matrix.cols()); 264 } 265 266 public: 267 EIGEN_DEVICE_FUNC 268 explicit inline VectorwiseOp(ExpressionType& matrix) : m_matrix(matrix) {} 269 270 /** \internal */ 271 EIGEN_DEVICE_FUNC 272 inline const ExpressionType& _expression() const { return m_matrix; } 273 274 #ifdef EIGEN_PARSED_BY_DOXYGEN 275 /** STL-like <a href="https://en.cppreference.com/w/cpp/named_req/RandomAccessIterator">RandomAccessIterator</a> 276 * iterator type over the columns or rows as returned by the begin() and end() methods. 277 */ 278 random_access_iterator_type iterator; 279 /** This is the const version of iterator (aka read-only) */ 280 random_access_iterator_type const_iterator; 281 #else 282 typedef internal::subvector_stl_iterator<ExpressionType, DirectionType(Direction)> iterator; 283 typedef internal::subvector_stl_iterator<const ExpressionType, DirectionType(Direction)> const_iterator; 284 typedef internal::subvector_stl_reverse_iterator<ExpressionType, DirectionType(Direction)> reverse_iterator; 285 typedef internal::subvector_stl_reverse_iterator<const ExpressionType, DirectionType(Direction)> const_reverse_iterator; 286 #endif 287 288 /** returns an iterator to the first row (rowwise) or column (colwise) of the nested expression. 289 * \sa end(), cbegin() 290 */ 291 iterator begin() { return iterator (m_matrix, 0); } 292 /** const version of begin() */ 293 const_iterator begin() const { return const_iterator(m_matrix, 0); } 294 /** const version of begin() */ 295 const_iterator cbegin() const { return const_iterator(m_matrix, 0); } 296 297 /** returns a reverse iterator to the last row (rowwise) or column (colwise) of the nested expression. 298 * \sa rend(), crbegin() 299 */ 300 reverse_iterator rbegin() { return reverse_iterator (m_matrix, m_matrix.template subVectors<DirectionType(Direction)>()-1); } 301 /** const version of rbegin() */ 302 const_reverse_iterator rbegin() const { return const_reverse_iterator (m_matrix, m_matrix.template subVectors<DirectionType(Direction)>()-1); } 303 /** const version of rbegin() */ 304 const_reverse_iterator crbegin() const { return const_reverse_iterator (m_matrix, m_matrix.template subVectors<DirectionType(Direction)>()-1); } 305 306 /** returns an iterator to the row (resp. column) following the last row (resp. column) of the nested expression 307 * \sa begin(), cend() 308 */ 309 iterator end() { return iterator (m_matrix, m_matrix.template subVectors<DirectionType(Direction)>()); } 310 /** const version of end() */ 311 const_iterator end() const { return const_iterator(m_matrix, m_matrix.template subVectors<DirectionType(Direction)>()); } 312 /** const version of end() */ 313 const_iterator cend() const { return const_iterator(m_matrix, m_matrix.template subVectors<DirectionType(Direction)>()); } 314 315 /** returns a reverse iterator to the row (resp. column) before the first row (resp. column) of the nested expression 316 * \sa begin(), cend() 317 */ 318 reverse_iterator rend() { return reverse_iterator (m_matrix, -1); } 319 /** const version of rend() */ 320 const_reverse_iterator rend() const { return const_reverse_iterator (m_matrix, -1); } 321 /** const version of rend() */ 322 const_reverse_iterator crend() const { return const_reverse_iterator (m_matrix, -1); } 323 324 /** \returns a row or column vector expression of \c *this reduxed by \a func 325 * 326 * The template parameter \a BinaryOp is the type of the functor 327 * of the custom redux operator. Note that func must be an associative operator. 328 * 329 * \warning the size along the reduction direction must be strictly positive, 330 * otherwise an assertion is triggered. 331 * 332 * \sa class VectorwiseOp, DenseBase::colwise(), DenseBase::rowwise() 333 */ 334 template<typename BinaryOp> 335 EIGEN_DEVICE_FUNC 336 const typename ReduxReturnType<BinaryOp>::Type 337 redux(const BinaryOp& func = BinaryOp()) const 338 { 339 eigen_assert(redux_length()>0 && "you are using an empty matrix"); 340 return typename ReduxReturnType<BinaryOp>::Type(_expression(), internal::member_redux<BinaryOp,Scalar>(func)); 341 } 342 343 typedef typename ReturnType<internal::member_minCoeff>::Type MinCoeffReturnType; 344 typedef typename ReturnType<internal::member_maxCoeff>::Type MaxCoeffReturnType; 345 typedef PartialReduxExpr<const CwiseUnaryOp<internal::scalar_abs2_op<Scalar>, const ExpressionTypeNestedCleaned>,internal::member_sum<RealScalar,RealScalar>,Direction> SquaredNormReturnType; 346 typedef CwiseUnaryOp<internal::scalar_sqrt_op<RealScalar>, const SquaredNormReturnType> NormReturnType; 347 typedef typename ReturnType<internal::member_blueNorm,RealScalar>::Type BlueNormReturnType; 348 typedef typename ReturnType<internal::member_stableNorm,RealScalar>::Type StableNormReturnType; 349 typedef typename ReturnType<internal::member_hypotNorm,RealScalar>::Type HypotNormReturnType; 350 typedef typename ReturnType<internal::member_sum>::Type SumReturnType; 351 typedef EIGEN_EXPR_BINARYOP_SCALAR_RETURN_TYPE(SumReturnType,Scalar,quotient) MeanReturnType; 352 typedef typename ReturnType<internal::member_all>::Type AllReturnType; 353 typedef typename ReturnType<internal::member_any>::Type AnyReturnType; 354 typedef PartialReduxExpr<ExpressionType, internal::member_count<Index,Scalar>, Direction> CountReturnType; 355 typedef typename ReturnType<internal::member_prod>::Type ProdReturnType; 356 typedef Reverse<const ExpressionType, Direction> ConstReverseReturnType; 357 typedef Reverse<ExpressionType, Direction> ReverseReturnType; 358 359 template<int p> struct LpNormReturnType { 360 typedef PartialReduxExpr<ExpressionType, internal::member_lpnorm<p,RealScalar,Scalar>,Direction> Type; 361 }; 362 363 /** \returns a row (or column) vector expression of the smallest coefficient 364 * of each column (or row) of the referenced expression. 365 * 366 * \warning the size along the reduction direction must be strictly positive, 367 * otherwise an assertion is triggered. 368 * 369 * \warning the result is undefined if \c *this contains NaN. 370 * 371 * Example: \include PartialRedux_minCoeff.cpp 372 * Output: \verbinclude PartialRedux_minCoeff.out 373 * 374 * \sa DenseBase::minCoeff() */ 375 EIGEN_DEVICE_FUNC 376 const MinCoeffReturnType minCoeff() const 377 { 378 eigen_assert(redux_length()>0 && "you are using an empty matrix"); 379 return MinCoeffReturnType(_expression()); 380 } 381 382 /** \returns a row (or column) vector expression of the largest coefficient 383 * of each column (or row) of the referenced expression. 384 * 385 * \warning the size along the reduction direction must be strictly positive, 386 * otherwise an assertion is triggered. 387 * 388 * \warning the result is undefined if \c *this contains NaN. 389 * 390 * Example: \include PartialRedux_maxCoeff.cpp 391 * Output: \verbinclude PartialRedux_maxCoeff.out 392 * 393 * \sa DenseBase::maxCoeff() */ 394 EIGEN_DEVICE_FUNC 395 const MaxCoeffReturnType maxCoeff() const 396 { 397 eigen_assert(redux_length()>0 && "you are using an empty matrix"); 398 return MaxCoeffReturnType(_expression()); 399 } 400 401 /** \returns a row (or column) vector expression of the squared norm 402 * of each column (or row) of the referenced expression. 403 * This is a vector with real entries, even if the original matrix has complex entries. 404 * 405 * Example: \include PartialRedux_squaredNorm.cpp 406 * Output: \verbinclude PartialRedux_squaredNorm.out 407 * 408 * \sa DenseBase::squaredNorm() */ 409 EIGEN_DEVICE_FUNC 410 const SquaredNormReturnType squaredNorm() const 411 { return SquaredNormReturnType(m_matrix.cwiseAbs2()); } 412 413 /** \returns a row (or column) vector expression of the norm 414 * of each column (or row) of the referenced expression. 415 * This is a vector with real entries, even if the original matrix has complex entries. 416 * 417 * Example: \include PartialRedux_norm.cpp 418 * Output: \verbinclude PartialRedux_norm.out 419 * 420 * \sa DenseBase::norm() */ 421 EIGEN_DEVICE_FUNC 422 const NormReturnType norm() const 423 { return NormReturnType(squaredNorm()); } 424 425 /** \returns a row (or column) vector expression of the norm 426 * of each column (or row) of the referenced expression. 427 * This is a vector with real entries, even if the original matrix has complex entries. 428 * 429 * Example: \include PartialRedux_norm.cpp 430 * Output: \verbinclude PartialRedux_norm.out 431 * 432 * \sa DenseBase::norm() */ 433 template<int p> 434 EIGEN_DEVICE_FUNC 435 const typename LpNormReturnType<p>::Type lpNorm() const 436 { return typename LpNormReturnType<p>::Type(_expression()); } 437 438 439 /** \returns a row (or column) vector expression of the norm 440 * of each column (or row) of the referenced expression, using 441 * Blue's algorithm. 442 * This is a vector with real entries, even if the original matrix has complex entries. 443 * 444 * \sa DenseBase::blueNorm() */ 445 EIGEN_DEVICE_FUNC 446 const BlueNormReturnType blueNorm() const 447 { return BlueNormReturnType(_expression()); } 448 449 450 /** \returns a row (or column) vector expression of the norm 451 * of each column (or row) of the referenced expression, avoiding 452 * underflow and overflow. 453 * This is a vector with real entries, even if the original matrix has complex entries. 454 * 455 * \sa DenseBase::stableNorm() */ 456 EIGEN_DEVICE_FUNC 457 const StableNormReturnType stableNorm() const 458 { return StableNormReturnType(_expression()); } 459 460 461 /** \returns a row (or column) vector expression of the norm 462 * of each column (or row) of the referenced expression, avoiding 463 * underflow and overflow using a concatenation of hypot() calls. 464 * This is a vector with real entries, even if the original matrix has complex entries. 465 * 466 * \sa DenseBase::hypotNorm() */ 467 EIGEN_DEVICE_FUNC 468 const HypotNormReturnType hypotNorm() const 469 { return HypotNormReturnType(_expression()); } 470 471 /** \returns a row (or column) vector expression of the sum 472 * of each column (or row) of the referenced expression. 473 * 474 * Example: \include PartialRedux_sum.cpp 475 * Output: \verbinclude PartialRedux_sum.out 476 * 477 * \sa DenseBase::sum() */ 478 EIGEN_DEVICE_FUNC 479 const SumReturnType sum() const 480 { return SumReturnType(_expression()); } 481 482 /** \returns a row (or column) vector expression of the mean 483 * of each column (or row) of the referenced expression. 484 * 485 * \sa DenseBase::mean() */ 486 EIGEN_DEVICE_FUNC 487 const MeanReturnType mean() const 488 { return sum() / Scalar(Direction==Vertical?m_matrix.rows():m_matrix.cols()); } 489 490 /** \returns a row (or column) vector expression representing 491 * whether \b all coefficients of each respective column (or row) are \c true. 492 * This expression can be assigned to a vector with entries of type \c bool. 493 * 494 * \sa DenseBase::all() */ 495 EIGEN_DEVICE_FUNC 496 const AllReturnType all() const 497 { return AllReturnType(_expression()); } 498 499 /** \returns a row (or column) vector expression representing 500 * whether \b at \b least one coefficient of each respective column (or row) is \c true. 501 * This expression can be assigned to a vector with entries of type \c bool. 502 * 503 * \sa DenseBase::any() */ 504 EIGEN_DEVICE_FUNC 505 const AnyReturnType any() const 506 { return AnyReturnType(_expression()); } 507 508 /** \returns a row (or column) vector expression representing 509 * the number of \c true coefficients of each respective column (or row). 510 * This expression can be assigned to a vector whose entries have the same type as is used to 511 * index entries of the original matrix; for dense matrices, this is \c std::ptrdiff_t . 512 * 513 * Example: \include PartialRedux_count.cpp 514 * Output: \verbinclude PartialRedux_count.out 515 * 516 * \sa DenseBase::count() */ 517 EIGEN_DEVICE_FUNC 518 const CountReturnType count() const 519 { return CountReturnType(_expression()); } 520 521 /** \returns a row (or column) vector expression of the product 522 * of each column (or row) of the referenced expression. 523 * 524 * Example: \include PartialRedux_prod.cpp 525 * Output: \verbinclude PartialRedux_prod.out 526 * 527 * \sa DenseBase::prod() */ 528 EIGEN_DEVICE_FUNC 529 const ProdReturnType prod() const 530 { return ProdReturnType(_expression()); } 531 532 533 /** \returns a matrix expression 534 * where each column (or row) are reversed. 535 * 536 * Example: \include Vectorwise_reverse.cpp 537 * Output: \verbinclude Vectorwise_reverse.out 538 * 539 * \sa DenseBase::reverse() */ 540 EIGEN_DEVICE_FUNC 541 const ConstReverseReturnType reverse() const 542 { return ConstReverseReturnType( _expression() ); } 543 544 /** \returns a writable matrix expression 545 * where each column (or row) are reversed. 546 * 547 * \sa reverse() const */ 548 EIGEN_DEVICE_FUNC 549 ReverseReturnType reverse() 550 { return ReverseReturnType( _expression() ); } 551 552 typedef Replicate<ExpressionType,(isVertical?Dynamic:1),(isHorizontal?Dynamic:1)> ReplicateReturnType; 553 EIGEN_DEVICE_FUNC 554 const ReplicateReturnType replicate(Index factor) const; 555 556 /** 557 * \return an expression of the replication of each column (or row) of \c *this 558 * 559 * Example: \include DirectionWise_replicate.cpp 560 * Output: \verbinclude DirectionWise_replicate.out 561 * 562 * \sa VectorwiseOp::replicate(Index), DenseBase::replicate(), class Replicate 563 */ 564 // NOTE implemented here because of sunstudio's compilation errors 565 // isVertical*Factor+isHorizontal instead of (isVertical?Factor:1) to handle CUDA bug with ternary operator 566 template<int Factor> const Replicate<ExpressionType,isVertical*Factor+isHorizontal,isHorizontal*Factor+isVertical> 567 EIGEN_DEVICE_FUNC 568 replicate(Index factor = Factor) const 569 { 570 return Replicate<ExpressionType,(isVertical?Factor:1),(isHorizontal?Factor:1)> 571 (_expression(),isVertical?factor:1,isHorizontal?factor:1); 572 } 573 574 /////////// Artithmetic operators /////////// 575 576 /** Copies the vector \a other to each subvector of \c *this */ 577 template<typename OtherDerived> 578 EIGEN_DEVICE_FUNC 579 ExpressionType& operator=(const DenseBase<OtherDerived>& other) 580 { 581 EIGEN_STATIC_ASSERT_VECTOR_ONLY(OtherDerived) 582 EIGEN_STATIC_ASSERT_SAME_XPR_KIND(ExpressionType, OtherDerived) 583 //eigen_assert((m_matrix.isNull()) == (other.isNull())); FIXME 584 return m_matrix = extendedTo(other.derived()); 585 } 586 587 /** Adds the vector \a other to each subvector of \c *this */ 588 template<typename OtherDerived> 589 EIGEN_DEVICE_FUNC 590 ExpressionType& operator+=(const DenseBase<OtherDerived>& other) 591 { 592 EIGEN_STATIC_ASSERT_VECTOR_ONLY(OtherDerived) 593 EIGEN_STATIC_ASSERT_SAME_XPR_KIND(ExpressionType, OtherDerived) 594 return m_matrix += extendedTo(other.derived()); 595 } 596 597 /** Substracts the vector \a other to each subvector of \c *this */ 598 template<typename OtherDerived> 599 EIGEN_DEVICE_FUNC 600 ExpressionType& operator-=(const DenseBase<OtherDerived>& other) 601 { 602 EIGEN_STATIC_ASSERT_VECTOR_ONLY(OtherDerived) 603 EIGEN_STATIC_ASSERT_SAME_XPR_KIND(ExpressionType, OtherDerived) 604 return m_matrix -= extendedTo(other.derived()); 605 } 606 607 /** Multiples each subvector of \c *this by the vector \a other */ 608 template<typename OtherDerived> 609 EIGEN_DEVICE_FUNC 610 ExpressionType& operator*=(const DenseBase<OtherDerived>& other) 611 { 612 EIGEN_STATIC_ASSERT_VECTOR_ONLY(OtherDerived) 613 EIGEN_STATIC_ASSERT_ARRAYXPR(ExpressionType) 614 EIGEN_STATIC_ASSERT_SAME_XPR_KIND(ExpressionType, OtherDerived) 615 m_matrix *= extendedTo(other.derived()); 616 return m_matrix; 617 } 618 619 /** Divides each subvector of \c *this by the vector \a other */ 620 template<typename OtherDerived> 621 EIGEN_DEVICE_FUNC 622 ExpressionType& operator/=(const DenseBase<OtherDerived>& other) 623 { 624 EIGEN_STATIC_ASSERT_VECTOR_ONLY(OtherDerived) 625 EIGEN_STATIC_ASSERT_ARRAYXPR(ExpressionType) 626 EIGEN_STATIC_ASSERT_SAME_XPR_KIND(ExpressionType, OtherDerived) 627 m_matrix /= extendedTo(other.derived()); 628 return m_matrix; 629 } 630 631 /** Returns the expression of the sum of the vector \a other to each subvector of \c *this */ 632 template<typename OtherDerived> EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC 633 CwiseBinaryOp<internal::scalar_sum_op<Scalar,typename OtherDerived::Scalar>, 634 const ExpressionTypeNestedCleaned, 635 const typename ExtendedType<OtherDerived>::Type> 636 operator+(const DenseBase<OtherDerived>& other) const 637 { 638 EIGEN_STATIC_ASSERT_VECTOR_ONLY(OtherDerived) 639 EIGEN_STATIC_ASSERT_SAME_XPR_KIND(ExpressionType, OtherDerived) 640 return m_matrix + extendedTo(other.derived()); 641 } 642 643 /** Returns the expression of the difference between each subvector of \c *this and the vector \a other */ 644 template<typename OtherDerived> 645 EIGEN_DEVICE_FUNC 646 CwiseBinaryOp<internal::scalar_difference_op<Scalar,typename OtherDerived::Scalar>, 647 const ExpressionTypeNestedCleaned, 648 const typename ExtendedType<OtherDerived>::Type> 649 operator-(const DenseBase<OtherDerived>& other) const 650 { 651 EIGEN_STATIC_ASSERT_VECTOR_ONLY(OtherDerived) 652 EIGEN_STATIC_ASSERT_SAME_XPR_KIND(ExpressionType, OtherDerived) 653 return m_matrix - extendedTo(other.derived()); 654 } 655 656 /** Returns the expression where each subvector is the product of the vector \a other 657 * by the corresponding subvector of \c *this */ 658 template<typename OtherDerived> EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC 659 CwiseBinaryOp<internal::scalar_product_op<Scalar>, 660 const ExpressionTypeNestedCleaned, 661 const typename ExtendedType<OtherDerived>::Type> 662 EIGEN_DEVICE_FUNC 663 operator*(const DenseBase<OtherDerived>& other) const 664 { 665 EIGEN_STATIC_ASSERT_VECTOR_ONLY(OtherDerived) 666 EIGEN_STATIC_ASSERT_ARRAYXPR(ExpressionType) 667 EIGEN_STATIC_ASSERT_SAME_XPR_KIND(ExpressionType, OtherDerived) 668 return m_matrix * extendedTo(other.derived()); 669 } 670 671 /** Returns the expression where each subvector is the quotient of the corresponding 672 * subvector of \c *this by the vector \a other */ 673 template<typename OtherDerived> 674 EIGEN_DEVICE_FUNC 675 CwiseBinaryOp<internal::scalar_quotient_op<Scalar>, 676 const ExpressionTypeNestedCleaned, 677 const typename ExtendedType<OtherDerived>::Type> 678 operator/(const DenseBase<OtherDerived>& other) const 679 { 680 EIGEN_STATIC_ASSERT_VECTOR_ONLY(OtherDerived) 681 EIGEN_STATIC_ASSERT_ARRAYXPR(ExpressionType) 682 EIGEN_STATIC_ASSERT_SAME_XPR_KIND(ExpressionType, OtherDerived) 683 return m_matrix / extendedTo(other.derived()); 684 } 685 686 /** \returns an expression where each column (or row) of the referenced matrix are normalized. 687 * The referenced matrix is \b not modified. 688 * \sa MatrixBase::normalized(), normalize() 689 */ 690 EIGEN_DEVICE_FUNC 691 CwiseBinaryOp<internal::scalar_quotient_op<Scalar>, 692 const ExpressionTypeNestedCleaned, 693 const typename OppositeExtendedType<NormReturnType>::Type> 694 normalized() const { return m_matrix.cwiseQuotient(extendedToOpposite(this->norm())); } 695 696 697 /** Normalize in-place each row or columns of the referenced matrix. 698 * \sa MatrixBase::normalize(), normalized() 699 */ 700 EIGEN_DEVICE_FUNC void normalize() { 701 m_matrix = this->normalized(); 702 } 703 704 EIGEN_DEVICE_FUNC inline void reverseInPlace(); 705 706 /////////// Geometry module /////////// 707 708 typedef Homogeneous<ExpressionType,Direction> HomogeneousReturnType; 709 EIGEN_DEVICE_FUNC 710 HomogeneousReturnType homogeneous() const; 711 712 typedef typename ExpressionType::PlainObject CrossReturnType; 713 template<typename OtherDerived> 714 EIGEN_DEVICE_FUNC 715 const CrossReturnType cross(const MatrixBase<OtherDerived>& other) const; 716 717 enum { 718 HNormalized_Size = Direction==Vertical ? internal::traits<ExpressionType>::RowsAtCompileTime 719 : internal::traits<ExpressionType>::ColsAtCompileTime, 720 HNormalized_SizeMinusOne = HNormalized_Size==Dynamic ? Dynamic : HNormalized_Size-1 721 }; 722 typedef Block<const ExpressionType, 723 Direction==Vertical ? int(HNormalized_SizeMinusOne) 724 : int(internal::traits<ExpressionType>::RowsAtCompileTime), 725 Direction==Horizontal ? int(HNormalized_SizeMinusOne) 726 : int(internal::traits<ExpressionType>::ColsAtCompileTime)> 727 HNormalized_Block; 728 typedef Block<const ExpressionType, 729 Direction==Vertical ? 1 : int(internal::traits<ExpressionType>::RowsAtCompileTime), 730 Direction==Horizontal ? 1 : int(internal::traits<ExpressionType>::ColsAtCompileTime)> 731 HNormalized_Factors; 732 typedef CwiseBinaryOp<internal::scalar_quotient_op<typename internal::traits<ExpressionType>::Scalar>, 733 const HNormalized_Block, 734 const Replicate<HNormalized_Factors, 735 Direction==Vertical ? HNormalized_SizeMinusOne : 1, 736 Direction==Horizontal ? HNormalized_SizeMinusOne : 1> > 737 HNormalizedReturnType; 738 739 EIGEN_DEVICE_FUNC 740 const HNormalizedReturnType hnormalized() const; 741 742 # ifdef EIGEN_VECTORWISEOP_PLUGIN 743 # include EIGEN_VECTORWISEOP_PLUGIN 744 # endif 745 746 protected: 747 Index redux_length() const 748 { 749 return Direction==Vertical ? m_matrix.rows() : m_matrix.cols(); 750 } 751 ExpressionTypeNested m_matrix; 752 }; 753 754 //const colwise moved to DenseBase.h due to CUDA compiler bug 755 756 757 /** \returns a writable VectorwiseOp wrapper of *this providing additional partial reduction operations 758 * 759 * \sa rowwise(), class VectorwiseOp, \ref TutorialReductionsVisitorsBroadcasting 760 */ 761 template<typename Derived> 762 EIGEN_DEVICE_FUNC inline typename DenseBase<Derived>::ColwiseReturnType 763 DenseBase<Derived>::colwise() 764 { 765 return ColwiseReturnType(derived()); 766 } 767 768 //const rowwise moved to DenseBase.h due to CUDA compiler bug 769 770 771 /** \returns a writable VectorwiseOp wrapper of *this providing additional partial reduction operations 772 * 773 * \sa colwise(), class VectorwiseOp, \ref TutorialReductionsVisitorsBroadcasting 774 */ 775 template<typename Derived> 776 EIGEN_DEVICE_FUNC inline typename DenseBase<Derived>::RowwiseReturnType 777 DenseBase<Derived>::rowwise() 778 { 779 return RowwiseReturnType(derived()); 780 } 781 782 } // end namespace Eigen 783 784 #endif // EIGEN_PARTIAL_REDUX_H 785