Re^5: OT: Finding Factor Closest To Square Root

[BrowserUk]

Yeah! I'm reaching a similar conclusion..now if only I could work out to use Algorithm::Loops, it'd be done in a trice :)

[Roy Johnson]

Here's a recursive solution to find the closest set of factors to a target (with wrapper script). I think it's solid.

use strict; 
use warnings;
use List::Util qw[ reduce ];

my @pfs = ( 2,3,3,3,7, 997 );
my $num = reduce { $a * $b } @pfs;
my $root = sqrt $num;
print "N=$num, R=$root\n";

sub closest {
  # Args are target and factor-list
  my ($target, $car, @cdr) = @_;
  @cdr or return $car;
  reduce { abs($target - $a) < abs($target - $b) ? $a : $b }
         ($car*closest($target/$car, @cdr), closest($target, @cdr), $c
+ar);
}

print closest($root, @pfs), "\n";

[BrowserUk]

It seems to work well as far as I can verify it, with one exception. Like one of my other attempts, if N is prime, it produces N as the nearest factor rather than 1.

However, being recursive, it suffers in performance one you really start making it work.

Try 8585937356. I get a nearest factor of 4 with my dumb crunch code. I gave up waiting for your's :)

[Roy Johnson]

When I plug in 1 and 5 as factors, I get 1. Does your factor list include 1?

[Limbic~Region]

BrowserUk,
The following code demonstrates that this is not actually true:
Read more... (2 kB)
• If a subset has been seen before, it is skipped prior to determining the product of the subset
• If a subset is the start of a progression of subsets previously seen, the entire progression is skipped
• If a factor is determined to be larger than the sqrt, determining if it is the closest is skipped
Unfortunately, I do not know if this is the fastest implementation. I have several variations ranging from bitstrings to avoiding push/slice. The trouble is that I have found an apparent bug in Math::Pari so I am not sure which solution should be fastest. The bug only seems to manifest itself when I add the sqrt() code. To see the bug for yourself using 24777695232, add print "$_\n" for @factor; right after the array is created. To see it go away, comment out the $sqrt definition and change the return statement to return $f.

[tall_man]

I got an entirely different error when I did a sqrt($N) before factoring:

  DB<1> n
main::(factsqrt.pl:5):    my $N = $_;
  DB<1> n
main::(factsqrt.pl:6):    my $sqrt2 = sqrt($N);
  DB<1> n
main::(factsqrt.pl:7):    my $a = factorint($N);
  DB<1> n
  DB<1> p $a
Mat([-1,1])

However, this is not a bug in Math::Pari, but a side effect of perl's conversion of string to numbers.

The number "24777695232" first comes in as a string from the command line. Applying the function "sqrt" to it converts it to a double. When you pass it to Math::Pari, it goes in as a double. Naturally, the module fails to do the right thing. If you call a Math::Pari function first, like factorint, it automatically converts the string into a long Pari integer and operates on that.

[Limbic~Region]

tall_man,
Interesting. That would imply I couldn't do this either:

my $a = factorint( 24777695232 );

It seems like a bug to me, but I am finding the documentation on Math::Pari lacking. Is there a good tutorial that covers all of the functions offered and examples on how to use them?

[BrowserUk]

The bug only seems to manifest itself when I add the sqrt() code. .... To see it go away, comment out the $sqrt definition and change the return statement to return $f.

Sorry Limbic~Region, I do not follow these instructions?

The only definition of $sqrt I can see is     my $sqrt = sqrt( $N );

But I do not understand what you mean by "change the return statement to return $f."?

[Limbic~Region]

BrowserUk,
Sorry - it was past my bedtime last night when I posted. The code in question (with a few extra lines for context looks like this:

    my $N = shift;
    my $sqrt = sqrt( $N );
    my @factor = prime_factors( $N );
    return 1 if @factor == 1;
        # ... rest of the code
        return $f < $sqrt ? $f : undef;

To see the problem, add the following line:

    my $N = shift;
    my $sqrt = sqrt( $N );
    my @factor = prime_factors( $N );
    # ADD THIS LINE
    print "$_\n" for @factor;
    return 1 if @factor == 1;
        # ... rest of the code
        return $f < $sqrt ? $f : undef;

And to see the problem go away:

    my $N = shift;    
    my @factor = prime_factors( $N );
    # KEEP THIS LINE
    print "$_\n" for @factor;
    return 1 if @factor == 1;
        # ... rest of the code
        return $f;

You can see the list of factors changes with the presence of the sqrt() code (at least it does on my machine).

[Roy Johnson]

Here's an Algorithm::Loops solution, then, based on lidden's nifty solution:

use strict; 
use warnings;
use List::Util qw[ reduce ];
use Algorithm::Loops 'NextPermuteNum';

sub closest{
    my $exact = sqrt(shift);
    my $best = 1;
    do {
        my $guess = 1;
        for (@_) {
        $guess *= $_;
            $best = $guess if (abs($exact - $guess) < abs($exact - $be
+st));
            last if $guess > $exact;
        }
    } while NextPermuteNum(@_);
    $best;
}

my @pfs = (2,2,3,3,5,5,7);
# Compute product of factors
our ($a, $b);
my $num = reduce { $a * $b } @pfs;
my $root = sqrt $num;
print "N=$num, R=$root\n";

print closest($num, @pfs), "\n";