| Subject: | Possible feature: prime count |
Consider adding a prime count function, so one can retrieve counts
without having to get a giant array back. Here is a simple
implementation for XS.xs.
void
xs_count_primes (number, base)
unsigned long number
unsigned long base
PROTOTYPE: $$
INIT:
unsigned char *composite = NULL;
unsigned long i, n;
unsigned long count = 0;
PPCODE:
const unsigned long square_root = sqrt (number); /* truncates */
if (base > number) XSRETURN_UV(0);
Newxz (composite, number/16 + 1, unsigned char);
for (n = 3; n <= square_root; n += 2) {
if ((composite[n/16] & (1 << ((n/2)%8))) == 0)
for (i = n*n; i <= number; i += 2*n)
composite[i/16] |= 1 << ((i/2)%8);
}
if (base <= 2) {
if (number >= 2)
count++;
base = 3;
}
base |= 1; /* make odd */
for (n = base; n <= number; n += 2)
if ((composite[n/16] & (1 << ((n/2)%8))) == 0)
count++;
Safefree (composite);
XSRETURN_UV(count);
There are faster solutions and this doesn't segment so isn't fast for
things like 10^14 to 10^14+1000, but those are shared by the existing
sieve code. For large ranges this is faster and much less memory than
getting an array back. It's about 4x faster than using the sieve from
v0.26.
The usual changes would need to additionally be made to
lib/Math/Prime/XS.pm for exporting and documentation.