    template<typename Dst, typename Src> void operator()(const Dst& dst, const Src& src) const {
        dst.const_cast_derived() += m_scale * src;
      }
    };
    
    template<typename Dest>
    inline void evalTo(Dest& dest) const {
      internal::outer_product_selector_run(*this, dest, set(), IsRowMajor<Dest>());
    }
    
    template<typename Dest>
    inline void addTo(Dest& dest) const {
      internal::outer_product_selector_run(*this, dest, add(), IsRowMajor<Dest>());
    }

    template<typename Dest>
    inline void subTo(Dest& dest) const {
      internal::outer_product_selector_run(*this, dest, sub(), IsRowMajor<Dest>());
    }

    template<typename Dest> void scaleAndAddTo(Dest& dest, const Scalar& alpha) const
    {
      internal::outer_product_selector_run(*this, dest, adds(alpha), IsRowMajor<Dest>());
    }
};

/***********************************************************************
*  Implementation of General Matrix Vector Product
***********************************************************************/

/*  According to the shape/flags of the matrix we have to distinghish 3 different cases:
 *   1 - the matrix is col-major, BLAS compatible and M is large => call fast BLAS-like colmajor routine
 *   2 - the matrix is row-major, BLAS compatible and N is large => call fast BLAS-like rowmajor routine
 *   3 - all other cases are handled using a simple loop along the outer-storage direction.
 *  Therefore we need a lower level meta selector.
 *  Furthermore, if the matrix is the rhs, then the product has to be transposed.
 */
namespace internal {

template<typename Lhs, typename Rhs>
struct traits<GeneralProduct<Lhs,Rhs,GemvProduct> >
 : traits<ProductBase<GeneralProduct<Lhs,Rhs,GemvProduct>, Lhs, Rhs> >
{};

template<int Side, int StorageOrder, bool BlasCompatible>
struct gemv_selector;

} // end namespace internal

template<typename Lhs, typename Rhs>
class GeneralProduct<Lhs, Rhs, GemvProduct>
  : public ProductBase<GeneralProduct<Lhs,Rhs,GemvProduct>, Lhs, Rhs>
{
  public:
    EIGEN_PROD