[libmath-prime-util-perl] 31/59: Update for older Perls and bigint
Partha P. Mukherjee
ppm-guest at moszumanska.debian.org
Thu May 21 18:44:57 UTC 2015
This is an automated email from the git hooks/post-receive script.
ppm-guest pushed a commit to annotated tag v0.10
in repository libmath-prime-util-perl.
commit c1a3e12cf3ebb8c4842818bf42ff01fdb6058032
Author: Dana Jacobsen <dana at acm.org>
Date: Fri Jul 6 04:05:11 2012 -0600
Update for older Perls and bigint
---
examples/bench-is-prime.pl | 2 +-
examples/bench-miller-rabin.pl | 2 +-
examples/bench-random-prime.pl | 2 +-
examples/sophie_germain.pl | 2 +-
examples/twin_primes.pl | 2 +-
lib/Math/Prime/Util.pm | 86 +++++++++++++++++++++++++++---------------
lib/Math/Prime/Util/PP.pm | 19 ++++------
t/12-nextprime.t | 6 ++-
t/14-nthprime.t | 21 ++++++-----
t/80-pp.t | 7 +++-
t/81-bignum.t | 20 ++++------
11 files changed, 98 insertions(+), 71 deletions(-)
diff --git a/examples/bench-is-prime.pl b/examples/bench-is-prime.pl
index 83fc71d..7930a18 100755
--- a/examples/bench-is-prime.pl
+++ b/examples/bench-is-prime.pl
@@ -3,7 +3,7 @@ use strict;
use warnings;
#use Math::Primality;
use Math::Prime::XS;
-use Math::Prime::Util;
+use Math::Prime::Util '-nobigint';
#use Math::Pari;
#use Math::Prime::FastSieve;
use Benchmark qw/:all/;
diff --git a/examples/bench-miller-rabin.pl b/examples/bench-miller-rabin.pl
index b79e121..02c4478 100755
--- a/examples/bench-miller-rabin.pl
+++ b/examples/bench-miller-rabin.pl
@@ -3,7 +3,7 @@ use strict;
use warnings;
#use Math::Primality;
use Math::Prime::XS;
-use Math::Prime::Util;
+use Math::Prime::Util '-nobigint';
#use Math::Prime::FastSieve;
use Benchmark qw/:all/;
use List::Util qw/min max/;
diff --git a/examples/bench-random-prime.pl b/examples/bench-random-prime.pl
index e60a221..13ea69e 100755
--- a/examples/bench-random-prime.pl
+++ b/examples/bench-random-prime.pl
@@ -2,7 +2,7 @@
use strict;
use warnings;
-use Math::Prime::Util qw/random_prime random_ndigit_prime/;
+use Math::Prime::Util qw/-nobigint random_prime random_ndigit_prime/;
use Benchmark qw/:all/;
use List::Util qw/min max/;
my $count = shift || -3;
diff --git a/examples/sophie_germain.pl b/examples/sophie_germain.pl
index ed7dde3..38dd3bd 100755
--- a/examples/sophie_germain.pl
+++ b/examples/sophie_germain.pl
@@ -2,7 +2,7 @@
use strict;
use warnings;
-use Math::Prime::Util qw/next_prime is_prime/;
+use Math::Prime::Util qw/-nobigint next_prime is_prime/;
my $count = shift || 20;
diff --git a/examples/twin_primes.pl b/examples/twin_primes.pl
index 2b3fbab..aa1ec5a 100755
--- a/examples/twin_primes.pl
+++ b/examples/twin_primes.pl
@@ -2,7 +2,7 @@
use strict;
use warnings;
-use Math::Prime::Util qw/next_prime is_prime/;
+use Math::Prime::Util qw/-nobigint next_prime is_prime/;
my $count = shift || 20;
diff --git a/lib/Math/Prime/Util.pm b/lib/Math/Prime/Util.pm
index 844ba75..a7a814e 100644
--- a/lib/Math/Prime/Util.pm
+++ b/lib/Math/Prime/Util.pm
@@ -16,6 +16,7 @@ our @EXPORT_OK = qw(
prime_precalc prime_memfree
is_prime is_prob_prime
is_strong_pseudoprime is_strong_lucas_pseudoprime
+ miller_rabin
primes
next_prime prev_prime
prime_count prime_count_lower prime_count_upper prime_count_approx
@@ -35,13 +36,15 @@ sub import {
}
sub _import_nobigint {
- undef *factor; *factor = \&_XS_factor;
- undef *is_prime; *is_prime = \&_XS_is_prime;
- undef *next_prime; *next_prime = \&_XS_next_prime;
- undef *prev_prime; *prev_prime = \&_XS_prev_prime;
- undef *prime_count; *prime_count = \&_XS_prime_count;
- undef *nth_prime; *nth_prime = \&_XS_nth_prime;
+ undef *factor; *factor = \&_XS_factor;
+ undef *is_prime; *is_prime = \&_XS_is_prime;
+ undef *is_prob_prime; *is_prob_prime = \&_XS_is_prob_prime;
+ undef *next_prime; *next_prime = \&_XS_next_prime;
+ undef *prev_prime; *prev_prime = \&_XS_prev_prime;
+ undef *prime_count; *prime_count = \&_XS_prime_count;
+ undef *nth_prime; *nth_prime = \&_XS_nth_prime;
undef *is_strong_pseudoprime; *is_strong_pseudoprime = \&_XS_miller_rabin;
+ undef *miller_rabin; *miller_rabin = \&_XS_miller_rabin;
}
my %_Config;
@@ -116,6 +119,12 @@ my $_XS_MAXVAL = $_Config{'xs'} ? $_Config{'maxparam'} : -1;
#
# $n = $n->numify if $n < ~0 && ref($n) =~ /^Math::Big/;
# get us out of big math if we can
+#
+# Sadly, non-modern versions of bignum (5.12.4 and earlier) completely make a
+# mess of things like BigInt::numify and int(BigFloat). Using int($x->bstr)
+# seems to work.
+# E.g.:
+# $n = 33662485846146713; $n->numify; $n is now 3.36624858461467e+16
sub prime_get_config {
@@ -137,7 +146,7 @@ sub _validate_positive_integer {
croak "Parameter '$n' must be <= $max" if defined $max && $n > $max;
if ($n <= $_Config{'maxparam'}) {
$_[0] = $n->as_number() if ref($n) eq 'Math::BigFloat';
- $_[0] = $n->numify() if ref($n) eq 'Math::BigInt';
+ $_[0] = int($n->bstr) if ref($n) eq 'Math::BigInt';
} elsif (ref($n) ne 'Math::BigInt') {
croak "Parameter '$n' outside of integer range" if !defined $bigint::VERSION;
$_[0] = Math::BigInt->new("$n"); # Make $n a proper bigint object
@@ -321,6 +330,8 @@ sub primes {
$randbase = $randbase << 31;
# Now loop looking for a prime. There are lots of ways we could speed
# this up, especially for special cases.
+ # TODO: This has to be updated for 5.6.2. It keeps turning these numbers
+ # into floating point.
while (1) {
my $rand = $randbase + $irandf->(2147483648);
$prime = $low + ($rand % $range);
@@ -413,7 +424,7 @@ sub all_factors {
# to make sure we're using bigints when we calculate the product.
foreach my $f2 (@all) {
my $product = Math::BigInt->new("$f1") * Math::BigInt->new("$f2");
- $product = $product->numify if $product <= ~0;
+ $product = int($product->bstr) if $product <= ~0;
$all_factors{$product} = 1 if $product < $n;
}
}
@@ -498,11 +509,11 @@ sub euler_phi {
sub is_prime {
my($n) = @_;
- return 0 if $n <= 0; # everything else below 7 is composite
+ return 0 if $n <= 0;
_validate_positive_integer($n);
return _XS_is_prime($n) if $n <= $_XS_MAXVAL;
- return Math::Prime::Util::PP::is_prime($n);
+ return is_prob_prime($n);
}
sub next_prime {
@@ -565,7 +576,7 @@ sub is_strong_lucas_pseudoprime {
}
sub miller_rabin {
- warn "Use of miller_rabin is deprecated. Use is_strong_pseudoprime instead.";
+ #warn "miller_rabin() is deprecated. Use is_strong_pseudoprime instead.";
return is_strong_pseudoprime(@_);
}
@@ -651,7 +662,10 @@ sub prime_count_approx {
# return int( LogarithmicIntegral($x) - LogarithmicIntegral(sqrt($x))/2 );
- return int(RiemannR($x)+0.5);
+ my $result = RiemannR($x) + 0.5;
+
+ $result = Math::BigInt->new($result->bfloor->bstr()) if ref($result) eq 'Math::BigFloat';
+ return int($result);
}
sub prime_count_lower {
@@ -686,7 +700,9 @@ sub prime_count_lower {
elsif ($x <60000000000) { $a = 2.15; }
else { $a = 1.80; } # Dusart 1999, page 14
- return int( ($x/$flogx) * (1.0 + 1.0/$flogx + $a/($flogx*$flogx)) );
+ my $result = ($x/$flogx) * (1.0 + 1.0/$flogx + $a/($flogx*$flogx));
+ $result = Math::BigInt->new($result->bfloor->bstr()) if ref($result) eq 'Math::BigFloat';
+ return int($result);
}
sub prime_count_upper {
@@ -736,7 +752,12 @@ sub prime_count_upper {
elsif ($x <60000000000) { $a = 2.362; }
else { $a = 2.51; }
- return int( ($x/$flogx) * (1.0 + 1.0/$flogx + $a/($flogx*$flogx)) + 1.0 );
+ # Old versions of Math::BigFloat will do the Wrong Thing with this.
+ #return int( ($x/$flogx) * (1.0 + 1.0/$flogx + $a/($flogx*$flogx)) + 1.0 );
+ my $result = ($x/$flogx) * (1.0 + 1.0/$flogx + $a/($flogx*$flogx)) + 1.0;
+ $result = Math::BigInt->new($result->bfloor->bstr()) if ref($result) eq 'Math::BigFloat';
+ return int($result);
+
}
#############################################################################
@@ -881,7 +902,7 @@ sub LogarithmicIntegral {
return 0 if $n == 0;
croak("Invalid input to LogarithmicIntegral: x must be >= 0") if $n <= 0;
- if ( (defined $bignum::VERSION && (!defined &bignum::in_effect || bignum::in_effect())) || (ref($n) eq 'Math::BigFloat') ) {
+ if ( defined $bignum::VERSION || ref($n) eq 'Math::BigFloat' ) {
return Math::BigFloat->binf('-') if $n == 1;
return Math::BigFloat->new('1.045163780117492784844588889194613136522615578151201575832909144075013205210359530172717405626383356306') if $n == 2;
} else {
@@ -1206,26 +1227,26 @@ generate any primes. Uses the Cipolla 1902 approximation with two
polynomials, plus a correction term for small values to reduce the error.
-=head2 miller_rabin
+=head2 is_strong_pseudoprime
- my $maybe_prime = miller_rabin($n, 2);
- my $probably_prime = miller_rabin($n, 2, 3, 5, 7, 11, 13, 17);
+ my $maybe_prime = is_strong_pseudoprime($n, 2);
+ my $probably_prime = is_strong_pseudoprime($n, 2, 3, 5, 7, 11, 13, 17);
Takes a positive number as input and one or more bases. The bases must be
-between C<2> and C<n - 2>. Returns 2 is C<n> is definitely prime, 1 if C<n>
-is probably prime, and 0 if C<n> is definitely composite. A value of 2 will
-only be returned for the inputs of 2 and 3, which are shortcut.
+between C<2> and C<n - 2>. Returns 1 is C<n> is a prime or a strong
+pseudoprime to all of the bases, and 0 if not.
If 0 is returned, then the number really is a composite. If 1 is returned,
then it is either a prime or a strong pseudoprime to all the given bases.
-Given enough bases, the chances become very, very strong that the number is
-actually prime.
+Given enough distinct bases, the chances become very, very strong that the
+number is actually prime.
This is usually used in combination with other tests to make either stronger
tests (e.g. the strong BPSW test) or deterministic results for numbers less
-than some verified limit (such as the C<is_prob_prime> function in this module).
-However, given the chances of passing multiple bases, there are some math
-packages that just use multiple MR tests for primality testing.
+than some verified limit (e.g. it has long been known that no more than three
+selected bases are required to give correct primality test results for any
+32-bit number). Given the small chances of passing multiple bases, there
+are some math packages that just use multiple MR tests for primality testing.
Even numbers other than 2 will always return 0 (composite). While the
algorithm does run with even input, most sources define it only on odd input.
@@ -1233,13 +1254,18 @@ Returning composite for all non-2 even input makes the function match most
other implementations including L<Math::Primality>'s C<is_strong_pseudoprime>
function.
+=head2 miller_rabin
+
+An alias for C<is_strong_pseudoprime>. This name is being deprecated.
+
=head2 is_strong_lucas_pseudoprime
Takes a positive number as input, and returns 1 if the input is a strong
-Lucas pseudoprime using the Selfridge method of choosing D, P, and Q (hence
-some sources call this a strong Lucas-Selfridge pseudoprime). This is one
-half of the BPSW primality test (the Miller-Rabin test being the other).
+Lucas pseudoprime using the Selfridge method of choosing D, P, and Q (some
+sources call this a strong Lucas-Selfridge pseudoprime). This is one half
+of the BPSW primality test (the Miller-Rabin strong pseudoprime test with
+base 2 being the other half).
=head2 is_prob_prime
@@ -1560,7 +1586,7 @@ Accuracy should be at least 14 digits.
Print pseudoprimes base 17:
- perl -MMath::Prime::Util=:all -E 'my $n=$base|1; while(1) { print "$n " if miller_rabin($n,$base) && !is_prime($n); $n+=2; } BEGIN {$|=1; $base=17}'
+ perl -MMath::Prime::Util=:all -E 'my $n=$base|1; while(1) { print "$n " if is_strong_pseudoprime($n,$base) && !is_prime($n); $n+=2; } BEGIN {$|=1; $base=17}'
Print some primes above 64-bit range:
diff --git a/lib/Math/Prime/Util/PP.pm b/lib/Math/Prime/Util/PP.pm
index ff39e71..18ef19a 100644
--- a/lib/Math/Prime/Util/PP.pm
+++ b/lib/Math/Prime/Util/PP.pm
@@ -61,7 +61,7 @@ sub _validate_positive_integer {
croak "Parameter '$n' must be <= $max" if defined $max && $n > $max;
if ($n <= ~0) {
$_[0] = $n->as_number() if ref($n) eq 'Math::BigFloat';
- $_[0] = $n->numify() if ref($n) eq 'Math::BigInt';
+ $_[0] = int($n->bstr) if ref($n) eq 'Math::BigInt';
} elsif (ref($n) ne 'Math::BigInt') {
croak "Parameter '$n' outside of integer range" if !defined $bigint::VERSION;
$_[0] = Math::BigInt->new("$n"); # Make $n a proper bigint object
@@ -643,7 +643,7 @@ sub is_strong_lucas_pseudoprime {
# Check for perfect square
if (ref($n) eq 'Math::BigInt') {
- my $mcheck = ($n & 127)->numify;
+ my $mcheck = int(($n & 127)->bstr);
if (($mcheck*0x8bc40d7d) & ($mcheck*0xa1e2f5d1) & 0x14020a) {
# ~82% of non-squares were rejected by the bloom filter
my $sq = int($n->copy()->bsqrt());
@@ -659,8 +659,6 @@ sub is_strong_lucas_pseudoprime {
# Determine Selfridge D, P, and Q parameters
my $D = _find_jacobi_d_sequence($n);
- #die if ref($D) eq 'Math::BigInt';
- #$D = $D->numify if $D <= ~0 && ref($D) eq 'Math::BigInt';
return 0 if $D == 0; # We found a divisor in the sequence
my $P = 1;
my $Q = int( (1 - $D) / 4 );
@@ -694,7 +692,7 @@ sub is_strong_lucas_pseudoprime {
my $Q_m2 = ($ZERO + $Q) * 2;
my $Qkd = $ZERO + $Q;
$d >>= 1;
- $d = $d->numify if $d <= ~0 && ref($d) eq 'Math::BigInt';
+ $d = int($d->bstr) if $d <= ~0 && ref($d) eq 'Math::BigInt';
#my $i = 1;
while ($d != 0) {
#warn "U=$U V=$V Qm=$Q_m Qm2=$Q_m2\n";
@@ -734,9 +732,6 @@ sub is_strong_lucas_pseudoprime {
$Qkd2 = 2 * $Qkd;
}
}
- #warn "Math::BigInt is loaded\n" if defined $Math::BigInt::VERSION;
- #warn "bigint is loaded\n" if defined $bigint::VERSION;
- #warn "bigint is in effect\n" if bigint::in_effect();
return 0;
}
@@ -751,7 +746,7 @@ sub _basic_factor {
while ( ($_[0] % 5) == 0 ) { push @factors, 5; $_[0] = int($_[0] / 5); }
# Stop using bignum if we can
- $_[0] = $_[0]->numify if ref($_[0]) eq 'Math::BigInt' && $_[0] <= ~0;
+ $_[0] = int($_[0]->bstr) if ref($_[0]) eq 'Math::BigInt' && $_[0] <= ~0;
if ( ($_[0] > 1) && _is_prime7($_[0]) ) {
push @factors, $_[0];
@@ -809,7 +804,7 @@ sub factor {
while (@nstack) {
$n = pop @nstack;
# Don't use bignum on $n if it has gotten small enough.
- $n = $n->numify if ref($n) eq 'Math::BigInt' && $n <= ~0;
+ $n = int($n->bstr) if ref($n) eq 'Math::BigInt' && $n <= ~0;
#print "Looking at $n with stack ", join(",", at nstack), "\n";
while ( ($n >= (31*31)) && !_is_prime7($n) ) {
my @ftry;
@@ -1299,9 +1294,9 @@ the prime_count and nth prime, next_prime and prev_prime, factoring utilities,
and more.
All routines should work with native integers or multi-precision numbers. To
-enable big numbers, use bignum:
+enable big numbers, use bigint or bignum:
- use bignum;
+ use bigint;
say prime_count_approx(1000000000000000000000000)'
# says 18435599767347543283712
diff --git a/t/12-nextprime.t b/t/12-nextprime.t
index 5f15a4b..9db936d 100644
--- a/t/12-nextprime.t
+++ b/t/12-nextprime.t
@@ -72,7 +72,11 @@ is( prev_prime(19660), 19609, "prev prime of 19660 is 19609" );
is( prev_prime(19610), 19609, "prev prime of 19610 is 19609" );
is( prev_prime(2), 0, "Previous prime of 2 returns 0" );
-is( next_prime(~0-4), 0, "Next prime of ~0-4 returns 0" );
+if ($use64) {
+ is( next_prime(18446744073709551611), 0, "Next prime of ~0-4 returns 0" );
+} else {
+ is( next_prime(4294967291), 0, "Next prime of ~0-4 returns 0" );
+}
# Turns out the testing of prev/next from 0-3572 still misses some cases.
foreach my $n (2010733 .. 2010880) {
diff --git a/t/14-nthprime.t b/t/14-nthprime.t
index a5c1840..6979b9d 100644
--- a/t/14-nthprime.t
+++ b/t/14-nthprime.t
@@ -96,20 +96,21 @@ if ($use64) {
}
}
-my $maxindex = $use64 ? 425656284035217743 : 203280221;
-my $maxprime = $use64 ? 18446744073709551557 : 4294967291;
+my $maxindex = $use64 ? 425656284035217743 : 203280221;
+my $maxindexp1 = $use64 ? 425656284035217744 : 203280222;
+my $maxprime = $use64 ? 18446744073709551557 : 4294967291;
cmp_ok( nth_prime_lower($maxindex), '<=', $maxprime, "nth_prime_lower(maxindex) <= maxprime");
cmp_ok( nth_prime_upper($maxindex), '>=', $maxprime, "nth_prime_upper(maxindex) >= maxprime");
cmp_ok( nth_prime_approx($maxindex), '==', $maxprime, "nth_prime_approx(maxindex) == maxprime");
-cmp_ok( nth_prime_lower($maxindex+1), '>=', nth_prime_lower($maxindex), "nth_prime_lower(maxindex+1) >= nth_prime_lower(maxindex)");
+cmp_ok( nth_prime_lower($maxindexp1), '>=', nth_prime_lower($maxindex), "nth_prime_lower(maxindex+1) >= nth_prime_lower(maxindex)");
-my $add = ($broken64) ? 100 : 1;
+my $overindex = ($broken64) ? 425656284035217843 : $maxindexp1;
-eval { nth_prime_upper($maxindex+$add); };
-like($@, qr/overflow/, "nth_prime_upper(maxindex+$add) overflows");
+eval { nth_prime_upper($overindex); };
+like($@, qr/overflow/, "nth_prime_upper($overindex) overflows");
-eval { nth_prime_approx($maxindex+$add); };
-like($@, qr/overflow/, "nth_prime_approx(maxindex+$add) overflows");
+eval { nth_prime_approx($overindex); };
+like($@, qr/overflow/, "nth_prime_approx($overindex) overflows");
-eval { nth_prime($maxindex+$add); };
-like($@, qr/overflow/, "nth_prime(maxindex+$add) overflows");
+eval { nth_prime($overindex); };
+like($@, qr/overflow/, "nth_prime($overindex) overflows");
diff --git a/t/80-pp.t b/t/80-pp.t
index 79d4f52..4f05556 100644
--- a/t/80-pp.t
+++ b/t/80-pp.t
@@ -5,6 +5,7 @@ use warnings;
# This is a subset of our tests. You really should run the whole test suite
# on the PP code. What this will do is basic regression testing.
my $extra = defined $ENV{RELEASE_TESTING} && $ENV{RELEASE_TESTING};
+my $use64 = ~0 > 4294967295;
use Test::More;
my @small_primes = qw/
@@ -294,7 +295,11 @@ is( prev_prime(19660), 19609, "prev prime of 19660 is 19609" );
is( prev_prime(19610), 19609, "prev prime of 19610 is 19609" );
is( prev_prime(2), 0, "Previous prime of 2 returns 0" );
-is( next_prime(~0-4), 0, "Next prime of ~0-4 returns 0" );
+if ($use64) {
+ is( next_prime(18446744073709551611), 0, "Next prime of ~0-4 returns 0" );
+} else {
+ is( next_prime(4294967291), 0, "Next prime of ~0-4 returns 0" );
+}
foreach my $n (2010733 .. 2010880) {
is(next_prime($n), 2010881, "next_prime($n) == 2010881");
diff --git a/t/81-bignum.t b/t/81-bignum.t
index 25cc484..802d79f 100644
--- a/t/81-bignum.t
+++ b/t/81-bignum.t
@@ -77,19 +77,15 @@ use Math::Prime::Util qw/
ExponentialIntegral
LogarithmicIntegral
RiemannR
+ primes
+ prime_count
+ nth_prime
+ is_prime
+ next_prime
+ prev_prime
+ is_strong_pseudoprime
/;
- *primes = \&Math::Prime::Util::PP::primes;
-
- *prime_count = \&Math::Prime::Util::PP::prime_count;
- #*nth_prime = \&Math::Prime::Util::PP::nth_prime;
-
- *is_prime = \&Math::Prime::Util::PP::is_prime;
- *next_prime = \&Math::Prime::Util::PP::next_prime;
- *prev_prime = \&Math::Prime::Util::PP::prev_prime;
-
- *miller_rabin = \&Math::Prime::Util::PP::miller_rabin;
-
###############################################################################
foreach my $n (@primes) {
@@ -125,7 +121,7 @@ is( prime_count(877777777777777777777752, 877777777777777777777872), 2, "prime_c
while (my($psrp, $baseref) = each (%pseudoprimes)) {
foreach my $base (@$baseref) {
- ok( miller_rabin($psrp, $base), "$psrp is a strong pseudoprime to base $base" );
+ ok( is_strong_pseudoprime($psrp, $base), "$psrp is a strong pseudoprime to base $base" );
}
}
--
Alioth's /usr/local/bin/git-commit-notice on /srv/git.debian.org/git/pkg-perl/packages/libmath-prime-util-perl.git
More information about the Pkg-perl-cvs-commits
mailing list