1#!/usr/bin/env perl 2use strict; 3use warnings; 4use Math::Prime::Util qw/factor/; 5# Compare to Math::Factor::XS, which uses trial division. 6use Math::Factor::XS qw/prime_factors/; 7 8use Benchmark qw/:all/; 9use List::Util qw/min max reduce/; 10my $count = shift || -2; 11my $is64bit = (~0 > 4294967295); 12my $maxdigits = ($is64bit) ? 20 : 10; # Noting the range is limited for max. 13my $semiprimes = 0; 14my $howmany = 1000; 15 16for my $d ( 3 .. $maxdigits ) { 17 print "Factor $howmany $d-digit numbers\n"; 18 test_at_digits($d, $howmany); 19} 20 21sub test_at_digits { 22 my $digits = shift; 23 die "Digits has to be >= 1" unless $digits >= 1; 24 die "Digits has to be <= $maxdigits" if $digits > $maxdigits; 25 my $quantity = shift; 26 27 my @rnd = ndigit_rand($digits, $quantity); 28 my @smp = genrough($digits, $quantity); 29 30 # verify (can be _really_ slow for 18+ digits) 31 foreach my $p (@rnd, @smp) { 32 next if $p < 2; 33 verify_factor($p, [prime_factors($p)], [factor($p)], "Math::Prime::Util $Math::Prime::Util::VERSION"); 34 } 35 36 #my $min_num = min @nums; 37 #my $max_num = max @nums; 38 #my $whatstr = "$digits-digit ", $semiprimes ? "semiprime" : "random"; 39 #print "factoring 1000 $digits-digit ", 40 # $semiprimes ? "semiprimes" : "random numbers", 41 # " ($min_num - $max_num)\n"; 42 43 my $lref = { 44 "MPU random" => sub { my@a=factor($_) for @rnd }, 45 "MPU nonsmooth" => sub { my@a=factor($_) for @smp }, 46 "MFXS random" => sub { my@a=prime_factors($_) for @rnd }, 47 "MFXS nonsmooth" => sub { my@a=prime_factors($_) for @smp }, 48 }; 49 cmpthese($count, $lref); 50} 51 52sub verify_factor { 53 my ($n, $aref1, $aref2, $name) = @_; 54 55 return 1 if "@$aref1" eq "@$aref2"; 56 57 my @master = @$aref1; 58 my @check = @$aref2; 59 die "Factor $n master fail!" unless $n == reduce { $a * $b } @master; 60 die "Factor $n fail: $name" unless $#check == $#master; 61 die "Factor $n fail: $name" unless $n == reduce { $a * $b } @check; 62 for (0 .. $#master) { 63 die "Factor $n fail: $name" unless $master[$_] == $check[$_]; 64 } 65 1; 66} 67 68sub genrough { 69 my ($digits, $num) = @_; 70 71 my @min_factors_by_digit = (2,2,3,5,7,13,23,47,97); 72 my $smallest_factor = $min_factors_by_digit[$digits]; 73 $smallest_factor = $min_factors_by_digit[-1] unless defined $smallest_factor; 74 75 my @semiprimes; 76 foreach my $i (1 .. $num) { 77 my $n; 78 my @facn; 79 do { 80 $n = ndigit_rand($digits, 1); 81 @facn = Math::Prime::Util::trial_factor($n,$smallest_factor); 82 } while scalar(@facn) > 1; 83 push @semiprimes, $n; 84 } 85 return @semiprimes; 86} 87 88use Bytes::Random::Secure qw/random_string_from/; 89sub ndigit_rand { 90 my($digits, $howmany) = @_; 91 die "digits must be > 0" if $digits < 1; 92 $howmany = 1 unless defined $howmany; 93 # TODO: need to skip things larger than ~0 for this module 94 my @nums = map { random_string_from("123456789",1) . random_string_from("0123456789",$digits-1) } 1 .. $howmany; 95 if (10**$digits > ~0) { @nums = map { Math::BigInt->new($_) } @nums; } 96 else { @nums = map { int($_) } @nums; } 97 return wantarray ? @nums : $nums[0]; 98} 99