Speeding up FactorInteger for product of two primes

Question

I have a large integer $N$ of size about 10^150. I know that $N$ is a product of two primes $p$ and $q$. I also know that both $p$ and $q$ are of roughly equal size, so one of them is not, for example, 2 or 3.

Can I use this information to somehow speed up FactorInteger ? Or is factoring such a large number already beyond the scope of a average desktop computer. I also noticed that FactorInteger uses one of my cpu cores only, is it possible to parallelize this?

Answer 1

It all depends on "p and q are of roughly equal size". If they differ no more than say 10^6, the following works rather fast. E.g.:

fac[big_] := Module[{ne = NextPrime[Sqrt[big] - 1]},
  While[Mod[big, ne] > 0, ne = NextPrime[ne]];
  {ne,big/ne}
  ]
n = 80;
fac[NextPrime[10^n] NextPrime[10^n + 10^6] ]
(* {100000000000000000000000000000000000000000000000000000000000000000000

000001000507, 

1000000000000000000000000000000000000000000000000000000000000000000000

00000000129} )*

Answer 2

If x is around 10^150, let y = floor(sqrt(x)), and if y^2≠x then solve (y+k+1)(y-k)=x, (y+k+2)(y-k) = x, (y + k + 3)(y - k) = x etc. until you find a k that is an integer. This will cover a huge range, much larger than 10^6.

We know 0 <= d = x - y^2 <= 2y ≈ 2 x 10^75.

Take for example (y+k+2)(y-k) = y^2 + yk + 2y - ky - k^2 - 2k. We solve 2y - k^2 - 2k = d or k^2 + 2k + 1 = 2y +1 - d. k will be between sqrt(y) and sqrt(2y), that is between 4.5 x 10^37 and 6.3 x 10^37. You can halve the number of calculations by not calculating k if the solution would have to consist of an even and an odd number.

(y+k+c)(y-k) = x will have a solution around sqrt (cy) = sqrt(2c) * x^(1/4). j calculations let you cover even c up to 2c, with factors up to 4 sqrt(j) x^(1/4) apart.

Answer 3

Translating some of @gnasher729's suggestions to Mathematica gives indeed a superfast recipe:

$MaxExtraPrecision = 10^3;
f[x_Integer] := Module[{y, i, z, s},
  y = Floor[Sqrt[x]];
  Catch[
    i = 0;
    While[True,
      z = (2 y + i)^2 - 4 x;
      If[z >= 0 && IntegerQ[Sqrt[z]],
        s = i/2 + y + {1, -1} Sqrt[z]/2;
        If[IntegerQ[s[[1]]] && IntegerQ[s[[2]]], Throw[s]]];
      i++]]]
f[15]
(*    {5, 3}    *)
f[16]
(*    {4, 4}    *)
f[17]
(*    {17, 1}    *)
n = 80;
f[NextPrime[10^n] NextPrime[10^n + 10^6]] // AbsoluteTiming
(*    {0.010925,
       {100000000000000000000000000000000000000000000000000000000000000000000000001000507,
        100000000000000000000000000000000000000000000000000000000000000000000000000000129}}    *)
f[NextPrime[10^n] NextPrime[10^n + 10^12]] // AbsoluteTiming
(*    {0.007086,
       {100000000000000000000000000000000000000000000000000000000000000000001000000000191,
        100000000000000000000000000000000000000000000000000000000000000000000000000000129}}    *)

Things get a lot slower once the difference between the two numbers becomes larger than their square root:

f[NextPrime[10^n] NextPrime[10^n + 10^43]] // AbsoluteTiming
(*    {132.132,
       {100000000000000000000000000000000000010000000000000000000000000000000000000000357,
        100000000000000000000000000000000000000000000000000000000000000000000000000000129}}    *)

The remaining point is a fast method to determine IntegerQ[Sqrt[z]] for very large values of z. Using @MichaelE2's fast perfect-square test gives a 150-fold speedup, which is however no match for the exponentially increasing difficulty with factor-spacing:

$MaxExtraPrecision = 10^3;
f[x_Integer] := Module[{y, i, z, s},
  y = Floor[Sqrt[x]];
  Catch[
    i = 0;
    While[True,
    z = (2 y + i)^2 - 4 x;
    If[z >= 0 && FractionalPart[Sqrt[z + 0``1]] == 0,
      s = i/2 + y + {1, -1} Sqrt[z]/2;
      If[IntegerQ[s[[1]]], Throw[s]]];
    i++]]]
n = 80;
f[NextPrime[10^n] NextPrime[10^n + 10^43]] // AbsoluteTiming
(*    {0.894308,
       {100000000000000000000000000000000000010000000000000000000000000000000000000000357,
        100000000000000000000000000000000000000000000000000000000000000000000000000000129}}    *)
f[NextPrime[10^n] NextPrime[10^n + 10^44]] // AbsoluteTiming
(*    {92.1891,
       {100000000000000000000000000000000000100000000000000000000000000000000000000000179,
        100000000000000000000000000000000000000000000000000000000000000000000000000000129}}