I'm looking for the highest-performance method of calculating integer square roots in Mathematica of very big arbitrary-precision numbers.
As an example testcase, I use:
n = 10^1000000 - 3 ^ 2095903
On my Core i7 linux machine, calculating the integer (floor) square root using the straightforward method takes 2.68 seconds:
In:= n = 10^1000000 - 3^2095903; In:= First@Timing@Floor[Sqrt[n]] Out= 2.68984
However, the same machine can calculate the integer square root much faster using the GMP library. Here's an example python program doing just that:
from gmpy2 import mpz, isqrt from time import time n = mpz(10**1000000 - 3**2095903) t0 = time() s = isqrt(n) t1 = time() print(t1 - t0)
On the same machine, this is approximately 160 times faster:
$ python3 benchmark.py 0.01685619354248047
The frustrating thing is that Mathematica already depends on the GMP library for its arbitrary-precision integers, but I cannot seem to make Mathematica use the GMP implementation of the integer square root.
Instead I suspect that its Floor[Sqrt[...]] function falls back on a generic algorithm that must calculate its argument to sufficient precision using arbitrary-precision floating-point approximation.
Is there any way to speed up the calculation of integer square roots, preferably approaching GMP's raw performance?