Re: [eigen] log10 support?

On 7 June 2010 18:06, Benoit Jacob <jacob.benoit.1@xxxxxxxxx> wrote:

2010/6/7 Benoit Jacob <jacob..benoit.1@xxxxxxxxx>:

> ok so time for another answer :)
>
> The only thing that I really don't want is for us to add
> SSE/AltiVec/whatever code for log10. We must reuse the existing code
> for log there.
>
> Then there are two ways:
> a) the way of your patch: add a new custom functor for log10
> b) the way i suggested: just let log10 return a product _expression_.
>
> Your way a) has the advantage of sparing one multiplication (at least
> in the non vectorized case), but b) means less extra code.
>
> Now Hauke raises a good point, that I overlooked: even with a) you can
> still reuse our vectorization code for log10. So actually, I'm OK with
> your approach, provided that we do that. You can do it generically for
> all types at once (Hauke showed the way for float but there's no
> reason to special-case it): edit Core/GenericPacketMath.h and
> implement ei_plog10 there like the other functions, returning
> ei_pmul(
> ei_pset1(log(Scalar(0.1))),

woooooops!!

replace log(Scalar(0.1)) by Scalar(1) / log(Scalar(10))

and sorry Hauke, I saw too late that you already wrote that. Though I
guess there's value in insisting that it's faster to multiply by a
constant, than to divide by a constant.

Benoit

> ei_plog(x)
> )
> The next problem is to cache the result of ei_pset1(log(Scalar(0.1))),
> probably a static constant, maybe have a look at what's done in
> ei_plog itself in SSE/. If that sounds obscure, leave that for later.
>
> In any case, we don't need a new HasLog10. Since we'll be using
> ei_plog, all we need is the existing HasLog.
>
> Benoit
>
> 2010/6/7 Trevor Irons <trevorirons@xxxxxxxxx>:
>> OK, thanks. I'm on a really old machine and can't check sse here. I'll test
>> your code and generate another patch. Unless I hear that is is unwanted.
>>
>> (This just goes in eigenNes/Eigen/src/Core/arch/SSE/MathFunctions.h) right?
>>
>> -Trevor
>>
>> On 7 June 2010 15:46, Hauke Heibel <hauke.heibel@xxxxxxxxxxxxxx> wrote:
>>>
>>> Let me try to answer for Benoit. :)
>>>
>>> On Mon, Jun 7, 2010 at 11:20 PM, Trevor Irons <trevorirons@xxxxxxxxx>
>>> wrote:
>>> > Sure, I can do this. I prefer store a minimum of constants. But if it
>>> > will
>>> > be faster, then sure.
>>> >
>>> > Why won't calls to std::log10 be vectorized? No sse instructions?
>>>
>>> This is because std::log is not right away available as an SSE
>>> intrinsic and we have a special version based on the code from here:
>>> http://gruntthepeon.free.fr/ssemath/
>>>
>>> Nonetheless I think your patch is the right way to go. We could
>>> achieve easily what Benoit suggests by offering a version of ei_plog10
>>> such as
>>>
>>> Packet4f ei_plog10(Packet4f x)
>>> {
>>> const Packet4f ei_p4f_log10 = ei_plog(ei_pset1(10.0f));
>>> return ei_pdiv(ei_plog(x), ei_p4f_log10);
>>> }
>>>
>>> which uses "Packet4f ei_plog(Packet4f x)". The constant would probably
>>> be taken care of by a decent compiler, i.e. I don't think we need
>>> manual caching.
>>>
>>> > Also, is there any desire to have log10 defined your way in Eigen, for
>>> > convenience? Or are you just suggesting I do this locally?
>>>
>>> I am pretty sure he meant that it would be nice to really have this in
>>> Eigen and not just in a local copy... though I am wildly guessing.
>>>
>>> - Hauke
>>>
>>> p.s. I am pretty tired an hope that I picked the right conversion factor
>>> ...
>>>
>>>
>>
>>
>

