Home | History | Annotate | Download | only in Core
      1 // This file is part of Eigen, a lightweight C++ template library
      2 // for linear algebra.
      3 //
      4 // Copyright (C) 2008-2010 Gael Guennebaud <gael.guennebaud (at) inria.fr>
      5 //
      6 // This Source Code Form is subject to the terms of the Mozilla
      7 // Public License v. 2.0. If a copy of the MPL was not distributed
      8 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
      9 
     10 #ifndef EIGEN_FUNCTORS_H
     11 #define EIGEN_FUNCTORS_H
     12 
     13 namespace Eigen {
     14 
     15 namespace internal {
     16 
     17 // associative functors:
     18 
     19 /** \internal
     20   * \brief Template functor to compute the sum of two scalars
     21   *
     22   * \sa class CwiseBinaryOp, MatrixBase::operator+, class VectorwiseOp, MatrixBase::sum()
     23   */
     24 template<typename Scalar> struct scalar_sum_op {
     25   EIGEN_EMPTY_STRUCT_CTOR(scalar_sum_op)
     26   EIGEN_STRONG_INLINE const Scalar operator() (const Scalar& a, const Scalar& b) const { return a + b; }
     27   template<typename Packet>
     28   EIGEN_STRONG_INLINE const Packet packetOp(const Packet& a, const Packet& b) const
     29   { return internal::padd(a,b); }
     30   template<typename Packet>
     31   EIGEN_STRONG_INLINE const Scalar predux(const Packet& a) const
     32   { return internal::predux(a); }
     33 };
     34 template<typename Scalar>
     35 struct functor_traits<scalar_sum_op<Scalar> > {
     36   enum {
     37     Cost = NumTraits<Scalar>::AddCost,
     38     PacketAccess = packet_traits<Scalar>::HasAdd
     39   };
     40 };
     41 
     42 /** \internal
     43   * \brief Template functor to compute the product of two scalars
     44   *
     45   * \sa class CwiseBinaryOp, Cwise::operator*(), class VectorwiseOp, MatrixBase::redux()
     46   */
     47 template<typename LhsScalar,typename RhsScalar> struct scalar_product_op {
     48   enum {
     49     // TODO vectorize mixed product
     50     Vectorizable = is_same<LhsScalar,RhsScalar>::value && packet_traits<LhsScalar>::HasMul && packet_traits<RhsScalar>::HasMul
     51   };
     52   typedef typename scalar_product_traits<LhsScalar,RhsScalar>::ReturnType result_type;
     53   EIGEN_EMPTY_STRUCT_CTOR(scalar_product_op)
     54   EIGEN_STRONG_INLINE const result_type operator() (const LhsScalar& a, const RhsScalar& b) const { return a * b; }
     55   template<typename Packet>
     56   EIGEN_STRONG_INLINE const Packet packetOp(const Packet& a, const Packet& b) const
     57   { return internal::pmul(a,b); }
     58   template<typename Packet>
     59   EIGEN_STRONG_INLINE const result_type predux(const Packet& a) const
     60   { return internal::predux_mul(a); }
     61 };
     62 template<typename LhsScalar,typename RhsScalar>
     63 struct functor_traits<scalar_product_op<LhsScalar,RhsScalar> > {
     64   enum {
     65     Cost = (NumTraits<LhsScalar>::MulCost + NumTraits<RhsScalar>::MulCost)/2, // rough estimate!
     66     PacketAccess = scalar_product_op<LhsScalar,RhsScalar>::Vectorizable
     67   };
     68 };
     69 
     70 /** \internal
     71   * \brief Template functor to compute the conjugate product of two scalars
     72   *
     73   * This is a short cut for conj(x) * y which is needed for optimization purpose; in Eigen2 support mode, this becomes x * conj(y)
     74   */
     75 template<typename LhsScalar,typename RhsScalar> struct scalar_conj_product_op {
     76 
     77   enum {
     78     Conj = NumTraits<LhsScalar>::IsComplex
     79   };
     80 
     81   typedef typename scalar_product_traits<LhsScalar,RhsScalar>::ReturnType result_type;
     82 
     83   EIGEN_EMPTY_STRUCT_CTOR(scalar_conj_product_op)
     84   EIGEN_STRONG_INLINE const result_type operator() (const LhsScalar& a, const RhsScalar& b) const
     85   { return conj_helper<LhsScalar,RhsScalar,Conj,false>().pmul(a,b); }
     86 
     87   template<typename Packet>
     88   EIGEN_STRONG_INLINE const Packet packetOp(const Packet& a, const Packet& b) const
     89   { return conj_helper<Packet,Packet,Conj,false>().pmul(a,b); }
     90 };
     91 template<typename LhsScalar,typename RhsScalar>
     92 struct functor_traits<scalar_conj_product_op<LhsScalar,RhsScalar> > {
     93   enum {
     94     Cost = NumTraits<LhsScalar>::MulCost,
     95     PacketAccess = internal::is_same<LhsScalar, RhsScalar>::value && packet_traits<LhsScalar>::HasMul
     96   };
     97 };
     98 
     99 /** \internal
    100   * \brief Template functor to compute the min of two scalars
    101   *
    102   * \sa class CwiseBinaryOp, MatrixBase::cwiseMin, class VectorwiseOp, MatrixBase::minCoeff()
    103   */
    104 template<typename Scalar> struct scalar_min_op {
    105   EIGEN_EMPTY_STRUCT_CTOR(scalar_min_op)
    106   EIGEN_STRONG_INLINE const Scalar operator() (const Scalar& a, const Scalar& b) const { using std::min; return (min)(a, b); }
    107   template<typename Packet>
    108   EIGEN_STRONG_INLINE const Packet packetOp(const Packet& a, const Packet& b) const
    109   { return internal::pmin(a,b); }
    110   template<typename Packet>
    111   EIGEN_STRONG_INLINE const Scalar predux(const Packet& a) const
    112   { return internal::predux_min(a); }
    113 };
    114 template<typename Scalar>
    115 struct functor_traits<scalar_min_op<Scalar> > {
    116   enum {
    117     Cost = NumTraits<Scalar>::AddCost,
    118     PacketAccess = packet_traits<Scalar>::HasMin
    119   };
    120 };
    121 
    122 /** \internal
    123   * \brief Template functor to compute the max of two scalars
    124   *
    125   * \sa class CwiseBinaryOp, MatrixBase::cwiseMax, class VectorwiseOp, MatrixBase::maxCoeff()
    126   */
    127 template<typename Scalar> struct scalar_max_op {
    128   EIGEN_EMPTY_STRUCT_CTOR(scalar_max_op)
    129   EIGEN_STRONG_INLINE const Scalar operator() (const Scalar& a, const Scalar& b) const { using std::max; return (max)(a, b); }
    130   template<typename Packet>
    131   EIGEN_STRONG_INLINE const Packet packetOp(const Packet& a, const Packet& b) const
    132   { return internal::pmax(a,b); }
    133   template<typename Packet>
    134   EIGEN_STRONG_INLINE const Scalar predux(const Packet& a) const
    135   { return internal::predux_max(a); }
    136 };
    137 template<typename Scalar>
    138 struct functor_traits<scalar_max_op<Scalar> > {
    139   enum {
    140     Cost = NumTraits<Scalar>::AddCost,
    141     PacketAccess = packet_traits<Scalar>::HasMax
    142   };
    143 };
    144 
    145 /** \internal
    146   * \brief Template functor to compute the hypot of two scalars
    147   *
    148   * \sa MatrixBase::stableNorm(), class Redux
    149   */
    150 template<typename Scalar> struct scalar_hypot_op {
    151   EIGEN_EMPTY_STRUCT_CTOR(scalar_hypot_op)
    152 //   typedef typename NumTraits<Scalar>::Real result_type;
    153   EIGEN_STRONG_INLINE const Scalar operator() (const Scalar& _x, const Scalar& _y) const
    154   {
    155     using std::max;
    156     using std::min;
    157     using std::sqrt;
    158     Scalar p = (max)(_x, _y);
    159     Scalar q = (min)(_x, _y);
    160     Scalar qp = q/p;
    161     return p * sqrt(Scalar(1) + qp*qp);
    162   }
    163 };
    164 template<typename Scalar>
    165 struct functor_traits<scalar_hypot_op<Scalar> > {
    166   enum { Cost = 5 * NumTraits<Scalar>::MulCost, PacketAccess=0 };
    167 };
    168 
    169 /** \internal
    170   * \brief Template functor to compute the pow of two scalars
    171   */
    172 template<typename Scalar, typename OtherScalar> struct scalar_binary_pow_op {
    173   EIGEN_EMPTY_STRUCT_CTOR(scalar_binary_pow_op)
    174   inline Scalar operator() (const Scalar& a, const OtherScalar& b) const { return numext::pow(a, b); }
    175 };
    176 template<typename Scalar, typename OtherScalar>
    177 struct functor_traits<scalar_binary_pow_op<Scalar,OtherScalar> > {
    178   enum { Cost = 5 * NumTraits<Scalar>::MulCost, PacketAccess = false };
    179 };
    180 
    181 // other binary functors:
    182 
    183 /** \internal
    184   * \brief Template functor to compute the difference of two scalars
    185   *
    186   * \sa class CwiseBinaryOp, MatrixBase::operator-
    187   */
    188 template<typename Scalar> struct scalar_difference_op {
    189   EIGEN_EMPTY_STRUCT_CTOR(scalar_difference_op)
    190   EIGEN_STRONG_INLINE const Scalar operator() (const Scalar& a, const Scalar& b) const { return a - b; }
    191   template<typename Packet>
    192   EIGEN_STRONG_INLINE const Packet packetOp(const Packet& a, const Packet& b) const
    193   { return internal::psub(a,b); }
    194 };
    195 template<typename Scalar>
    196 struct functor_traits<scalar_difference_op<Scalar> > {
    197   enum {
    198     Cost = NumTraits<Scalar>::AddCost,
    199     PacketAccess = packet_traits<Scalar>::HasSub
    200   };
    201 };
    202 
    203 /** \internal
    204   * \brief Template functor to compute the quotient of two scalars
    205   *
    206   * \sa class CwiseBinaryOp, Cwise::operator/()
    207   */
    208 template<typename LhsScalar,typename RhsScalar> struct scalar_quotient_op {
    209   enum {
    210     // TODO vectorize mixed product
    211     Vectorizable = is_same<LhsScalar,RhsScalar>::value && packet_traits<LhsScalar>::HasDiv && packet_traits<RhsScalar>::HasDiv
    212   };
    213   typedef typename scalar_product_traits<LhsScalar,RhsScalar>::ReturnType result_type;
    214   EIGEN_EMPTY_STRUCT_CTOR(scalar_quotient_op)
    215   EIGEN_STRONG_INLINE const result_type operator() (const LhsScalar& a, const RhsScalar& b) const { return a / b; }
    216   template<typename Packet>
    217   EIGEN_STRONG_INLINE const Packet packetOp(const Packet& a, const Packet& b) const
    218   { return internal::pdiv(a,b); }
    219 };
    220 template<typename LhsScalar,typename RhsScalar>
    221 struct functor_traits<scalar_quotient_op<LhsScalar,RhsScalar> > {
    222   enum {
    223     Cost = (NumTraits<LhsScalar>::MulCost + NumTraits<RhsScalar>::MulCost), // rough estimate!
    224     PacketAccess = scalar_quotient_op<LhsScalar,RhsScalar>::Vectorizable
    225   };
    226 };
    227 
    228 
    229 
    230 /** \internal
    231   * \brief Template functor to compute the and of two booleans
    232   *
    233   * \sa class CwiseBinaryOp, ArrayBase::operator&&
    234   */
    235 struct scalar_boolean_and_op {
    236   EIGEN_EMPTY_STRUCT_CTOR(scalar_boolean_and_op)
    237   EIGEN_STRONG_INLINE bool operator() (const bool& a, const bool& b) const { return a && b; }
    238 };
    239 template<> struct functor_traits<scalar_boolean_and_op> {
    240   enum {
    241     Cost = NumTraits<bool>::AddCost,
    242     PacketAccess = false
    243   };
    244 };
    245 
    246 /** \internal
    247   * \brief Template functor to compute the or of two booleans
    248   *
    249   * \sa class CwiseBinaryOp, ArrayBase::operator||
    250   */
    251 struct scalar_boolean_or_op {
    252   EIGEN_EMPTY_STRUCT_CTOR(scalar_boolean_or_op)
    253   EIGEN_STRONG_INLINE bool operator() (const bool& a, const bool& b) const { return a || b; }
    254 };
    255 template<> struct functor_traits<scalar_boolean_or_op> {
    256   enum {
    257     Cost = NumTraits<bool>::AddCost,
    258     PacketAccess = false
    259   };
    260 };
    261 
    262 // unary functors:
    263 
    264 /** \internal
    265   * \brief Template functor to compute the opposite of a scalar
    266   *
    267   * \sa class CwiseUnaryOp, MatrixBase::operator-
    268   */
    269 template<typename Scalar> struct scalar_opposite_op {
    270   EIGEN_EMPTY_STRUCT_CTOR(scalar_opposite_op)
    271   EIGEN_STRONG_INLINE const Scalar operator() (const Scalar& a) const { return -a; }
    272   template<typename Packet>
    273   EIGEN_STRONG_INLINE const Packet packetOp(const Packet& a) const
    274   { return internal::pnegate(a); }
    275 };
    276 template<typename Scalar>
    277 struct functor_traits<scalar_opposite_op<Scalar> >
    278 { enum {
    279     Cost = NumTraits<Scalar>::AddCost,
    280     PacketAccess = packet_traits<Scalar>::HasNegate };
    281 };
    282 
    283 /** \internal
    284   * \brief Template functor to compute the absolute value of a scalar
    285   *
    286   * \sa class CwiseUnaryOp, Cwise::abs
    287   */
    288 template<typename Scalar> struct scalar_abs_op {
    289   EIGEN_EMPTY_STRUCT_CTOR(scalar_abs_op)
    290   typedef typename NumTraits<Scalar>::Real result_type;
    291   EIGEN_STRONG_INLINE const result_type operator() (const Scalar& a) const { using std::abs; return abs(a); }
    292   template<typename Packet>
    293   EIGEN_STRONG_INLINE const Packet packetOp(const Packet& a) const
    294   { return internal::pabs(a); }
    295 };
    296 template<typename Scalar>
    297 struct functor_traits<scalar_abs_op<Scalar> >
    298 {
    299   enum {
    300     Cost = NumTraits<Scalar>::AddCost,
    301     PacketAccess = packet_traits<Scalar>::HasAbs
    302   };
    303 };
    304 
    305 /** \internal
    306   * \brief Template functor to compute the squared absolute value of a scalar
    307   *
    308   * \sa class CwiseUnaryOp, Cwise::abs2
    309   */
    310 template<typename Scalar> struct scalar_abs2_op {
    311   EIGEN_EMPTY_STRUCT_CTOR(scalar_abs2_op)
    312   typedef typename NumTraits<Scalar>::Real result_type;
    313   EIGEN_STRONG_INLINE const result_type operator() (const Scalar& a) const { return numext::abs2(a); }
    314   template<typename Packet>
    315   EIGEN_STRONG_INLINE const Packet packetOp(const Packet& a) const
    316   { return internal::pmul(a,a); }
    317 };
    318 template<typename Scalar>
    319 struct functor_traits<scalar_abs2_op<Scalar> >
    320 { enum { Cost = NumTraits<Scalar>::MulCost, PacketAccess = packet_traits<Scalar>::HasAbs2 }; };
    321 
    322 /** \internal
    323   * \brief Template functor to compute the conjugate of a complex value
    324   *
    325   * \sa class CwiseUnaryOp, MatrixBase::conjugate()
    326   */
    327 template<typename Scalar> struct scalar_conjugate_op {
    328   EIGEN_EMPTY_STRUCT_CTOR(scalar_conjugate_op)
    329   EIGEN_STRONG_INLINE const Scalar operator() (const Scalar& a) const { using numext::conj; return conj(a); }
    330   template<typename Packet>
    331   EIGEN_STRONG_INLINE const Packet packetOp(const Packet& a) const { return internal::pconj(a); }
    332 };
    333 template<typename Scalar>
    334 struct functor_traits<scalar_conjugate_op<Scalar> >
    335 {
    336   enum {
    337     Cost = NumTraits<Scalar>::IsComplex ? NumTraits<Scalar>::AddCost : 0,
    338     PacketAccess = packet_traits<Scalar>::HasConj
    339   };
    340 };
    341 
    342 /** \internal
    343   * \brief Template functor to cast a scalar to another type
    344   *
    345   * \sa class CwiseUnaryOp, MatrixBase::cast()
    346   */
    347 template<typename Scalar, typename NewType>
    348 struct scalar_cast_op {
    349   EIGEN_EMPTY_STRUCT_CTOR(scalar_cast_op)
    350   typedef NewType result_type;
    351   EIGEN_STRONG_INLINE const NewType operator() (const Scalar& a) const { return cast<Scalar, NewType>(a); }
    352 };
    353 template<typename Scalar, typename NewType>
    354 struct functor_traits<scalar_cast_op<Scalar,NewType> >
    355 { enum { Cost = is_same<Scalar, NewType>::value ? 0 : NumTraits<NewType>::AddCost, PacketAccess = false }; };
    356 
    357 /** \internal
    358   * \brief Template functor to extract the real part of a complex
    359   *
    360   * \sa class CwiseUnaryOp, MatrixBase::real()
    361   */
    362 template<typename Scalar>
    363 struct scalar_real_op {
    364   EIGEN_EMPTY_STRUCT_CTOR(scalar_real_op)
    365   typedef typename NumTraits<Scalar>::Real result_type;
    366   EIGEN_STRONG_INLINE result_type operator() (const Scalar& a) const { return numext::real(a); }
    367 };
    368 template<typename Scalar>
    369 struct functor_traits<scalar_real_op<Scalar> >
    370 { enum { Cost = 0, PacketAccess = false }; };
    371 
    372 /** \internal
    373   * \brief Template functor to extract the imaginary part of a complex
    374   *
    375   * \sa class CwiseUnaryOp, MatrixBase::imag()
    376   */
    377 template<typename Scalar>
    378 struct scalar_imag_op {
    379   EIGEN_EMPTY_STRUCT_CTOR(scalar_imag_op)
    380   typedef typename NumTraits<Scalar>::Real result_type;
    381   EIGEN_STRONG_INLINE result_type operator() (const Scalar& a) const { return numext::imag(a); }
    382 };
    383 template<typename Scalar>
    384 struct functor_traits<scalar_imag_op<Scalar> >
    385 { enum { Cost = 0, PacketAccess = false }; };
    386 
    387 /** \internal
    388   * \brief Template functor to extract the real part of a complex as a reference
    389   *
    390   * \sa class CwiseUnaryOp, MatrixBase::real()
    391   */
    392 template<typename Scalar>
    393 struct scalar_real_ref_op {
    394   EIGEN_EMPTY_STRUCT_CTOR(scalar_real_ref_op)
    395   typedef typename NumTraits<Scalar>::Real result_type;
    396   EIGEN_STRONG_INLINE result_type& operator() (const Scalar& a) const { return numext::real_ref(*const_cast<Scalar*>(&a)); }
    397 };
    398 template<typename Scalar>
    399 struct functor_traits<scalar_real_ref_op<Scalar> >
    400 { enum { Cost = 0, PacketAccess = false }; };
    401 
    402 /** \internal
    403   * \brief Template functor to extract the imaginary part of a complex as a reference
    404   *
    405   * \sa class CwiseUnaryOp, MatrixBase::imag()
    406   */
    407 template<typename Scalar>
    408 struct scalar_imag_ref_op {
    409   EIGEN_EMPTY_STRUCT_CTOR(scalar_imag_ref_op)
    410   typedef typename NumTraits<Scalar>::Real result_type;
    411   EIGEN_STRONG_INLINE result_type& operator() (const Scalar& a) const { return numext::imag_ref(*const_cast<Scalar*>(&a)); }
    412 };
    413 template<typename Scalar>
    414 struct functor_traits<scalar_imag_ref_op<Scalar> >
    415 { enum { Cost = 0, PacketAccess = false }; };
    416 
    417 /** \internal
    418   *
    419   * \brief Template functor to compute the exponential of a scalar
    420   *
    421   * \sa class CwiseUnaryOp, Cwise::exp()
    422   */
    423 template<typename Scalar> struct scalar_exp_op {
    424   EIGEN_EMPTY_STRUCT_CTOR(scalar_exp_op)
    425   inline const Scalar operator() (const Scalar& a) const { using std::exp; return exp(a); }
    426   typedef typename packet_traits<Scalar>::type Packet;
    427   inline Packet packetOp(const Packet& a) const { return internal::pexp(a); }
    428 };
    429 template<typename Scalar>
    430 struct functor_traits<scalar_exp_op<Scalar> >
    431 { enum { Cost = 5 * NumTraits<Scalar>::MulCost, PacketAccess = packet_traits<Scalar>::HasExp }; };
    432 
    433 /** \internal
    434   *
    435   * \brief Template functor to compute the logarithm of a scalar
    436   *
    437   * \sa class CwiseUnaryOp, Cwise::log()
    438   */
    439 template<typename Scalar> struct scalar_log_op {
    440   EIGEN_EMPTY_STRUCT_CTOR(scalar_log_op)
    441   inline const Scalar operator() (const Scalar& a) const { using std::log; return log(a); }
    442   typedef typename packet_traits<Scalar>::type Packet;
    443   inline Packet packetOp(const Packet& a) const { return internal::plog(a); }
    444 };
    445 template<typename Scalar>
    446 struct functor_traits<scalar_log_op<Scalar> >
    447 { enum { Cost = 5 * NumTraits<Scalar>::MulCost, PacketAccess = packet_traits<Scalar>::HasLog }; };
    448 
    449 /** \internal
    450   * \brief Template functor to multiply a scalar by a fixed other one
    451   *
    452   * \sa class CwiseUnaryOp, MatrixBase::operator*, MatrixBase::operator/
    453   */
    454 /* NOTE why doing the pset1() in packetOp *is* an optimization ?
    455  * indeed it seems better to declare m_other as a Packet and do the pset1() once
    456  * in the constructor. However, in practice:
    457  *  - GCC does not like m_other as a Packet and generate a load every time it needs it
    458  *  - on the other hand GCC is able to moves the pset1() outside the loop :)
    459  *  - simpler code ;)
    460  * (ICC and gcc 4.4 seems to perform well in both cases, the issue is visible with y = a*x + b*y)
    461  */
    462 template<typename Scalar>
    463 struct scalar_multiple_op {
    464   typedef typename packet_traits<Scalar>::type Packet;
    465   // FIXME default copy constructors seems bugged with std::complex<>
    466   EIGEN_STRONG_INLINE scalar_multiple_op(const scalar_multiple_op& other) : m_other(other.m_other) { }
    467   EIGEN_STRONG_INLINE scalar_multiple_op(const Scalar& other) : m_other(other) { }
    468   EIGEN_STRONG_INLINE Scalar operator() (const Scalar& a) const { return a * m_other; }
    469   EIGEN_STRONG_INLINE const Packet packetOp(const Packet& a) const
    470   { return internal::pmul(a, pset1<Packet>(m_other)); }
    471   typename add_const_on_value_type<typename NumTraits<Scalar>::Nested>::type m_other;
    472 };
    473 template<typename Scalar>
    474 struct functor_traits<scalar_multiple_op<Scalar> >
    475 { enum { Cost = NumTraits<Scalar>::MulCost, PacketAccess = packet_traits<Scalar>::HasMul }; };
    476 
    477 template<typename Scalar1, typename Scalar2>
    478 struct scalar_multiple2_op {
    479   typedef typename scalar_product_traits<Scalar1,Scalar2>::ReturnType result_type;
    480   EIGEN_STRONG_INLINE scalar_multiple2_op(const scalar_multiple2_op& other) : m_other(other.m_other) { }
    481   EIGEN_STRONG_INLINE scalar_multiple2_op(const Scalar2& other) : m_other(other) { }
    482   EIGEN_STRONG_INLINE result_type operator() (const Scalar1& a) const { return a * m_other; }
    483   typename add_const_on_value_type<typename NumTraits<Scalar2>::Nested>::type m_other;
    484 };
    485 template<typename Scalar1,typename Scalar2>
    486 struct functor_traits<scalar_multiple2_op<Scalar1,Scalar2> >
    487 { enum { Cost = NumTraits<Scalar1>::MulCost, PacketAccess = false }; };
    488 
    489 /** \internal
    490   * \brief Template functor to divide a scalar by a fixed other one
    491   *
    492   * This functor is used to implement the quotient of a matrix by
    493   * a scalar where the scalar type is not necessarily a floating point type.
    494   *
    495   * \sa class CwiseUnaryOp, MatrixBase::operator/
    496   */
    497 template<typename Scalar>
    498 struct scalar_quotient1_op {
    499   typedef typename packet_traits<Scalar>::type Packet;
    500   // FIXME default copy constructors seems bugged with std::complex<>
    501   EIGEN_STRONG_INLINE scalar_quotient1_op(const scalar_quotient1_op& other) : m_other(other.m_other) { }
    502   EIGEN_STRONG_INLINE scalar_quotient1_op(const Scalar& other) : m_other(other) {}
    503   EIGEN_STRONG_INLINE Scalar operator() (const Scalar& a) const { return a / m_other; }
    504   EIGEN_STRONG_INLINE const Packet packetOp(const Packet& a) const
    505   { return internal::pdiv(a, pset1<Packet>(m_other)); }
    506   typename add_const_on_value_type<typename NumTraits<Scalar>::Nested>::type m_other;
    507 };
    508 template<typename Scalar>
    509 struct functor_traits<scalar_quotient1_op<Scalar> >
    510 { enum { Cost = 2 * NumTraits<Scalar>::MulCost, PacketAccess = packet_traits<Scalar>::HasDiv }; };
    511 
    512 // nullary functors
    513 
    514 template<typename Scalar>
    515 struct scalar_constant_op {
    516   typedef typename packet_traits<Scalar>::type Packet;
    517   EIGEN_STRONG_INLINE scalar_constant_op(const scalar_constant_op& other) : m_other(other.m_other) { }
    518   EIGEN_STRONG_INLINE scalar_constant_op(const Scalar& other) : m_other(other) { }
    519   template<typename Index>
    520   EIGEN_STRONG_INLINE const Scalar operator() (Index, Index = 0) const { return m_other; }
    521   template<typename Index>
    522   EIGEN_STRONG_INLINE const Packet packetOp(Index, Index = 0) const { return internal::pset1<Packet>(m_other); }
    523   const Scalar m_other;
    524 };
    525 template<typename Scalar>
    526 struct functor_traits<scalar_constant_op<Scalar> >
    527 // FIXME replace this packet test by a safe one
    528 { enum { Cost = 1, PacketAccess = packet_traits<Scalar>::Vectorizable, IsRepeatable = true }; };
    529 
    530 template<typename Scalar> struct scalar_identity_op {
    531   EIGEN_EMPTY_STRUCT_CTOR(scalar_identity_op)
    532   template<typename Index>
    533   EIGEN_STRONG_INLINE const Scalar operator() (Index row, Index col) const { return row==col ? Scalar(1) : Scalar(0); }
    534 };
    535 template<typename Scalar>
    536 struct functor_traits<scalar_identity_op<Scalar> >
    537 { enum { Cost = NumTraits<Scalar>::AddCost, PacketAccess = false, IsRepeatable = true }; };
    538 
    539 template <typename Scalar, bool RandomAccess> struct linspaced_op_impl;
    540 
    541 // linear access for packet ops:
    542 // 1) initialization
    543 //   base = [low, ..., low] + ([step, ..., step] * [-size, ..., 0])
    544 // 2) each step (where size is 1 for coeff access or PacketSize for packet access)
    545 //   base += [size*step, ..., size*step]
    546 //
    547 // TODO: Perhaps it's better to initialize lazily (so not in the constructor but in packetOp)
    548 //       in order to avoid the padd() in operator() ?
    549 template <typename Scalar>
    550 struct linspaced_op_impl<Scalar,false>
    551 {
    552   typedef typename packet_traits<Scalar>::type Packet;
    553 
    554   linspaced_op_impl(const Scalar& low, const Scalar& step) :
    555   m_low(low), m_step(step),
    556   m_packetStep(pset1<Packet>(packet_traits<Scalar>::size*step)),
    557   m_base(padd(pset1<Packet>(low), pmul(pset1<Packet>(step),plset<Scalar>(-packet_traits<Scalar>::size)))) {}
    558 
    559   template<typename Index>
    560   EIGEN_STRONG_INLINE const Scalar operator() (Index i) const
    561   {
    562     m_base = padd(m_base, pset1<Packet>(m_step));
    563     return m_low+Scalar(i)*m_step;
    564   }
    565 
    566   template<typename Index>
    567   EIGEN_STRONG_INLINE const Packet packetOp(Index) const { return m_base = padd(m_base,m_packetStep); }
    568 
    569   const Scalar m_low;
    570   const Scalar m_step;
    571   const Packet m_packetStep;
    572   mutable Packet m_base;
    573 };
    574 
    575 // random access for packet ops:
    576 // 1) each step
    577 //   [low, ..., low] + ( [step, ..., step] * ( [i, ..., i] + [0, ..., size] ) )
    578 template <typename Scalar>
    579 struct linspaced_op_impl<Scalar,true>
    580 {
    581   typedef typename packet_traits<Scalar>::type Packet;
    582 
    583   linspaced_op_impl(const Scalar& low, const Scalar& step) :
    584   m_low(low), m_step(step),
    585   m_lowPacket(pset1<Packet>(m_low)), m_stepPacket(pset1<Packet>(m_step)), m_interPacket(plset<Scalar>(0)) {}
    586 
    587   template<typename Index>
    588   EIGEN_STRONG_INLINE const Scalar operator() (Index i) const { return m_low+i*m_step; }
    589 
    590   template<typename Index>
    591   EIGEN_STRONG_INLINE const Packet packetOp(Index i) const
    592   { return internal::padd(m_lowPacket, pmul(m_stepPacket, padd(pset1<Packet>(Scalar(i)),m_interPacket))); }
    593 
    594   const Scalar m_low;
    595   const Scalar m_step;
    596   const Packet m_lowPacket;
    597   const Packet m_stepPacket;
    598   const Packet m_interPacket;
    599 };
    600 
    601 // ----- Linspace functor ----------------------------------------------------------------
    602 
    603 // Forward declaration (we default to random access which does not really give
    604 // us a speed gain when using packet access but it allows to use the functor in
    605 // nested expressions).
    606 template <typename Scalar, bool RandomAccess = true> struct linspaced_op;
    607 template <typename Scalar, bool RandomAccess> struct functor_traits< linspaced_op<Scalar,RandomAccess> >
    608 { enum { Cost = 1, PacketAccess = packet_traits<Scalar>::HasSetLinear, IsRepeatable = true }; };
    609 template <typename Scalar, bool RandomAccess> struct linspaced_op
    610 {
    611   typedef typename packet_traits<Scalar>::type Packet;
    612   linspaced_op(const Scalar& low, const Scalar& high, DenseIndex num_steps) : impl((num_steps==1 ? high : low), (num_steps==1 ? Scalar() : (high-low)/Scalar(num_steps-1))) {}
    613 
    614   template<typename Index>
    615   EIGEN_STRONG_INLINE const Scalar operator() (Index i) const { return impl(i); }
    616 
    617   // We need this function when assigning e.g. a RowVectorXd to a MatrixXd since
    618   // there row==0 and col is used for the actual iteration.
    619   template<typename Index>
    620   EIGEN_STRONG_INLINE const Scalar operator() (Index row, Index col) const
    621   {
    622     eigen_assert(col==0 || row==0);
    623     return impl(col + row);
    624   }
    625 
    626   template<typename Index>
    627   EIGEN_STRONG_INLINE const Packet packetOp(Index i) const { return impl.packetOp(i); }
    628 
    629   // We need this function when assigning e.g. a RowVectorXd to a MatrixXd since
    630   // there row==0 and col is used for the actual iteration.
    631   template<typename Index>
    632   EIGEN_STRONG_INLINE const Packet packetOp(Index row, Index col) const
    633   {
    634     eigen_assert(col==0 || row==0);
    635     return impl.packetOp(col + row);
    636   }
    637 
    638   // This proxy object handles the actual required temporaries, the different
    639   // implementations (random vs. sequential access) as well as the
    640   // correct piping to size 2/4 packet operations.
    641   const linspaced_op_impl<Scalar,RandomAccess> impl;
    642 };
    643 
    644 // all functors allow linear access, except scalar_identity_op. So we fix here a quick meta
    645 // to indicate whether a functor allows linear access, just always answering 'yes' except for
    646 // scalar_identity_op.
    647 // FIXME move this to functor_traits adding a functor_default
    648 template<typename Functor> struct functor_has_linear_access { enum { ret = 1 }; };
    649 template<typename Scalar> struct functor_has_linear_access<scalar_identity_op<Scalar> > { enum { ret = 0 }; };
    650 
    651 // In Eigen, any binary op (Product, CwiseBinaryOp) require the Lhs and Rhs to have the same scalar type, except for multiplication
    652 // where the mixing of different types is handled by scalar_product_traits
    653 // In particular, real * complex<real> is allowed.
    654 // FIXME move this to functor_traits adding a functor_default
    655 template<typename Functor> struct functor_is_product_like { enum { ret = 0 }; };
    656 template<typename LhsScalar,typename RhsScalar> struct functor_is_product_like<scalar_product_op<LhsScalar,RhsScalar> > { enum { ret = 1 }; };
    657 template<typename LhsScalar,typename RhsScalar> struct functor_is_product_like<scalar_conj_product_op<LhsScalar,RhsScalar> > { enum { ret = 1 }; };
    658 template<typename LhsScalar,typename RhsScalar> struct functor_is_product_like<scalar_quotient_op<LhsScalar,RhsScalar> > { enum { ret = 1 }; };
    659 
    660 
    661 /** \internal
    662   * \brief Template functor to add a scalar to a fixed other one
    663   * \sa class CwiseUnaryOp, Array::operator+
    664   */
    665 /* If you wonder why doing the pset1() in packetOp() is an optimization check scalar_multiple_op */
    666 template<typename Scalar>
    667 struct scalar_add_op {
    668   typedef typename packet_traits<Scalar>::type Packet;
    669   // FIXME default copy constructors seems bugged with std::complex<>
    670   inline scalar_add_op(const scalar_add_op& other) : m_other(other.m_other) { }
    671   inline scalar_add_op(const Scalar& other) : m_other(other) { }
    672   inline Scalar operator() (const Scalar& a) const { return a + m_other; }
    673   inline const Packet packetOp(const Packet& a) const
    674   { return internal::padd(a, pset1<Packet>(m_other)); }
    675   const Scalar m_other;
    676 };
    677 template<typename Scalar>
    678 struct functor_traits<scalar_add_op<Scalar> >
    679 { enum { Cost = NumTraits<Scalar>::AddCost, PacketAccess = packet_traits<Scalar>::HasAdd }; };
    680 
    681 /** \internal
    682   * \brief Template functor to compute the square root of a scalar
    683   * \sa class CwiseUnaryOp, Cwise::sqrt()
    684   */
    685 template<typename Scalar> struct scalar_sqrt_op {
    686   EIGEN_EMPTY_STRUCT_CTOR(scalar_sqrt_op)
    687   inline const Scalar operator() (const Scalar& a) const { using std::sqrt; return sqrt(a); }
    688   typedef typename packet_traits<Scalar>::type Packet;
    689   inline Packet packetOp(const Packet& a) const { return internal::psqrt(a); }
    690 };
    691 template<typename Scalar>
    692 struct functor_traits<scalar_sqrt_op<Scalar> >
    693 { enum {
    694     Cost = 5 * NumTraits<Scalar>::MulCost,
    695     PacketAccess = packet_traits<Scalar>::HasSqrt
    696   };
    697 };
    698 
    699 /** \internal
    700   * \brief Template functor to compute the cosine of a scalar
    701   * \sa class CwiseUnaryOp, ArrayBase::cos()
    702   */
    703 template<typename Scalar> struct scalar_cos_op {
    704   EIGEN_EMPTY_STRUCT_CTOR(scalar_cos_op)
    705   inline Scalar operator() (const Scalar& a) const { using std::cos; return cos(a); }
    706   typedef typename packet_traits<Scalar>::type Packet;
    707   inline Packet packetOp(const Packet& a) const { return internal::pcos(a); }
    708 };
    709 template<typename Scalar>
    710 struct functor_traits<scalar_cos_op<Scalar> >
    711 {
    712   enum {
    713     Cost = 5 * NumTraits<Scalar>::MulCost,
    714     PacketAccess = packet_traits<Scalar>::HasCos
    715   };
    716 };
    717 
    718 /** \internal
    719   * \brief Template functor to compute the sine of a scalar
    720   * \sa class CwiseUnaryOp, ArrayBase::sin()
    721   */
    722 template<typename Scalar> struct scalar_sin_op {
    723   EIGEN_EMPTY_STRUCT_CTOR(scalar_sin_op)
    724   inline const Scalar operator() (const Scalar& a) const { using std::sin; return sin(a); }
    725   typedef typename packet_traits<Scalar>::type Packet;
    726   inline Packet packetOp(const Packet& a) const { return internal::psin(a); }
    727 };
    728 template<typename Scalar>
    729 struct functor_traits<scalar_sin_op<Scalar> >
    730 {
    731   enum {
    732     Cost = 5 * NumTraits<Scalar>::MulCost,
    733     PacketAccess = packet_traits<Scalar>::HasSin
    734   };
    735 };
    736 
    737 
    738 /** \internal
    739   * \brief Template functor to compute the tan of a scalar
    740   * \sa class CwiseUnaryOp, ArrayBase::tan()
    741   */
    742 template<typename Scalar> struct scalar_tan_op {
    743   EIGEN_EMPTY_STRUCT_CTOR(scalar_tan_op)
    744   inline const Scalar operator() (const Scalar& a) const { using std::tan; return tan(a); }
    745   typedef typename packet_traits<Scalar>::type Packet;
    746   inline Packet packetOp(const Packet& a) const { return internal::ptan(a); }
    747 };
    748 template<typename Scalar>
    749 struct functor_traits<scalar_tan_op<Scalar> >
    750 {
    751   enum {
    752     Cost = 5 * NumTraits<Scalar>::MulCost,
    753     PacketAccess = packet_traits<Scalar>::HasTan
    754   };
    755 };
    756 
    757 /** \internal
    758   * \brief Template functor to compute the arc cosine of a scalar
    759   * \sa class CwiseUnaryOp, ArrayBase::acos()
    760   */
    761 template<typename Scalar> struct scalar_acos_op {
    762   EIGEN_EMPTY_STRUCT_CTOR(scalar_acos_op)
    763   inline const Scalar operator() (const Scalar& a) const { using std::acos; return acos(a); }
    764   typedef typename packet_traits<Scalar>::type Packet;
    765   inline Packet packetOp(const Packet& a) const { return internal::pacos(a); }
    766 };
    767 template<typename Scalar>
    768 struct functor_traits<scalar_acos_op<Scalar> >
    769 {
    770   enum {
    771     Cost = 5 * NumTraits<Scalar>::MulCost,
    772     PacketAccess = packet_traits<Scalar>::HasACos
    773   };
    774 };
    775 
    776 /** \internal
    777   * \brief Template functor to compute the arc sine of a scalar
    778   * \sa class CwiseUnaryOp, ArrayBase::asin()
    779   */
    780 template<typename Scalar> struct scalar_asin_op {
    781   EIGEN_EMPTY_STRUCT_CTOR(scalar_asin_op)
    782   inline const Scalar operator() (const Scalar& a) const { using std::asin; return asin(a); }
    783   typedef typename packet_traits<Scalar>::type Packet;
    784   inline Packet packetOp(const Packet& a) const { return internal::pasin(a); }
    785 };
    786 template<typename Scalar>
    787 struct functor_traits<scalar_asin_op<Scalar> >
    788 {
    789   enum {
    790     Cost = 5 * NumTraits<Scalar>::MulCost,
    791     PacketAccess = packet_traits<Scalar>::HasASin
    792   };
    793 };
    794 
    795 /** \internal
    796   * \brief Template functor to raise a scalar to a power
    797   * \sa class CwiseUnaryOp, Cwise::pow
    798   */
    799 template<typename Scalar>
    800 struct scalar_pow_op {
    801   // FIXME default copy constructors seems bugged with std::complex<>
    802   inline scalar_pow_op(const scalar_pow_op& other) : m_exponent(other.m_exponent) { }
    803   inline scalar_pow_op(const Scalar& exponent) : m_exponent(exponent) {}
    804   inline Scalar operator() (const Scalar& a) const { return numext::pow(a, m_exponent); }
    805   const Scalar m_exponent;
    806 };
    807 template<typename Scalar>
    808 struct functor_traits<scalar_pow_op<Scalar> >
    809 { enum { Cost = 5 * NumTraits<Scalar>::MulCost, PacketAccess = false }; };
    810 
    811 /** \internal
    812   * \brief Template functor to compute the quotient between a scalar and array entries.
    813   * \sa class CwiseUnaryOp, Cwise::inverse()
    814   */
    815 template<typename Scalar>
    816 struct scalar_inverse_mult_op {
    817   scalar_inverse_mult_op(const Scalar& other) : m_other(other) {}
    818   inline Scalar operator() (const Scalar& a) const { return m_other / a; }
    819   template<typename Packet>
    820   inline const Packet packetOp(const Packet& a) const
    821   { return internal::pdiv(pset1<Packet>(m_other),a); }
    822   Scalar m_other;
    823 };
    824 
    825 /** \internal
    826   * \brief Template functor to compute the inverse of a scalar
    827   * \sa class CwiseUnaryOp, Cwise::inverse()
    828   */
    829 template<typename Scalar>
    830 struct scalar_inverse_op {
    831   EIGEN_EMPTY_STRUCT_CTOR(scalar_inverse_op)
    832   inline Scalar operator() (const Scalar& a) const { return Scalar(1)/a; }
    833   template<typename Packet>
    834   inline const Packet packetOp(const Packet& a) const
    835   { return internal::pdiv(pset1<Packet>(Scalar(1)),a); }
    836 };
    837 template<typename Scalar>
    838 struct functor_traits<scalar_inverse_op<Scalar> >
    839 { enum { Cost = NumTraits<Scalar>::MulCost, PacketAccess = packet_traits<Scalar>::HasDiv }; };
    840 
    841 /** \internal
    842   * \brief Template functor to compute the square of a scalar
    843   * \sa class CwiseUnaryOp, Cwise::square()
    844   */
    845 template<typename Scalar>
    846 struct scalar_square_op {
    847   EIGEN_EMPTY_STRUCT_CTOR(scalar_square_op)
    848   inline Scalar operator() (const Scalar& a) const { return a*a; }
    849   template<typename Packet>
    850   inline const Packet packetOp(const Packet& a) const
    851   { return internal::pmul(a,a); }
    852 };
    853 template<typename Scalar>
    854 struct functor_traits<scalar_square_op<Scalar> >
    855 { enum { Cost = NumTraits<Scalar>::MulCost, PacketAccess = packet_traits<Scalar>::HasMul }; };
    856 
    857 /** \internal
    858   * \brief Template functor to compute the cube of a scalar
    859   * \sa class CwiseUnaryOp, Cwise::cube()
    860   */
    861 template<typename Scalar>
    862 struct scalar_cube_op {
    863   EIGEN_EMPTY_STRUCT_CTOR(scalar_cube_op)
    864   inline Scalar operator() (const Scalar& a) const { return a*a*a; }
    865   template<typename Packet>
    866   inline const Packet packetOp(const Packet& a) const
    867   { return internal::pmul(a,pmul(a,a)); }
    868 };
    869 template<typename Scalar>
    870 struct functor_traits<scalar_cube_op<Scalar> >
    871 { enum { Cost = 2*NumTraits<Scalar>::MulCost, PacketAccess = packet_traits<Scalar>::HasMul }; };
    872 
    873 // default functor traits for STL functors:
    874 
    875 template<typename T>
    876 struct functor_traits<std::multiplies<T> >
    877 { enum { Cost = NumTraits<T>::MulCost, PacketAccess = false }; };
    878 
    879 template<typename T>
    880 struct functor_traits<std::divides<T> >
    881 { enum { Cost = NumTraits<T>::MulCost, PacketAccess = false }; };
    882 
    883 template<typename T>
    884 struct functor_traits<std::plus<T> >
    885 { enum { Cost = NumTraits<T>::AddCost, PacketAccess = false }; };
    886 
    887 template<typename T>
    888 struct functor_traits<std::minus<T> >
    889 { enum { Cost = NumTraits<T>::AddCost, PacketAccess = false }; };
    890 
    891 template<typename T>
    892 struct functor_traits<std::negate<T> >
    893 { enum { Cost = NumTraits<T>::AddCost, PacketAccess = false }; };
    894 
    895 template<typename T>
    896 struct functor_traits<std::logical_or<T> >
    897 { enum { Cost = 1, PacketAccess = false }; };
    898 
    899 template<typename T>
    900 struct functor_traits<std::logical_and<T> >
    901 { enum { Cost = 1, PacketAccess = false }; };
    902 
    903 template<typename T>
    904 struct functor_traits<std::logical_not<T> >
    905 { enum { Cost = 1, PacketAccess = false }; };
    906 
    907 template<typename T>
    908 struct functor_traits<std::greater<T> >
    909 { enum { Cost = 1, PacketAccess = false }; };
    910 
    911 template<typename T>
    912 struct functor_traits<std::less<T> >
    913 { enum { Cost = 1, PacketAccess = false }; };
    914 
    915 template<typename T>
    916 struct functor_traits<std::greater_equal<T> >
    917 { enum { Cost = 1, PacketAccess = false }; };
    918 
    919 template<typename T>
    920 struct functor_traits<std::less_equal<T> >
    921 { enum { Cost = 1, PacketAccess = false }; };
    922 
    923 template<typename T>
    924 struct functor_traits<std::equal_to<T> >
    925 { enum { Cost = 1, PacketAccess = false }; };
    926 
    927 template<typename T>
    928 struct functor_traits<std::not_equal_to<T> >
    929 { enum { Cost = 1, PacketAccess = false }; };
    930 
    931 template<typename T>
    932 struct functor_traits<std::binder2nd<T> >
    933 { enum { Cost = functor_traits<T>::Cost, PacketAccess = false }; };
    934 
    935 template<typename T>
    936 struct functor_traits<std::binder1st<T> >
    937 { enum { Cost = functor_traits<T>::Cost, PacketAccess = false }; };
    938 
    939 template<typename T>
    940 struct functor_traits<std::unary_negate<T> >
    941 { enum { Cost = 1 + functor_traits<T>::Cost, PacketAccess = false }; };
    942 
    943 template<typename T>
    944 struct functor_traits<std::binary_negate<T> >
    945 { enum { Cost = 1 + functor_traits<T>::Cost, PacketAccess = false }; };
    946 
    947 #ifdef EIGEN_STDEXT_SUPPORT
    948 
    949 template<typename T0,typename T1>
    950 struct functor_traits<std::project1st<T0,T1> >
    951 { enum { Cost = 0, PacketAccess = false }; };
    952 
    953 template<typename T0,typename T1>
    954 struct functor_traits<std::project2nd<T0,T1> >
    955 { enum { Cost = 0, PacketAccess = false }; };
    956 
    957 template<typename T0,typename T1>
    958 struct functor_traits<std::select2nd<std::pair<T0,T1> > >
    959 { enum { Cost = 0, PacketAccess = false }; };
    960 
    961 template<typename T0,typename T1>
    962 struct functor_traits<std::select1st<std::pair<T0,T1> > >
    963 { enum { Cost = 0, PacketAccess = false }; };
    964 
    965 template<typename T0,typename T1>
    966 struct functor_traits<std::unary_compose<T0,T1> >
    967 { enum { Cost = functor_traits<T0>::Cost + functor_traits<T1>::Cost, PacketAccess = false }; };
    968 
    969 template<typename T0,typename T1,typename T2>
    970 struct functor_traits<std::binary_compose<T0,T1,T2> >
    971 { enum { Cost = functor_traits<T0>::Cost + functor_traits<T1>::Cost + functor_traits<T2>::Cost, PacketAccess = false }; };
    972 
    973 #endif // EIGEN_STDEXT_SUPPORT
    974 
    975 // allow to add new functors and specializations of functor_traits from outside Eigen.
    976 // this macro is really needed because functor_traits must be specialized after it is declared but before it is used...
    977 #ifdef EIGEN_FUNCTORS_PLUGIN
    978 #include EIGEN_FUNCTORS_PLUGIN
    979 #endif
    980 
    981 } // end namespace internal
    982 
    983 } // end namespace Eigen
    984 
    985 #endif // EIGEN_FUNCTORS_H
    986