Algorithm to find nth root of a number

Viewed 37134

I am looking for an efficient algorithm to find nth root of a number. The answer must be an integer. I have found that newtons method and bisection method are popular methods. Are there any efficient and simple methods for integer output?

5 Answers

Here is an efficient general implementation in C, using a simplified version of the "shifting nth root algorithm" to compute the floor of the nth root of x:

uint64_t iroot(const uint64_t x, const unsigned n)
{
  if ((x == 0) || (n == 0)) return 0;
  if (n == 1) return x;

  uint64_t r = 1;
  for (int s = ((ilog2(x) / n) * n) - n; s >= 0; s -= n)
  {
    r <<= 1;
    r |= (ipow(r|1, n) <= (x >> s));
  }

  return r;
}

It needs this function to compute the nth power of x (using the method of exponentiation by squaring):

uint64_t ipow(uint64_t x, unsigned n)
{
  if (x <= 1) return x;

  uint64_t y = 1;
  for (; n != 0; n >>= 1, x *= x)
    if (n & 1)
      y *= x;
  return y;
}

and this function to compute the floor of base-2 logarithm of x:

int ilog2(uint64_t x)
{
  #if __has_builtin(__builtin_clzll)
    return 63 - ((x != 0) * (int)__builtin_clzll(x)) - ((x == 0) * 64);
  #else
    int y = -(x == 0);
    for (unsigned k = 64 / 2; k != 0; k /= 2)
      if ((x >> k) != 0)
        { x >>= k; y += k; }
    return y;
  #endif
}

Note: This assumes that your compiler understands GCC's __has_builtin test and that your compiler's uint64_t type is the same size as an unsigned long long.

You can try this C function to get the nth_root of an unsigned integer :

unsigned initial_guess_nth_root(unsigned n, unsigned nth){
    unsigned res = 1;
    for(; n >>= 1; ++res);
    return nth ? 1 << (res + nth - 1) / nth : 0 ;
}

// return a number that, when multiplied by itself nth times, makes N.
unsigned nth_root(const unsigned n, const unsigned nth) {
    unsigned a = initial_guess_nth_root(n , nth), b, c, r = nth ? a + (n > 0) : n == 1 ;
    for (; a < r; b = a + (nth - 1) * r, a = b / nth)
        for (r = a, a = n, c = nth - 1; c && (a /= r); --c);
    return r;
}

Example of output :

  24 == (int) pow(15625, 1.0/3)

  25 == nth_root(15625, 3)
   0 == nth_root(0, 0)
   1 == nth_root(1, 0)
   4 == nth_root(4096, 6)
  13 == nth_root(18446744073709551614, 17) // 64-bit 20 digits
  11 == nth_root(340282366920938463463374607431768211454, 37) // 128-bit 39 digits

Here is the github source.

Related