# HG changeset patch # User Trevor Irons <trevorirons@xxxxxxxxx> # Date 1275975037 21600 # Node ID f4d712b82883c12dbaed74bebd9258c78856ad3b # Parent a9d324abb4f1f08aff26daacc216904476f5342a Added log10 convenience support. For SSE instructions scaled versions of log are called. diff -r a9d324abb4f1 -r f4d712b82883 Eigen/src/Array/GlobalFunctions.h --- a/Eigen/src/Array/GlobalFunctions.h Tue Jun 08 00:41:33 2010 +0200 +++ b/Eigen/src/Array/GlobalFunctions.h Mon Jun 07 23:30:37 2010 -0600 @@ -58,6 +58,7 @@ EIGEN_ARRAY_DECLARE_GLOBAL_STD_UNARY(cos,ei_scalar_cos_op) EIGEN_ARRAY_DECLARE_GLOBAL_STD_UNARY(exp,ei_scalar_exp_op) EIGEN_ARRAY_DECLARE_GLOBAL_STD_UNARY(log,ei_scalar_log_op) + EIGEN_ARRAY_DECLARE_GLOBAL_STD_UNARY(log10,ei_scalar_log10_op) EIGEN_ARRAY_DECLARE_GLOBAL_STD_UNARY(abs,ei_scalar_abs_op) EIGEN_ARRAY_DECLARE_GLOBAL_STD_UNARY(sqrt,ei_scalar_sqrt_op) } @@ -70,6 +71,7 @@ EIGEN_ARRAY_DECLARE_GLOBAL_EIGEN_UNARY(ei_cos,ei_scalar_cos_op) EIGEN_ARRAY_DECLARE_GLOBAL_EIGEN_UNARY(ei_exp,ei_scalar_exp_op) EIGEN_ARRAY_DECLARE_GLOBAL_EIGEN_UNARY(ei_log,ei_scalar_log_op) + EIGEN_ARRAY_DECLARE_GLOBAL_EIGEN_UNARY(ei_log10,ei_scalar_log10_op) EIGEN_ARRAY_DECLARE_GLOBAL_EIGEN_UNARY(ei_abs,ei_scalar_abs_op) EIGEN_ARRAY_DECLARE_GLOBAL_EIGEN_UNARY(ei_abs2,ei_scalar_abs2_op) EIGEN_ARRAY_DECLARE_GLOBAL_EIGEN_UNARY(ei_sqrt,ei_scalar_sqrt_op) diff -r a9d324abb4f1 -r f4d712b82883 Eigen/src/Core/Functors.h --- a/Eigen/src/Core/Functors.h Tue Jun 08 00:41:33 2010 +0200 +++ b/Eigen/src/Core/Functors.h Mon Jun 07 23:30:37 2010 -0600 @@ -373,6 +373,22 @@ { enum { Cost = 5 * NumTraits<Scalar>::MulCost, PacketAccess = ei_packet_traits<Scalar>::HasLog }; }; /** \internal + * + * \brief Template functor to compute the base 10 logarithm of a scalar + * + * \sa class CwiseUnaryOp, Cwise::log10() + */ +template<typename Scalar> struct ei_scalar_log10_op { + EIGEN_EMPTY_STRUCT_CTOR(ei_scalar_log10_op) + inline const Scalar operator() (const Scalar& a) const { return ei_log10(a); } + typedef typename ei_packet_traits<Scalar>::type Packet; + inline Packet packetOp(const Packet& a) const { return ei_plog10(a); } +}; +template<typename Scalar> +struct ei_functor_traits<ei_scalar_log10_op<Scalar> > +{ enum { Cost = 6 * NumTraits<Scalar>::MulCost, PacketAccess = ei_packet_traits<Scalar>::HasLog }; }; + +/** \internal * \brief Template functor to multiply a scalar by a fixed other one * * \sa class CwiseUnaryOp, MatrixBase::operator*, MatrixBase::operator/ diff -r a9d324abb4f1 -r f4d712b82883 Eigen/src/Core/GenericPacketMath.h --- a/Eigen/src/Core/GenericPacketMath.h Tue Jun 08 00:41:33 2010 +0200 +++ b/Eigen/src/Core/GenericPacketMath.h Mon Jun 07 23:30:37 2010 -0600 @@ -221,6 +221,10 @@ /** \internal \returns the log of \a a (coeff-wise) */ template<typename Packet> inline static Packet ei_plog(Packet a) { return ei_log(a); } +/** \internal \returns the log10 of \a a (coeff-wise) */ +template<typename Scalar, typename Packet> inline static Packet ei_plog10(Packet a) +{ return ei_pmul(ei_pset1(Scalar(1)/log(Scalar(10))), ei_plog(a)); } + /** \internal \returns the square-root of \a a (coeff-wise) */ template<typename Packet> inline static Packet ei_psqrt(Packet a) { return ei_sqrt(a); } diff -r a9d324abb4f1 -r f4d712b82883 Eigen/src/Core/MathFunctions.h --- a/Eigen/src/Core/MathFunctions.h Tue Jun 08 00:41:33 2010 +0200 +++ b/Eigen/src/Core/MathFunctions.h Mon Jun 07 23:30:37 2010 -0600 @@ -600,6 +600,44 @@ } /**************************************************************************** +* Implementation of ei_log10 * +****************************************************************************/ + +template<typename Scalar, bool IsInteger> +struct ei_log10_default_impl +{ + static inline Scalar run(const Scalar& x) + { + return std::log10(x); + } +}; + +template<typename Scalar> +struct ei_log10_default_impl<Scalar, true> +{ + static inline Scalar run(const Scalar&) + { + EIGEN_STATIC_ASSERT_NON_INTEGER(Scalar) + return Scalar(0); + } +}; + +template<typename Scalar> +struct ei_log10_impl : ei_log10_default_impl<Scalar, NumTraits<Scalar>::IsInteger> {}; + +template<typename Scalar> +struct ei_log10_retval +{ + typedef Scalar type; +}; + +template<typename Scalar> +inline EIGEN_MATHFUNC_RETVAL(log10, Scalar) ei_log10(const Scalar& x) +{ + return EIGEN_MATHFUNC_IMPL(log10, Scalar)::run(x); +} + +/**************************************************************************** * Implementation of ei_atan2 * ****************************************************************************/ diff -r a9d324abb4f1 -r f4d712b82883 Eigen/src/Core/arch/SSE/MathFunctions.h --- a/Eigen/src/Core/arch/SSE/MathFunctions.h Tue Jun 08 00:41:33 2010 +0200 +++ b/Eigen/src/Core/arch/SSE/MathFunctions.h Mon Jun 07 23:30:37 2010 -0600 @@ -108,6 +108,12 @@ return _mm_or_ps(x, invalid_mask); // negative arg will be NAN } +static EIGEN_DONT_INLINE EIGEN_UNUSED Packet4f ei_plog10(Packet4f x) +{ + _EIGEN_DECLARE_CONST_Packet4f( invlog_0p1, 4.342944819032518721E-1f); + return ei_pmul(ei_p4f_invlog_0p1, ei_plog(x)); +} + static EIGEN_DONT_INLINE EIGEN_UNUSED Packet4f ei_pexp(Packet4f x) { _EIGEN_DECLARE_CONST_Packet4f(1 , 1.0f); diff -r a9d324abb4f1 -r f4d712b82883 Eigen/src/Core/util/ForwardDeclarations.h --- a/Eigen/src/Core/util/ForwardDeclarations.h Tue Jun 08 00:41:33 2010 +0200 +++ b/Eigen/src/Core/util/ForwardDeclarations.h Mon Jun 07 23:30:37 2010 -0600 @@ -116,6 +116,7 @@ template<typename Scalar> struct ei_scalar_sqrt_op; template<typename Scalar> struct ei_scalar_exp_op; template<typename Scalar> struct ei_scalar_log_op; +template<typename Scalar> struct ei_scalar_log10_op; template<typename Scalar> struct ei_scalar_cos_op; template<typename Scalar> struct ei_scalar_sin_op; template<typename Scalar> struct ei_scalar_pow_op; diff -r a9d324abb4f1 -r f4d712b82883 Eigen/src/plugins/ArrayCwiseUnaryOps.h --- a/Eigen/src/plugins/ArrayCwiseUnaryOps.h Tue Jun 08 00:41:33 2010 +0200 +++ b/Eigen/src/plugins/ArrayCwiseUnaryOps.h Mon Jun 07 23:30:37 2010 -0600 @@ -31,7 +31,7 @@ * Example: \include Cwise_exp.cpp * Output: \verbinclude Cwise_exp.out * - * \sa pow(), log(), sin(), cos() + * \sa pow(), log(), log10(), sin(), cos() */ inline const CwiseUnaryOp<ei_scalar_exp_op<Scalar>, Derived> exp() const @@ -52,6 +52,16 @@ return derived(); } +/** \returns an expression of the coefficient-wise base 10 logarithm of *this. + * + * \sa exp() + */ +inline const CwiseUnaryOp<ei_scalar_log10_op<Scalar>, Derived> +log10() const +{ + return derived(); +} + /** \returns an expression of the coefficient-wise square root of *this. * * Example: \include Cwise_sqrt.cpp