Patch attached. Hopefully clarifies the documentation. I will see if I will find time to work on better solution.
Petr On Tue, 2019-11-05 at 17:11 +0100, Christoph Hertzberg wrote: > Adding an additional template method would not be a problem (minimal > example: https://godbolt.org/z/ZB9WDI), though I don't see any real > advantage vs writing > > m.cast<int>().sum(); > > Sure, `m.sum<int>();` would be shorter, but also more confusing. > > Your "another template parameter for the matrices" formulation > sounded > to me like you wanted to add another template parameter to > `Eigen::Matrix` (which is not an option). > > Any patches clarifying the documentation are welcome! > > Christoph > > > On 05/11/2019 16.44, Petr Kubánek wrote: > > Hi, > > > > I know what's the problem. I am just looking either to document it > > or > > for a solution. > > > > Why would: > > > > template <typename dt> dt sum() > > > > break API/ABI? You will be able to call either sum() or > > sum<int32_t>(), > > shouldn't you? I will try that and submit a patch. > > > > Having sum<dt>() working, one can hopefully create mean<dt_sum, > > dt_mean>(), so one can code: > > > > mean<int32_t,double>() > > > > Petr > > > > On Tue, 2019-11-05 at 16:36 +0100, Christoph Hertzberg wrote: > > > On 05/11/2019 16.11, Peter wrote: > > > > [...] > > > > > > > > actually, this would also be interesting for the scalar > > > > products > > > > in > > > > general, namely a different > > > > type for accumulating the sums within a scalar product, e.g. as > > > > yet > > > > another template parameter for the matrices. > > > > > > Adding another template parameter to basic types is not an > > > option. > > > This > > > would break ABI and API compatibility (even if the parameter has > > > a > > > default). > > > > > > You could create your own custom type `my_int16` for which > > > `my_int16 > > > * > > > my_int16` results in a `my_int32` (this needs to be told to > > > Eigen, > > > similar to how real*complex products are handled). > > > > > > > > [...] > > > > > > > > I think it's more subtle than that. > > > > Even > > > > > > > > int16_t Two = 2; > > > > int16_t Max = INT16_MAX; > > > > int16_t mean = ( Max/2 + Max + Two + Two ) / int16_t(4); > > > > > > > > doesn't produce an overflow. > > > > > > Yes, because `Max/2` gets implicitly converted to `int`. > > > Actually, > > > even > > > adding two `int16_t` get implicitly converted to `int` (search > > > for > > > integer promotion rules -- I don't know them entirely either). > > > And > > > dividing an `int` by an `int16_t` results in an `int` which is > > > only > > > afterwards converted to `int16_t`. > > > > > > In contrast, Eigen::DenseBase::sum() does something (more or > > > less) > > > equivalent to: > > > > > > T sum = T(0); > > > for(Index i=0; i< size(); ++i) > > > sum+=coeff(i); > > > return sum; > > > > > > i.e., after each addition the result gets reduced to the scalar > > > type > > > of > > > the matrix, thus it will immediately overflow. > > > > > > And mean() essentially just takes the return value of `sum()` and > > > divides it by the size. > > > > > > > > > Christoph > > > > > > > > > > > > > >
# HG changeset patch # User Petr Kubánek <[email protected]> # Date 1573065305 25200 # Wed Nov 06 11:35:05 2019 -0700 # Node ID 42bca2d654395d5356bc053f03f5e7e558273cf7 # Parent afc120bc03bdc4265858d6f86218eb1fed51b1b9 improves reductions documentation diff --git a/Eigen/src/Core/Redux.h b/Eigen/src/Core/Redux.h --- a/Eigen/src/Core/Redux.h +++ b/Eigen/src/Core/Redux.h @@ -436,65 +436,69 @@ DenseBase<Derived>::minCoeff() const template<typename Derived> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE typename internal::traits<Derived>::Scalar DenseBase<Derived>::maxCoeff() const { return derived().redux(Eigen::internal::scalar_max_op<Scalar,Scalar>()); } /** \returns the sum of all coefficients of \c *this + * \warning internal accumulator type and returned type matches base type. Overflow might occur for e.g. int values. Use \link MatrixBase::cast() .cast() \endlink to cast your variable to proper type. * * If \c *this is empty, then the value 0 is returned. * * \sa trace(), prod(), mean() */ template<typename Derived> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE typename internal::traits<Derived>::Scalar DenseBase<Derived>::sum() const { if(SizeAtCompileTime==0 || (SizeAtCompileTime==Dynamic && size()==0)) return Scalar(0); return derived().redux(Eigen::internal::scalar_sum_op<Scalar,Scalar>()); } /** \returns the mean of all coefficients of *this -* -* \sa trace(), prod(), sum() -*/ + * \warning internal accumulator type and returned type matches base type. Overflow might occur for e.g. int values. Use \link MatrixBase::cast() .cast() \endlink to cast your variable to proper type. + * + * \sa trace(), prod(), sum() + */ template<typename Derived> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE typename internal::traits<Derived>::Scalar DenseBase<Derived>::mean() const { #ifdef __INTEL_COMPILER #pragma warning push #pragma warning ( disable : 2259 ) #endif return Scalar(derived().redux(Eigen::internal::scalar_sum_op<Scalar,Scalar>())) / Scalar(this->size()); #ifdef __INTEL_COMPILER #pragma warning pop #endif } /** \returns the product of all coefficients of *this + * \warning internal accumulator type and returned type matches base type. Overflow might occur for e.g. int values. Use \link MatrixBase::cast() .cast() \endlink to cast your variable to proper type. * * Example: \include MatrixBase_prod.cpp * Output: \verbinclude MatrixBase_prod.out * * \sa sum(), mean(), trace() */ template<typename Derived> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE typename internal::traits<Derived>::Scalar DenseBase<Derived>::prod() const { if(SizeAtCompileTime==0 || (SizeAtCompileTime==Dynamic && size()==0)) return Scalar(1); return derived().redux(Eigen::internal::scalar_product_op<Scalar>()); } /** \returns the trace of \c *this, i.e. the sum of the coefficients on the main diagonal. + * \warning internal accumulator type and returned type matches base type. Overflow might occur for e.g. int values. Use \link MatrixBase::cast() .cast() \endlink to cast your variable to proper type. * * \c *this can be any matrix, not necessarily square. * * \sa diagonal(), sum() */ template<typename Derived> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE typename internal::traits<Derived>::Scalar MatrixBase<Derived>::trace() const diff --git a/doc/TutorialReductionsVisitorsBroadcasting.dox b/doc/TutorialReductionsVisitorsBroadcasting.dox --- a/doc/TutorialReductionsVisitorsBroadcasting.dox +++ b/doc/TutorialReductionsVisitorsBroadcasting.dox @@ -16,16 +16,18 @@ returning the sum of all the coefficient <tr><th>Example:</th><th>Output:</th></tr> <tr><td> \include tut_arithmetic_redux_basic.cpp </td> <td> \verbinclude tut_arithmetic_redux_basic.out </td></tr></table> +Please note \link DenseBase::sum() .sum() \endlink calculates sum using same storage type as the matrix or array. + The \em trace of a matrix, as returned by the function \c trace(), is the sum of the diagonal coefficients and can equivalently be computed <tt>a.diagonal().sum()</tt>. \subsection TutorialReductionsVisitorsBroadcastingReductionsNorm Norm computations The (Euclidean a.k.a. \f$\ell^2\f$) squared norm of a vector can be obtained \link MatrixBase::squaredNorm() squaredNorm() \endlink. It is equal to the dot product of the vector by itself, and equivalently to the sum of squared absolute values of its coefficients. Eigen also provides the \link MatrixBase::norm() norm() \endlink method, which returns the square root of \link MatrixBase::squaredNorm() squaredNorm() \endlink.
