pow for Integral Values
Asked Answered
M

2

6

I need a version of pow for integers. I have 2 problems that I need to solve with pow:

  1. If the result is larger than my integral type I need to clamp to numeric_limits::max()
  2. I need to be able to deal with 41.99999 rounding up to 42, rather than down to 41

Does C++ provide me some kind of inline solution here, or am I stuck writing my own function:

template <typename T>
enable_if_t<is_integral_v<T>, T> mypow(const T base, unsigned int exp) {
    T result = exp == 0U ? base : 1;

    while(exp-- > 1U) {
        if(numeric_limits<T>::max() / result <= base) return numeric_limits<T>::max();
        result *= base;
    }
    return result;
}
Moron answered 15/6, 2017 at 15:8 Comment(3)
Can't you use std::pow and truncate the result if needed?Chancemedley
@Chancemedley Yeah... if there was a way to guarantee 1 and 2 wouldn't be violated by that. Is there a way to accomplish that?Moron
You can use tail recursion and some numeric trickery to improve runtime performance. Further, you can make this method constexpr. #9349433Petition
R
3

Does C++ provide me some kind of inline solution here

No, there is no integer pow in the standard library.

or am I stuck writing my own function

Yes, you get to write your own function. Note that your shown multiplication loop may be slower than using std::pow to implement the function, especially since you also have a branch and division in the loop:

template<class I>
I int_pow_no_overflow(I base, I exp)
{
    double max = std::numeric_limits<I>::max();
    double result = std::round(std::pow(base, exp));
    return result >= max
        ? max
        : result;
}

For a more generic approach, you may want to consider underflow as well.

There are as well other, faster algorithms (see for example Exponentiation by squaring) for integer exponentiation than the linear one you showed, but I'm unsure whether it would be worth considering them unless you deal with arbitrary precision arithmetic, or an embedded system with no floating point unit.

Radiometer answered 15/6, 2017 at 15:28 Comment(15)
If you have a better way to avoid the truncation and to clamp using pow I'm all ears. What I've got was the best I could come up with, and it is very expensive.Moron
Why do you need the generic programming?Vermis
@JonathanMee as a clue, N^8 = N^4 * N^4 so you don't need to calculate N^4 twice.Reasoning
@JonathanMee how about the one I added to the answer?Radiometer
@KlitosKyriacou Hmmm... I see what you're saying, I could have a logrithmic reduction in my loop count if I could just figure out how to do this...Moron
@user2079303 Probably is as good as we can do :/ We could inline too I suppose: static_cast<int>(min<double>(numeric_limits<int>::max(), round(pow(base, exp))))Moron
Problem with wrapping std::pow is that it works with double, and (at least on some platforms) a long long has some values that cannot be precisely represented by a double.Reasoning
@KlitosKyriacou good point. If long long is needed, I would write a specialization with explicit conversions to long double to use that overload of std::pow.Radiometer
@klitos-kyriacou That's a very good point. This could make a very interesting example as a recursive function.Huan
@user2079303 long double isn't necessarily able to store long long, because on platforms other than x86 it's still mostly the same as doubleQuintain
It's actually certain that it can't. not enough mantissa bits available.Huan
@LưuVĩnhPhúc well, you could check numeric_limits whether there are enough significand digits to represent the target values. If you really need to support long longs on all platforms, then perhaps fall-back to purely integer based algorithm.Radiometer
@MichaëlRoy not at all certain. 80-bit long doubles are very typical, and have just the appropriate number of mantissa bits to represent signed 64 bit integers.Radiometer
MS still doesn't support anything longer than the 64 bit double.Huan
@MichaëlRoy or rather: no longer supports :) They used to.Radiometer
R
2

Your code did not compile, you should if possible check your code compiles first, use your compiler, or check it on compiler explorer first.

Also, you forgot to take into account negative values. That's a very important feature of integral powers. Code below is for regular int type. I'll let you explore how you could extend it for other integral types.

#include <type_traits>
#include <iostream>
#include <cmath>
#include <limits>

using namespace std;

template <typename T>
enable_if_t< is_integral<T>::value, T> 
mypow(T base, unsigned int exp) 
{
    T result = T(1);
    bool sign = (base < 0);

    if (sign) base = -base;

    T temp = result;
    while(exp-- != 0) 
    {
        temp *= base;
        if (temp < result)
        {
            return (sign) ? numeric_limits<T>::min() 
                          : numeric_limits<T>::max();
        }
        result = temp;
    }
    return (sign && (exp & 1)) ? -result : result;
}

template <typename T>
enable_if_t< !is_integral<T>::value, int> 
mypow(const T& base, unsigned int exp) 
{
    T result = T(1);
    int i_base = int(floor(base + .5));
    bool sign = (i_base < 0);

    if (sign) i_base = -i_base;

    int temp = result;
    while(exp-- != 0) 
    {
        temp *= i_base;
        if (temp < result)
        {
            return (sign) ? numeric_limits<int>::min() : numeric_limits<int>::max();
        }
        result = temp;
    }
    return (sign && (exp & 1)) ? -result : result;
}

In real life, I would do this note the use of floor, even in the integral case.

  template<typename T>
   enable_if_t< is_integral<T>::value, T> 
   mypow(T x, unsigned int y) { return T(floor(pow(x, y) + .5)); }

   template<typename T>
   enable_if_t< !is_integral<T>::value, int> 
   mypow(T x, unsigned int y) { return int(floor(pow(floor(x + .5), y) + .5)); }
Rolon answered 15/6, 2017 at 15:50 Comment(5)
I only intended this as an example, but yeah, this is one more reason that I don't want to go with writing my own function :(Moron
Using floating point pow will generally be faster on CPU that support floating point ops. Some modern CPUs have a floating-point pow instruction.Huan
@JonathanMee: Maybe you should check this interesting discussion on the subject: #101939Huan
If x is_integral, then why are you doing floor(x+.5) rather than double(x)? are you worried about 64bit integers not representable as double? They screw your code anyway.Enlistment
@Enlistment sorry, it's a typo in second function, which is for floating points. I'm rounding the value, since that was part of the original question. Thanks for pointing this out..Huan

© 2022 - 2024 — McMap. All rights reserved.