[libmath-prime-util-perl] 11/16: Allow random_prime to use overridden rand

[Partha P. Mukherjee]

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

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

commit 2c3cc18d8ca3123954859f6bb637890df51724eb
Author: Dana Jacobsen <dana at acm.org>
Date:   Mon Jun 11 17:11:37 2012 -0600

    Allow random_prime to use overridden rand
---
 Changes                |  1 +
 lib/Math/Prime/Util.pm | 24 +++++++++++++++++++++---
 t/16-randomprime.t     |  3 +++
 util.c                 |  2 +-
 4 files changed, 26 insertions(+), 4 deletions(-)

diff --git a/Changes b/Changes
index 8e6a7e6..8fb3bc0 100644
--- a/Changes
+++ b/Changes
@@ -10,6 +10,7 @@ Revision history for Perl extension Math::Prime::Util.
             say prime_count( 10**16,  10**16 + 2**20 )
     - Add Ei(x), li(x), and R(x) functions.
     - prime_count_approx uses R(x), making it vastly more accurate.
+    - Let user override rand for random_prime.
 
 0.04  7 June 2012
     - Didn't do tests on 32-bit machine before release.  Test suite caught
diff --git a/lib/Math/Prime/Util.pm b/lib/Math/Prime/Util.pm
index 1801149..55619f9 100644
--- a/lib/Math/Prime/Util.pm
+++ b/lib/Math/Prime/Util.pm
@@ -124,23 +124,31 @@ sub random_prime {
   my $range = $high - $low + 1;
   my $prime;
 
+  # Note:  I was using rand($range), but Math::Random::MT ignores the argument
+  #        instead of following its documentation.
+  my $irandf = (exists &::rand) ? sub { return int(::rand()*shift); }
+                                : sub { return int(rand()*shift); };
+
   if ($high < 30000) {
     # nice deterministic solution, but gets very costly with large values.
     my $li = ($low == 2) ? 1 : prime_count($low);
     my $hi = prime_count($high);
     my $irange = $hi - $li + 1;
-    my $rand = int(rand($irange));
+    my $rand = $irandf->($irange);
     $prime = nth_prime($li + $rand);
   } else {
     # random loop
+    my $loop_limit = 2000 * 1000;  # To protect against broken rand
     if ($range <= 4294967295) {
       do {
-        $prime = $low + int(rand($range));
+        $prime = $low + $irandf->($range);
+        croak "Random function broken?" if $loop_limit-- < 0;
       }  while ( !($prime%2) || !($prime%3) || !is_prime($prime) );
     } else {
       do {
-        my $rand = ( (int(rand(4294967295)) << 32) + int(rand(4294967295)) ) % $range;
+        my $rand = ( ($irandf->(4294967295) << 32) + $irandf->(4294967295) ) % $range;
         $prime = $low + $rand;
+        croak "Random function broken?" if $loop_limit-- < 0;
       }  while ( !($prime%2) || !($prime%3) || !is_prime($prime) );
     }
   }
@@ -505,6 +513,16 @@ L<Crypt::Primes> or something similar.  The current L<Math::Prime::Util>
 module does not use strong randomness, and its primes are ridiculously small
 by cryptographic standards.
 
+Perl's L<rand> function is normally called, but if the sub C<main::rand>
+exists, it will be used instead.  When called with no arguments it should
+return a float value between 0 and 1-epsilon, with 32 bits of randomness.
+Examples:
+
+  # Use Mersenne Twister
+  use Math::Random::MT::Auto qw/rand/;
+
+  # Use my custom random function
+  sub rand { ... }
 
 =head2 random_ndigit_prime
 
diff --git a/t/16-randomprime.t b/t/16-randomprime.t
index bb51aba..19c3b44 100644
--- a/t/16-randomprime.t
+++ b/t/16-randomprime.t
@@ -3,6 +3,9 @@ use strict;
 use warnings;
 
 use Test::More;
+#use Math::Random::MT qw/rand/;
+#use Math::Random::MT::Auto qw/rand/;
+#sub rand { return 0.5; }
 use Math::Prime::Util qw/random_prime random_ndigit_prime is_prime/;
 
 my $use64 = Math::Prime::Util::_maxbits > 32;
diff --git a/util.c b/util.c
index 39f5413..b59ec1d 100644
--- a/util.c
+++ b/util.c
@@ -713,7 +713,7 @@ UV nth_prime_approx(UV n)
   else if (n < 12000) approx += 3.0 * order;
   else if (n <150000) approx += 2.1 * order;
   else if (n <200000000) approx += 0.0 * order;
-  else                approx += -0.023 * order;
+  else                approx += -0.010 * order; /* -0.025 is better */
 
   /* For all three analytical functions, it is possible that for a given valid
    * input, we will not be able to return an output that fits in the UV type.

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