Another optimisation possible if N is odd, is to step back by 2 each time.

1000 == 2**3 * 5**3
[download]

125 == 5**3
[download]

5 * 2**2 == 20
[download]

I moved away from using M::B::F. The "best" I came up with is:

#! perl -slw
use strict;

my $NUM = $ARGV[ 0 ] || die 'No arg';
$NUM = eval $NUM if $NUM =~ m[\D];

my $root = sqrt( $NUM );
my $near = int $root;
my $step = $NUM % 2 ? 2 : 1;

$near-= $step while $near and  $NUM % $near;
$near = 1 unless $near;

printf "%10.f (%10f) %.f\n", $NUM, $root, $near;
[download]

It's a simple brute force search for a factor < root.

There may be a way to use the prime factors to speed the search, but the time they take to produce, a linear search down from the root is quicker.

Even with a largish prime like 988041964007, which means a linear search all the way to 1, this takes less than a second, where Math::Big::Factors::Factors_wheel() churns for hours trying to factorise it.

Maybe there's a qucker factorising module out there somewhere? Even so, from what I've tried and seen from other peoples attempts, having the prime factors doesn't give you any obvious way to avoid what amounts to a linear search.

---

Examine what is said, not who speaks.

Silence betokens consent.

Love the truth but pardon error.

The original question started "given N and N's prime factorization", so it seems inappropriate to take into account the time taken to factorize.

In this case though, note that 988041964007 = 7 x 11087 x 12731023 - you are losing accuracy by throwing out the BigInts.

In general if the factorization of N is taking longer than it takes to do sqrt(N) trial divisions then the factorization algorithm is pretty bad, since "trial division up to the square root" is the classical example of a slow algorithm.

Many of the high precision maths libraries out there do factorization - for example the pari library (Math::Pari) has a factorint() function:

  use Math::Pari qw/ factorint /;
  my $n = Math::Pari->new("988041964007");
  my $f = factorint($n);
  my $count = $#{ $f->[0] };
  print join(" * ", map {
    my($prime, $power) = ($f->[0][$_], $f->[1][$_]);
    $power > 1 ? "$prime ^ $power" : $prime
  } 0 .. $count), "\n";
7 * 11087 * 12731023
[download]

Note that this algorithm does not guarantee primality: for large factors you need to verify with an isprime() check.

Hugo