[libmath-prime-util-perl] 17/43: Add nth_twin_prime
Partha P. Mukherjee
ppm-guest at moszumanska.debian.org
Thu May 21 18:53:06 UTC 2015
This is an automated email from the git hooks/post-receive script.
ppm-guest pushed a commit to annotated tag v0.40
in repository libmath-prime-util-perl.
commit 9e1d2afcd67c24306896cc96f00111d95c0d113a
Author: Dana Jacobsen <dana at acm.org>
Date: Thu Mar 27 18:13:25 2014 -0700
Add nth_twin_prime
---
Changes | 1 +
XS.xs | 27 +++++++++++++++------------
lib/Math/Prime/Util.pm | 8 +++++++-
t/14-nthprime.t | 12 ++++++++++--
util.c | 41 +++++++++++++++++++++++++++++++++++++++++
util.h | 1 +
6 files changed, 75 insertions(+), 15 deletions(-)
diff --git a/Changes b/Changes
index af0b19a..713c556 100644
--- a/Changes
+++ b/Changes
@@ -8,6 +8,7 @@ Revision history for Perl module Math::Prime::Util
- random_shawe_taylor_prime_with_cert as above with certificate
- twin_prime_count counts twin primes in range
- twin_prime_count_approx fast approximation to Pi_2(n)
+ - nth_twin_prime returns the nth twin prime
[FUNCTIONALITY AND PERFORMANCE]
diff --git a/XS.xs b/XS.xs
index cd9bfb8..1cbabd2 100644
--- a/XS.xs
+++ b/XS.xs
@@ -603,10 +603,11 @@ next_prime(IN SV* svn)
nth_prime_upper = 3
nth_prime_lower = 4
nth_prime_approx = 5
- prime_count_upper = 6
- prime_count_lower = 7
- prime_count_approx = 8
- twin_prime_count_approx = 9
+ nth_twin_prime = 6
+ prime_count_upper = 7
+ prime_count_lower = 8
+ prime_count_approx = 9
+ twin_prime_count_approx = 10
PPCODE:
if (_validate_int(aTHX_ svn, 0)) {
UV n = my_svuv(svn);
@@ -622,10 +623,11 @@ next_prime(IN SV* svn)
case 3: ret = nth_prime_upper(n); break;
case 4: ret = nth_prime_lower(n); break;
case 5: ret = nth_prime_approx(n); break;
- case 6: ret = prime_count_upper(n); break;
- case 7: ret = prime_count_lower(n); break;
- case 8 :ret = prime_count_approx(n); break;
- case 9:
+ case 6: ret = nth_twin_prime(n); break;
+ case 7: ret = prime_count_upper(n); break;
+ case 8: ret = prime_count_lower(n); break;
+ case 9 :ret = prime_count_approx(n); break;
+ case 10:
default:ret = twin_prime_count_approx(n); break;
}
XSRETURN_UV(ret);
@@ -642,10 +644,11 @@ next_prime(IN SV* svn)
case 3: _vcallsub_with_pp("nth_prime_upper"); break;
case 4: _vcallsub_with_pp("nth_prime_lower"); break;
case 5: _vcallsub_with_pp("nth_prime_approx"); break;
- case 6: _vcallsub_with_pp("prime_count_upper"); break;
- case 7: _vcallsub_with_pp("prime_count_lower"); break;
- case 8: _vcallsub_with_pp("prime_count_approx"); break;
- case 9:
+ case 6: _vcallsub_with_pp("nth_twin_prime"); break;
+ case 7: _vcallsub_with_pp("prime_count_upper"); break;
+ case 8: _vcallsub_with_pp("prime_count_lower"); break;
+ case 9: _vcallsub_with_pp("prime_count_approx"); break;
+ case 10:
default: _vcallsub_with_pp("twin_prime_count_approx"); break;
}
return; /* skip implicit PUTBACK */
diff --git a/lib/Math/Prime/Util.pm b/lib/Math/Prime/Util.pm
index f3785b6..4172777 100644
--- a/lib/Math/Prime/Util.pm
+++ b/lib/Math/Prime/Util.pm
@@ -33,7 +33,7 @@ our @EXPORT_OK =
prime_count
prime_count_lower prime_count_upper prime_count_approx
nth_prime nth_prime_lower nth_prime_upper nth_prime_approx
- twin_prime_count twin_prime_count_approx
+ twin_prime_count twin_prime_count_approx nth_twin_prime
random_prime random_ndigit_prime random_nbit_prime random_strong_prime
random_proven_prime random_proven_prime_with_cert
random_maurer_prime random_maurer_prime_with_cert
@@ -1368,6 +1368,12 @@ generate any primes. Uses the Cipolla 1902 approximation with two
polynomials, plus a correction for small values to reduce the error.
+=head2 nth_twin_prime
+
+Returns the Nth twin prime. This is done via sieving and counting, so
+is not very fast for large values.
+
+
=head2 is_pseudoprime
Takes a positive number C<n> and a base C<a> as input, and returns 1 if
diff --git a/t/14-nthprime.t b/t/14-nthprime.t
index 11632bb..c026157 100644
--- a/t/14-nthprime.t
+++ b/t/14-nthprime.t
@@ -3,7 +3,8 @@ use strict;
use warnings;
use Test::More;
-use Math::Prime::Util qw/primes nth_prime nth_prime_lower nth_prime_upper nth_prime_approx/;
+use Math::Prime::Util qw/primes nth_prime nth_twin_prime
+ nth_prime_lower nth_prime_upper nth_prime_approx/;
my $use64 = Math::Prime::Util::prime_get_config->{'maxbits'} > 32;
my $usexs = Math::Prime::Util::prime_get_config->{'xs'};
@@ -66,7 +67,8 @@ plan tests => 0 + 2*scalar(keys %pivals32)
+ 3*scalar(keys %nthprimes32)
+ scalar(keys %nthprimes_small)
+ $use64 * 3 * scalar(keys %nthprimes64)
- + 3
+ + 3 # nth_prime_lower with max index
+ + 3 # nth_twin_prime
+ (($extra && $use64 && $usexs) ? 1 : 0);
@@ -117,3 +119,9 @@ if ($extra && $use64 && $usexs) {
# Test an nth prime value that uses the binary-search-on-R(n) algorithm
is( nth_prime(21234567890), 551990503367, "nth_prime(21234567890)" );
}
+
+####################################3
+
+is( nth_twin_prime(0), 0, "nth_twin_prime(0) = 0" );
+is( nth_twin_prime(17), 239, "239 = 17th twin prime" );
+is( nth_twin_prime(1234), 101207, "101207 = 1234'th twin prime" );
diff --git a/util.c b/util.c
index b03c4e2..888d98c 100644
--- a/util.c
+++ b/util.c
@@ -938,6 +938,47 @@ UV twin_prime_count(UV beg, UV end)
return sum;
}
+UV nth_twin_prime(UV n)
+{
+ unsigned char* segment;
+ UV nth = 0;
+ UV beg, end;
+
+ if (n < 6) {
+ switch (n) {
+ case 0: nth = 0; break;
+ case 1: nth = 3; break;
+ case 2: nth = 5; break;
+ case 3: nth =11; break;
+ case 4: nth =17; break;
+ case 5:
+ default: nth =29; break;
+ }
+ return nth;
+ }
+ n -= 5;
+ beg = 31;
+ end = UV_MAX - 16;
+ {
+ UV seg_base, seg_low, seg_high;
+ void* ctx = start_segment_primes(beg, end, &segment);
+ while (next_segment_primes(ctx, &seg_base, &seg_low, &seg_high)) {
+ UV p, bytes = (seg_high-seg_low+29)/30;
+ for (p = 0; p < bytes; p++) {
+ UV s = segment[p];
+ int twin1 = !(s & 0x0C);
+ int twin2 = !(s & 0x30);
+ int twin3 = !(s & 0x80) && ((p+1 < bytes) ? !(segment[p+1] & 0x01) : _XS_is_prime(seg_high+2));
+ if (twin1 && !--n) { nth=seg_base+p*30+11; break; }
+ if (twin2 && !--n) { nth=seg_base+p*30+17; break; }
+ if (twin3 && !--n) { nth=seg_base+p*30+29; break; }
+ }
+ if (n == 0) break;
+ }
+ end_segment_primes(ctx);
+ }
+ return nth;
+}
/* Return a char array with lo-hi+1 elements. mu[k-lo] = µ(k) for k = lo .. hi.
diff --git a/util.h b/util.h
index e730a4b..2130671 100644
--- a/util.h
+++ b/util.h
@@ -22,6 +22,7 @@ extern UV prime_count_upper(UV x);
extern UV prime_count_approx(UV x);
extern UV twin_prime_count(UV low, UV high);
extern UV twin_prime_count_approx(UV n);
+extern UV nth_twin_prime(UV n);
extern int powerof(UV n);
extern int is_power(UV n, UV a);
--
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