9

1

**Goal**

The goal is to efficiently check whether a number is a perfect power.

**Attempts**

It is possible to check whether a number is a perfect power using `FactorInteger`

, as seen in (42251), but that method is very slow, especially for large numbers.

Since the upper bound of my list of interest is around `10^18`

, I realized that if a number in my list is a perfect power, its power would be at most `59`

(= `Floor[Log2[10^18]]`

).

Thus I wrote the following code:

```
powerq[in_] := Evaluate[Or @@ Array[FractionalPart[in^(1/Prime[#])] == 0 &, PrimePi[59]]];
read = OpenRead["numbers"];
write = OpenWrite["powers"];
Do[If[powerq[#], Write[write, #], , Break[]]&@Read[read], Infinity]
```

(I am using `Read`

and `Write`

because my list has approximately 2 billion numbers (~40 GB) and cannot be stored in RAM; I used `FractionalPart`

instead of `IntegerQ`

because `IntegerQ`

does not work symbolically)

This was slow, so I decided to use `Compile`

```
powerqcomp = Compile[{{in, _Integer}},
Evaluate[
Or @@ Array[FractionalPart[in^(1/Prime[#])] == 0 &,
PrimePi[59]]], CompilationTarget -> "C"]
```

The problem was that the `Compile`

d function uses machine real, so it sometimes returns fallacious answer, such as:

```
powerqcomp[2^15]
(* False *)
```

So, I used a different approach.

If a number $n$ can be expressed as $a^b$ ($a$ and $b$ are positive integers and $b > 1$), then $b \leq \log_{2}{n}$. Then, I could check whether there exists $a$ such that $n = a^b$ for at least one $b$.

Since my list consists of large numbers, trying all $a$ from $1$ to $n^{(1/b)}$ is out of the question. Perhaps I could use `BinarySearch`

, but then I would need to generate `Range[n^(1/b)]`

; I would need to make my own function that performs a binary search.

This is the code I made that uses binary search. It does work, but it is still too slow (and gives `CompiledFunction::cfn`

errors sometimes).

```
powerqcomp2 =
Compile[{{in, _Integer}},
Module[{up = 10^9, down = 1,
blist = Select[{2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41,
43, 47, 53, 59}, # <= Log[2, in] &]},
Catch[Do[
While[True,
Which[#^b == in, Throw[True],
#^b > in, up = #,
True, down = #
] &[Ceiling[(up + down)/2]];
If[up - down == 1, Break[]]
];
up = 10^9;
down = 1,
{b, blist}]; False]]]
```

**Question**

How could one create code that efficiently checks whether a number is a perfect power?

Have you seen this?

– J. M.'s ennui – 2016-07-23T20:43:34.640@J.M. I read the article, but it seems like

Mathematicais not designed to perform operations in there. E.g. separating the mantissa and the exponent from a floating-point number.`MantissaExponent`

calculates the mantissa and exponent, which is slower than simple extraction. Plus, bitwise operations are not permitted on floating-point reals in Mma... How could that algorithm be implemented in Mma? – JungHwan Min – 2016-07-24T03:18:46.243As there are only a million cubes in the range of interest (and less of higher powers) perhaps you should build a list of all these and only check for exact squares. – mikado – 2016-07-24T10:02:55.823