* Copyright (c) 1989 The Regents of the University of California.
* This code is derived from software contributed to Berkeley by
* Redistribution and use in source and binary forms are permitted
* provided that the above copyright notice and this paragraph are
* duplicated in all such forms and that any documentation,
* advertising materials, and other materials related to such
* distribution and use acknowledge that the software was developed
* by the University of California, Berkeley. The name of the
* University may not be used to endorse or promote products derived
* from this software without specific prior written permission.
* THIS SOFTWARE IS PROVIDED ``AS IS'' AND WITHOUT ANY EXPRESS OR
* IMPLIED WARRANTIES, INCLUDING, WITHOUT LIMITATION, THE IMPLIED
* WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
"@(#) Copyright (c) 1989 The Regents of the University of California.\n\
static char sccsid
[] = "@(#)primes.c 5.2 (Berkeley) %G%";
* primes - generate a table of primes between two values
* By: Landon Curt Noll chongo@toad.com, ...!{sun,tolsoft}!hoptoad!chongo
* chongo <for a good prime call: 391581 * 2^216193 - 1> /\oo/\
* Print primes >= start and <= stop. If stop is omitted,
* the value 4294967295 (2^32-1) is assumed. If start is
* omitted, start is read from standard input.
* Prints "ouch" if start or stop is > 4294967295 (2^32-1)
* or if a non-numeric is read on stdin.
* Validation check: there are 664579 primes between 0 and 10^7
* Eratosthenes sieve table
* We only sieve the odd numbers. The base of our sieve windows are always
* odd. If the base of table is 1, table[i] represents 2*i-1. After the
* sieve, table[i] == 1 if and only iff 2*i-1 is prime.
* We make TABSIZE large to reduce the overhead of inner loop setup.
char table
[TABSIZE
]; /* Eratosthenes sieve of odd numbers */
* prime[i] is the (i-1)th prime.
* We are able to sieve 2^32-1 because this byte table yields all primes
* up to 65537 and 65537^2 > 2^32-1.
extern ubig
*pr_limit
; /* largest prime in the prime array */
* To avoid excessive sieves for small factors, we use the table below to
* setup our sieve blocks. Each element represents a odd number starting
* with 1. All non-zero elements are factors of 3, 5, 7, 11 and 13.
extern int pattern_size
; /* length of pattern array */
#define MAX_LINE 255 /* max line allowed on stdin */
char *read_num_buf(); /* read a number buffer */
void primes(); /* print the primes in range */
char *program
; /* our name */
int argc
; /* arg count */
char buf
[MAX_LINE
+1]; /* input buffer */
char *ret
; /* return result */
ubig start
; /* where to start generating */
ubig stop
; /* don't generate at or above this value */
/* convert low and high args */
if (read_num_buf(NULL
, argv
[1]) == NULL
) {
fprintf(stderr
, "%s: ouch\n", program
);
if (read_num_buf(NULL
, argv
[2]) == NULL
) {
fprintf(stderr
, "%s: ouch\n", program
);
if (sscanf(argv
[1], "%ld", &start
) != 1) {
fprintf(stderr
, "%s: ouch\n", program
);
if (sscanf(argv
[2], "%ld", &stop
) != 1) {
fprintf(stderr
, "%s: ouch\n", program
);
if (read_num_buf(NULL
, argv
[1]) == NULL
) {
fprintf(stderr
, "%s: ouch\n", program
);
if (sscanf(argv
[1], "%ld", &start
) != 1) {
fprintf(stderr
, "%s: ouch\n", program
);
/* read input until we get a good line */
if (read_num_buf(stdin
, buf
) != NULL
) {
if (sscanf(buf
, "%ld", &start
) != 1) {
fprintf(stderr
, "%s: ouch\n", program
);
fprintf(stderr
, "%s: ouch\n", program
);
* read_num_buf - read a number buffer from a stream
* Read a number on a line of the form:
* ^[ \t]*\(+?[0-9][0-9]\)*.*$
* where ? is a 1-or-0 operator and the number is within \( \).
* If does not match the above pattern, it is ignored and a new
* line is read. If the number is too large or small, we will
* print ouch and read a new line.
* We have to be very careful on how we check the magnitude of the
* input. We can not use numeric checks because of the need to
* check values against maximum numeric values.
* This routine will return a line containing a ascii number between
* 0 and BIG, or it will return NULL.
* If the stream is NULL then buf will be processed as if were
* char * pointer to leading digit or +
FILE *input
; /* input stream or NULL */
char *buf
; /* input buffer */
static char limit
[MAX_LINE
+1]; /* ascii value of BIG */
static int limit_len
; /* digit count of limit */
int len
; /* digits in input (excluding +/-) */
char *s
; /* line start marker */
char *d
; /* first digit, skip +/- */
char *p
; /* scan pointer */
char *z
; /* zero scan pointer */
/* form the ascii value of SEMIBIG if needed */
if (!isascii(limit
[0]) || !isdigit(limit
[0])) {
sprintf(limit
, "%ld", SEMIBIG
);
limit_len
= strlen(limit
);
* the search for a good line
if (input
!= NULL
&& fgets(buf
, MAX_LINE
, input
) == NULL
) {
/* ignore leading whitespace */
for (s
=buf
; *s
&& s
< buf
+MAX_LINE
; ++s
) {
if (!isascii(*s
) || !isspace(*s
)) {
fprintf(stderr
, "%s: ouch\n", program
);
/* skip over any leading + */
for (z
=d
; *z
&& z
< buf
+MAX_LINE
; ++z
) {
/* scan for the first non-digit/non-plus/non-minus */
for (p
=d
; *p
&& p
< buf
+MAX_LINE
; ++p
) {
if (!isascii(*p
) || !isdigit(*p
)) {
/* object if too many digits */
len
= (len
<=0) ? 1 : len
;
/* accept if digit count is below limit */
/* reject very large numbers */
} else if (len
> limit_len
) {
fprintf(stderr
, "%s: ouch\n", program
);
/* carefully check against near limit numbers */
} else if (strcmp(z
, limit
) > 0) {
fprintf(stderr
, "%s: ouch\n", program
);
/* number is near limit, but is under it */
} while (input
!= NULL
&& fgets(buf
, MAX_LINE
, input
) != NULL
);
* primes - sieve and print primes from start up to and but not including stop
ubig start
; /* where to start generating */
ubig stop
; /* don't generate at or above this value */
register char *q
; /* sieve spot */
register ubig factor
; /* index and factor */
register char *tab_lim
; /* the limit to sieve on the table */
register ubig
*p
; /* prime table pointer */
register ubig fact_lim
; /* highest prime for current block */
* A number of systems can not convert double values
* into unsigned longs when the values are larger than
* the largest signed value. Thus we take case when
* the double is larger than the value SEMIBIG. *sigh*
* be sure that the values are odd, or 2
if (start
!= 2 && (start
&0x1) == 0) {
if (stop
!= 2 && (stop
&0x1) == 0) {
* quick list of primes <= pr_limit
if (start
<= *pr_limit
) {
/* skip primes up to the start value */
for (p
= &prime
[0], factor
= prime
[0];
factor
< stop
&& p
<= pr_limit
;
/* return early if we are done */
* we shall sieve a bytemap window, note primes and move the window
* upward until we pass the stop point
* factor out 3, 5, 7, 11 and 13
/* initial pattern copy */
factor
= (start
%(2*3*5*7*11*13))/2; /* starting copy spot */
memcpy(table
, &pattern
[factor
], pattern_size
-factor
);
/* main block pattern copies */
for (fact_lim
=pattern_size
-factor
;
fact_lim
+pattern_size
<=TABSIZE
;
fact_lim
+=pattern_size
) {
memcpy(&table
[fact_lim
], pattern
, pattern_size
);
/* final block pattern copy */
memcpy(&table
[fact_lim
], pattern
, TABSIZE
-fact_lim
);
* sieve for primes 17 and higher
/* note highest useful factor and sieve spot */
if (stop
-start
> TABSIZE
+TABSIZE
) {
tab_lim
= &table
[TABSIZE
]; /* sieve it all */
(double)(start
)+TABSIZE
+TABSIZE
+1.0);
tab_lim
= &table
[(stop
-start
)/2]; /* partial sieve */
fact_lim
= (int)sqrt((double)(stop
)+1.0);
/* sieve for factors >= 17 */
factor
= 17; /* 17 is first prime to use */
p
= &prime
[7]; /* 19 is next prime, pi(19)=7 */
/* determine the factor's initial sieve point */
q
= (char *)(start
%factor
); /* temp storage for mod */
q
= &table
[(factor
-(int)q
)/2];
q
= &table
[q
? factor
-((int)q
/2) : 0];
/* sive for our current factor */
for ( ; q
< tab_lim
; q
+= factor
) {
*q
= '\0'; /* sieve out a spot */
} while ((factor
=(ubig
)(*(p
++))) <= fact_lim
);
for (q
= table
; q
< tab_lim
; ++q
, start
+=2) {