[libmath-prime-util-perl] 02/04: Change PP perfect square test

[Partha P. Mukherjee]

This is an automated email from the git hooks/post-receive script.

ppm-guest pushed a commit to annotated tag v0.20
in repository libmath-prime-util-perl.

commit dc5a392af4c0fcb42482981cf5357b9a6f40fd71
Author: Dana Jacobsen <dana at acm.org>
Date:   Sat Feb 2 19:04:37 2013 -0800

    Change PP perfect square test
---
 Changes                          |  3 +++
 MANIFEST                         |  1 +
 examples/bench-mp-prime_count.pl | 21 +++++++++++++++++++++
 lib/Math/Prime/Util/PP.pm        | 19 ++++++++-----------
 4 files changed, 33 insertions(+), 11 deletions(-)

diff --git a/Changes b/Changes
index 904b316..11d6ef6 100644
--- a/Changes
+++ b/Changes
@@ -4,6 +4,9 @@ Revision history for Perl extension Math::Prime::Util.
 
     - Little speedup for PP AKS, and turn off test on 32-bit machines.
 
+    - Replaced fast sqrt detection in PP.pm with a slightly slower version.
+      The bloom filter doesn't work right in 32-bit Perl.
+
 0.19  1 February 2012
 
     - Update MR bases with newest from http://miller-rabin.appspot.com/.
diff --git a/MANIFEST b/MANIFEST
index 6f59e04..6945e23 100644
--- a/MANIFEST
+++ b/MANIFEST
@@ -41,6 +41,7 @@ examples/bench-pp-isprime.pl
 examples/bench-pp-sieve.pl
 examples/bench-mp-nextprime.pl
 examples/bench-mp-psrp.pl
+examples/bench-mp-prime_count.pl
 examples/test-factor-yafu.pl
 examples/test-factor-mpxs.pl
 examples/test-nextprime-yafu.pl
diff --git a/examples/bench-mp-prime_count.pl b/examples/bench-mp-prime_count.pl
new file mode 100644
index 0000000..d19797f
--- /dev/null
+++ b/examples/bench-mp-prime_count.pl
@@ -0,0 +1,21 @@
+#!/usr/bin/env perl
+use strict;
+use warnings;
+use Math::Prime::Util;
+use Math::Prime::Util::GMP;
+use Math::Primality;
+use Benchmark qw/:all/;
+my $count = shift || -2;
+
+#my($n, $exp) = (100000,9592);
+my($n, $exp) = (1000000,78498);
+#my($n, $exp) = (10000000,664579);
+cmpthese($count,{
+  'MP'           =>sub { die unless $exp == Math::Primality::prime_count($n); },
+  'MPU default'  =>sub { die unless $exp == Math::Prime::Util::prime_count($n); },
+  'MPU XS Sieve' =>sub { die unless $exp == Math::Prime::Util::_XS_prime_count($n); },
+  'MPU XS Lehmer'=>sub { die unless $exp == Math::Prime::Util::_XS_lehmer_pi($n); },
+  'MPU PP Sieve' =>sub { die unless $exp == Math::Prime::Util::PP::_sieve_prime_count($n); },
+  'MPU PP Lehmer'=>sub { die unless $exp == Math::Prime::Util::PP::_lehmer_pi($n); },
+  'MPU GMP Trial'=>sub { die unless $exp == Math::Prime::Util::GMP::prime_count(2,$n); },
+});
diff --git a/lib/Math/Prime/Util/PP.pm b/lib/Math/Prime/Util/PP.pm
index 051422e..54c8b10 100644
--- a/lib/Math/Prime/Util/PP.pm
+++ b/lib/Math/Prime/Util/PP.pm
@@ -852,16 +852,15 @@ sub is_strong_lucas_pseudoprime {
 
   # Check for perfect square
   if (ref($n) eq 'Math::BigInt') {
-    my $mcheck = int(($n & 127)->bstr);
-    if (($mcheck*0x8bc40d7d) & ($mcheck*0xa1e2f5d1) & 0x14020a) {
-      # ~82% of non-squares were rejected by the bloom filter
+    my $mc = int(($n & 31)->bstr);
+    if ($mc==0||$mc==1||$mc==4||$mc==9||$mc==16||$mc==17||$mc==25) {
       my $sq = $n->copy->bsqrt->bfloor;
       $sq->bmul($sq);
       return 0 if $sq == $n;
     }
   } else {
-    my $mcheck = $n & 127;
-    if (($mcheck*0x8bc40d7d) & ($mcheck*0xa1e2f5d1) & 0x14020a) {
+    my $mc = $n & 31;
+    if ($mc==0||$mc==1||$mc==4||$mc==9||$mc==16||$mc==17||$mc==25) {
       my $sq = int(sqrt($n));
       return 0 if ($sq*$sq) == $n;
     }
@@ -1412,9 +1411,8 @@ sub holf_factor {
       $s->binc if ($s * $s) != $ni;
       my $m = $s->copy->bmul($s)->bmod($n);
       # Check for perfect square
-      my $mcheck = int(($m & 127)->bstr);
-      next if (($mcheck*0x8bc40d7d) & ($mcheck*0xa1e2f5d1) & 0x14020a);
-      # ... 82% of non-squares were rejected by the bloom filter
+      my $mc = int(($m & 31)->bstr);
+      next unless $mc==0||$mc==1||$mc==4||$mc==9||$mc==16||$mc==17||$mc==25;
       my $f = $m->copy->bsqrt->bfloor->as_int;
       next unless ($f*$f) == $m;
       $f = Math::BigInt::bgcd( ($s > $f) ? $s-$f : $f-$s,  $n);
@@ -1432,9 +1430,8 @@ sub holf_factor {
       $s++ if ($s * $s) != ($n * $i);
       my $m = ($s < $_half_word) ? ($s*$s) % $n : _mulmod($s, $s, $n);
       # Check for perfect square
-      my $mcheck = $m & 127;
-      next if (($mcheck*0x8bc40d7d) & ($mcheck*0xa1e2f5d1) & 0x14020a);
-      # ... 82% of non-squares were rejected by the bloom filter
+      my $mc = $m & 31;
+      next unless $mc==0||$mc==1||$mc==4||$mc==9||$mc==16||$mc==17||$mc==25;
       my $f = int(sqrt($m));
       next unless $f*$f == $m;
       $f = _gcd_ui( ($s > $f)  ?  $s - $f  :  $f - $s,  $n);

-- 
Alioth's /usr/local/bin/git-commit-notice on /srv/git.debian.org/git/pkg-perl/packages/libmath-prime-util-perl.git