Take the 2-minute tour ×
MathOverflow is a question and answer site for professional mathematicians. It's 100% free, no registration required.

We have a natural number $n>1$. We want to determine whether exists any natural numbers $a, k>1$ such that $n = a^k$.

Please suggest a polynomial-time algorithm.

share|improve this question
3  
Naive answer: approximate the logarithm to about the same number of places as digits in n (there exist algorithms polynomial time in the number of digits), then check if kth roots are integers for k < log n. Each step takes polynomial time, so the the algorithm terminates in polynomial time. –  S. Carnahan Feb 2 '10 at 18:12
2  
    
Simply compute $n^{1/k}$ for $k=2,\ldots,\lfloor \log\_2 n \rfloor$. Arbitrary roots can be computed in polynomial time (directly or by using logarithms, as Scott Carnahan said), so this is a polynomial-time algorithm. –  Darsh Ranjan Feb 2 '10 at 18:26
    
I am ignorant in CS .. so this may be a dumb question. Expressing a number to the base p is not a polynomial time problem? –  Anweshi Feb 2 '10 at 18:56
1  
Oops, sorry, then one would need to determine the prime factors of a number too.. Which complicates the issue. Anyway I am totally ignorant of this type of things, as I said. –  Anweshi Feb 2 '10 at 19:34
show 2 more comments

4 Answers 4

up vote 15 down vote accepted

This can be done in "essentially linear time." Check out Daniel Bernstein's website: http://cr.yp.to/arith.html

Especially note his papers labeled [powers] and [powers2].

share|improve this answer
add comment

In order to test whether or not a natural number $n$ is a perfect power, we can conduct a binary search of the integers {1,2,...,n} for a number a number $m$ such that $n = m^b$ for some $b>1$. Let $b>1$. If a solution $m$ to $m^b =n$ exists,then it must lie in some interval $[c_i,d_i]$. When $i = 0$ we may take $[c_0,d_0] = [1,n]$. To define $[c_{i+1},d_{i+1}]$, consider $\alpha:= \left\lfloor \frac{(ci+di)}{2}\right\rfloor$. If $\alpha^b = n$ then we’re done. If $\alpha^b > n$, let $[c_{i+1}, d_{i+1}] = [c_i, \alpha]$; otherwise $\alpha^b < n$ and we let $[c_{i+1}, d_{i+1}] = [\alpha, d_i]$. We continue in this manner until $|c_i − d_i| \leq 1$. We then increase the value stored in variable $b$ and start the loop again. Performing this loop for all $b \leq log(n)$ completes the algorithm.

A pseudocode implementation of this algorithm can be found on page 21 of Dietzelbinger's Primality Testing in Polynoial Time. Its complexity is about $O(log^3(n))$.

share|improve this answer
add comment

For each $k \le \log n/\log 2$, compute an approximation to the positive real $k$-th root of $n$ using Newton's method to enough precision to check if it is an integer. Alternatively, use $p$-adic roots for a suitable $p$, with Newton turning into Hensel.

share|improve this answer
add comment

The computer algebra system GAP performs this test and determines a smallest root $a$ of a given integer $n$ quite efficiently. The following is copied directly from its source code (file gap4r6/lib/integer.gi), and should be self-explaining:

#############################################################################
##
#F  SmallestRootInt( <n> )  . . . . . . . . . . . smallest root of an integer
##
InstallGlobalFunction(SmallestRootInt,

  function ( n )

    local   k, r, s, p, l, q;

    # check the argument
    if   n > 0  then k := 2;  s :=  1;
    elif n < 0  then k := 3;  s := -1;  n := -n;
    else return 0;
    fi;

    # exclude small divisors, and thereby large exponents
    if n mod 2 = 0  then
        p := 2;
    else
        p := 3;  while p < 100  and n mod p <> 0  do p := p+2;  od;
    fi;
    l := LogInt( n, p );

    # loop over the possible prime divisors of exponents
    # use Euler's criterion to cast out impossible ones
    while k <= l  do
        q := 2*k+1;  while not IsPrimeInt(q)  do q := q+2*k;  od;
        if PowerModInt( n, (q-1)/k, q ) <= 1  then
            r := RootInt( n, k );
            if r ^ k = n  then
                n := r;
                l := QuoInt( l, k );
            else
                k := NextPrimeInt( k );
            fi;
        else
            k := NextPrimeInt( k );
        fi;
    od;

    return s * n;
end);
share|improve this answer
add comment

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.