I tried something different. I gave the problem to GPT and asked it to optimise the code. The first time around it gave code that didn't work. The next time the code did appear to do what it is supposed to do. On my computer at least the optimised version is vastly faster than the posted code. Here is the discussion with GPT.
So the original code is:
p = 2564855351
Monitor[For[i = 3, i <= p, i += 2,
If[(Sqrt[p^3/(p*i^2 + i)] - p) < 0.5, x = i, Print[i], i = p + 1]];
While[x <= p, x + 2 If[Divisible[p, x], Print[x], p = 1];], I]
The optimised code, which interrupts once it finds the first divisor, is:
p = 2564855351;
x = 3;
Monitor[While[x <= p, If[(Sqrt[p^3/(p*x^2 + x)] - p) < 0.5, Print[x];
Break[];];
x += 2;];
If[x <= p, While[x <= p, If[Divisible[p, x], Print[x];
Break[];];
x += 2;];], x]
Which I did not really test, but it gives
41227
If you run the built in Mathematica function
Divisors[p]
you get
{1, 41227, 62213, 2564855351}
This suggests that the code performs ok in this case. As you know the number is relatively prime, you can easily divide p/41227 to obtain 62213. I did not analyse the code or look at any other cases in more detail. It would also be possible to parallelise the entire code, which could give you another speed boost.
Best wishes,
Marco
