[libmath-prime-util-perl] 11/13: Change to unified method for small is_prime, next_prime, prev_prime, prime_count
Partha P. Mukherjee
ppm-guest at moszumanska.debian.org
Thu May 21 18:47:42 UTC 2015
This is an automated email from the git hooks/post-receive script.
ppm-guest pushed a commit to annotated tag v0.24
in repository libmath-prime-util-perl.
commit 49512f3d4bf90c7cc70a847a18975976ba5b9e06
Author: Dana Jacobsen <dana at acm.org>
Date: Sun Mar 10 19:04:01 2013 -0700
Change to unified method for small is_prime, next_prime, prev_prime, prime_count
---
util.c | 103 +++++++++++++++++++++++++++++++++++++++++++----------------------
1 file changed, 69 insertions(+), 34 deletions(-)
diff --git a/util.c b/util.c
index 123349b..3eb4929 100644
--- a/util.c
+++ b/util.c
@@ -128,17 +128,36 @@ static int _is_prime7(UV n)
}
-/* Marked bits for each n, indicating if the number is prime */
-static const unsigned char prime_is_small[] =
- {0xac,0x28,0x8a,0xa0,0x20,0x8a,0x20,0x28,0x88,0x82,0x08,0x02,0xa2,0x28,0x02,
- 0x80,0x08,0x0a,0xa0,0x20,0x88,0x20,0x28,0x80,0xa2,0x00,0x08,0x80,0x28,0x82,
- 0x02,0x08,0x82,0xa0,0x20,0x0a,0x20,0x00,0x88,0x22,0x00,0x08,0x02,0x28,0x82,
- 0x80,0x20,0x88,0x20,0x20,0x02,0x02,0x28,0x80,0x82,0x08,0x02,0xa2,0x08,0x80,
- 0x80,0x08,0x88,0x20,0x00,0x0a,0x00,0x20,0x08,0x20,0x08,0x0a,0x02,0x08,0x82,
- 0x82,0x20,0x0a,0x80,0x00,0x8a,0x20,0x28,0x00,0x22,0x08,0x08,0x20,0x20,0x80,
- 0x80,0x20,0x88,0x80,0x20,0x02,0x22,0x00,0x08,0x20,0x00,0x0a,0xa0,0x28,0x80,
- 0x00,0x20,0x8a,0x00,0x20,0x8a,0x00,0x00,0x88,0x80,0x00,0x02,0x22,0x08,0x02};
-#define NPRIME_IS_SMALL (sizeof(prime_is_small)/sizeof(prime_is_small[0]))
+/* We'll use this little static sieve to quickly answer small values of
+ * is_prime, next_prime, prev_prime, prime_count
+ * for non-threaded Perl it's basically the same as getting the primary
+ * cache. It guarantees we'll have an answer with no waiting on any version.
+ */
+static const unsigned char prime_sieve30[] =
+ {0x01,0x20,0x10,0x81,0x49,0x24,0xc2,0x06,0x2a,0xb0,0xe1,0x0c,0x15,0x59,0x12,
+ 0x61,0x19,0xf3,0x2c,0x2c,0xc4,0x22,0xa6,0x5a,0x95,0x98,0x6d,0x42,0x87,0xe1,
+ 0x59,0xa9,0xa9,0x1c,0x52,0xd2,0x21,0xd5,0xb3,0xaa,0x26,0x5c,0x0f,0x60,0xfc,
+ 0xab,0x5e,0x07,0xd1,0x02,0xbb,0x16,0x99,0x09,0xec,0xc5,0x47,0xb3,0xd4,0xc5,
+ 0xba,0xee,0x40,0xab,0x73,0x3e,0x85,0x4c,0x37,0x43,0x73,0xb0,0xde,0xa7,0x8e,
+ 0x8e,0x64,0x3e,0xe8,0x10,0xab,0x69,0xe5,0xf7,0x1a,0x7c,0x73,0xb9,0x8d,0x04,
+ 0x51,0x9a,0x6d,0x70,0xa7,0x78,0x2d,0x6d,0x27,0x7e,0x9a,0xd9,0x1c,0x5f,0xee,
+ 0xc7,0x38,0xd9,0xc3,0x7e,0x14,0x66,0x72,0xae,0x77,0xc1,0xdb,0x0c,0xcc,0xb2,
+ 0xa5,0x74,0xe3,0x58,0xd5,0x4b,0xa7,0xb3,0xb1,0xd9,0x09,0xe6,0x7d,0x23,0x7c,
+ 0x3c,0xd3,0x0e,0xc7,0xfd,0x4a,0x32,0x32,0xfd,0x4d,0xb5,0x6b,0xf3,0xa8,0xb3,
+ 0x85,0xcf,0xbc,0xf4,0x0e,0x34,0xbb,0x93,0xdb,0x07,0xe6,0xfe,0x6a,0x57,0xa3,
+ 0x8c,0x15,0x72,0xdb,0x69,0xd4,0xaf,0x59,0xdd,0xe1,0x3b,0x2e,0xb7,0xf9,0x2b,
+ 0xc5,0xd0,0x8b,0x63,0xf8,0x95,0xfa,0x77,0x40,0x97,0xea,0xd1,0x9f,0xaa,0x1c,
+ 0x48,0xae,0x67,0xf7,0xeb,0x79,0xa5,0x55,0xba,0xb2,0xb6,0x8f,0xd8,0x2d,0x6c,
+ 0x2a,0x35,0x54,0xfd,0x7c,0x9e,0xfa,0xdb,0x31,0x78,0xdd,0x3d,0x56,0x52,0xe7,
+ 0x73,0xb2,0x87,0x2e,0x76,0xe9,0x4f,0xa8,0x38,0x9d,0x5d,0x3f,0xcb,0xdb,0xad,
+ 0x51,0xa5,0xbf,0xcd,0x72,0xde,0xf7,0xbc,0xcb,0x49,0x2d,0x49,0x26,0xe6,0x1e,
+ 0x9f,0x98,0xe5,0xc6,0x9f,0x2f,0xbb,0x85,0x6b,0x65,0xf6,0x77,0x7c,0x57,0x8b,
+ 0xaa,0xef,0xd8,0x5e,0xa2,0x97,0xe1,0xdc,0x37,0xcd,0x1f,0xe6,0xfc,0xbb,0x8c,
+ 0xb7,0x4e,0xc7,0x3c,0x19,0xd5,0xa8,0x9e,0x67,0x4a,0xe3,0xf5,0x97,0x3a,0x7e,
+ 0x70,0x53,0xfd,0xd6,0xe5,0xb8,0x1c,0x6b,0xee,0xb1,0x9b,0xd1,0xeb,0x34,0xc2,
+ 0x23,0xeb,0x3a,0xf9,0xef,0x16,0xd6,0x4e,0x7d,0x16,0xcf,0xb8,0x1c,0xcb,0xe6,
+ 0x3c,0xda,0xf5,0xcf};
+#define NPRIME_SIEVE30 (sizeof(prime_sieve30)/sizeof(prime_sieve30[0]))
/* Return of 2 if n is prime, 0 if not. Do it fast. */
int _XS_is_prime(UV n)
@@ -148,9 +167,13 @@ int _XS_is_prime(UV n)
const unsigned char* sieve;
int isprime;
- if ( n < (NPRIME_IS_SMALL*8))
- return ((prime_is_small[n/8] >> (n%8)) & 1) ? 2 : 0;
-
+ if (n <= 10) {
+ switch (n) {
+ case 2: case 3: case 5: case 7: return 2; break;
+ default: break;
+ }
+ return 0;
+ }
d = n/30;
m = n - d*30;
mtab = masktab30[m]; /* Bitmask in mod30 wheel */
@@ -159,6 +182,9 @@ int _XS_is_prime(UV n)
if (mtab == 0)
return 0;
+ if (d < NPRIME_SIEVE30)
+ return (prime_sieve30[d] & mtab) ? 0 : 2;
+
isprime = (n <= get_prime_cache(0, &sieve))
? 2*((sieve[d] & mtab) == 0)
: -1;
@@ -168,18 +194,12 @@ int _XS_is_prime(UV n)
}
-static const unsigned char prime_next_small[] =
- {2,2,3,5,5,7,7,11,11,11,11,13,13,17,17,17,17,19,19,23,23,23,23,
- 29,29,29,29,29,29,31,31,37,37,37,37,37,37,41,41,41,41,43,43,47,
- 47,47,47,53,53,53,53,53,53,59,59,59,59,59,59,61,61,67,67,67,67,67,67,71};
-#define NPRIME_NEXT_SMALL (sizeof(prime_next_small)/sizeof(prime_next_small[0]))
-
UV next_trial_prime(UV n)
{
UV d,m;
- if (n < NPRIME_NEXT_SMALL)
- return prime_next_small[n];
+ if (n < 7)
+ return (n < 2) ? 2 : (n < 3) ? 3 : (n < 5) ? 5 : 7;
d = n/30;
m = n - d*30;
@@ -198,8 +218,20 @@ UV _XS_next_prime(UV n)
const unsigned char* sieve;
UV sieve_size;
- if (n < NPRIME_NEXT_SMALL)
- return prime_next_small[n];
+ if (n <= 10) {
+ switch (n) {
+ case 0: case 1: return 2; break;
+ case 2: return 3; break;
+ case 3: case 4: return 5; break;
+ case 5: case 6: return 7; break;
+ default: return 11; break;
+ }
+ }
+ if (n < 30*NPRIME_SIEVE30) {
+ START_DO_FOR_EACH_SIEVE_PRIME(prime_sieve30, n+1, 30*NPRIME_SIEVE30)
+ return p;
+ END_DO_FOR_EACH_SIEVE_PRIME;
+ }
/* Overflow */
#if BITS_PER_WORD == 32
@@ -236,13 +268,20 @@ UV _XS_prev_prime(UV n)
const unsigned char* sieve;
UV sieve_size;
- /* TODO: small prev prime */
if (n <= 7)
return (n <= 2) ? 0 : (n <= 3) ? 2 : (n <= 5) ? 3 : 5;
d = n/30;
m = n - d*30;
+ if (n < 30*NPRIME_SIEVE30) {
+ do {
+ m = prevwheel30[m];
+ if (m==29) { MPUassert(d>0, "d 0 in prev_prime"); d--; }
+ } while (prime_sieve30[d] & masktab30[m]);
+ return(d*30+m);
+ }
+
sieve_size = get_prime_cache(0, &sieve);
if (n < sieve_size) {
do {
@@ -406,12 +445,6 @@ static UV count_segment_ranged(const unsigned char* sieve, UV nbytes, UV lowp, U
*
*/
-static const unsigned char prime_count_small[] =
- {0,0,1,2,2,3,3,4,4,4,4,5,5,6,6,6,6,7,7,8,8,8,8,9,9,9,9,9,9,10,10,
- 11,11,11,11,11,11,12,12,12,12,13,13,14,14,14,14,15,15,15,15,15,15,
- 16,16,16,16,16,16,17,17,18,18,18,18,18,18,19};
-#define NPRIME_COUNT_SMALL (sizeof(prime_count_small)/sizeof(prime_count_small[0]))
-
#define USE_PC_TABLES 1
#if USE_PC_TABLES
/* These tables let us have fast answers up to 3000M for the cost of ~1.4k of
@@ -472,9 +505,6 @@ UV _XS_prime_count(UV low, UV high)
UV segment_size, low_d, high_d;
UV count = 0;
- if ( (low <= 2) && (high < NPRIME_COUNT_SMALL) )
- return prime_count_small[high];
-
if ((low <= 2) && (high >= 2)) count++;
if ((low <= 3) && (high >= 3)) count++;
if ((low <= 5) && (high >= 5)) count++;
@@ -482,6 +512,11 @@ UV _XS_prime_count(UV low, UV high)
if (low > high) return count;
+ if (low == 7 && high <= 30*NPRIME_SIEVE30) {
+ count += count_segment_ranged(prime_sieve30, NPRIME_SIEVE30, low, high);
+ return count;
+ }
+
#if USE_PC_TABLES
if (low == 7 && high >= 30000) {
UV i = 0;
--
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