/* prime.c
John Simpson 2003-02-02
quick-n-dirty prime number finder - uses 32-bit integers
does not use any arrays
does not use any floating point operations
only checks each number for divisibility by previously-known primes
this idea by Steve Litt .
2003-02-04 - adding an optimization which pre-loads the first two
primes (2 and 3) and then uses "add 2, add 4, add 2, add 4, ..."
instead of just incrementing the number we're checking.
2005-04-09 jms1 - changed copyright notice to indicate my intention
that this code is licensed under GPL version 2 only, rather than
GPL version 2 "or later".
****************************************
Copyright (C) 2004-2005 John Simpson.
This program is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License, version 2, as
published by the Free Software Foundation.
This program is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with this program; if not, write to the Free Software
Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
or visit http://www.gnu.org/licenses/gpl.txt
*/
#include
#include
#include
#define NUM unsigned long
#define FORMAT "%lu"
NUM *plist ;
NUM *pend ;
int pcount ;
int main ( int argc , char *argv[] )
{
int fcount , prime ;
NUM n , r , s , t ;
NUM *q ;
unsigned short j ;
struct timeval ts , te ;
struct timezone tz ;
if ( argc < 2 )
{
puts ( "Count not specified" ) ;
return 1 ;
}
fcount = atoi ( argv[1] ) ;
if ( fcount < 1 )
{
puts ( "Count must be 1 or higher" ) ;
return 1 ;
}
/* grab enough memory to hold all the primes we're gonna find */
plist = (NUM *) malloc ( fcount * sizeof ( NUM ) ) ;
if ( ! plist )
{
perror ( "malloc" ) ;
return 1 ;
}
/* preload primes */
pend = plist ;
pcount = 0 ;
*pend++ = 2UL ;
pcount ++ ;
if ( fcount > 1 )
{
*pend++ = 3UL ;
pcount ++ ;
}
/* start the timer */
gettimeofday ( &ts , &tz ) ;
/* let's do it
pcount how many primes we've found
plist pointer to beginning of that list
pend pointer to end of that list
n the number we're testing for prime-ness
r integer "following" n
s the square of r (makes comparisons faster)
t already-known prime that we're testing against
j how far to jump to find the next number to test
cycles between 2 and 4, which skips multiples
of two and multiples of three.
*/
r = s = 1 ;
j = 4 ;
for ( n = 5 ; pcount < fcount ; n += ( j ^= 6 ) )
{
/*
if n >= s, we know that sqrt(n) >= r
need to increment r until s > n
*/
while ( n >= s )
{
r ++ ;
s = r * r ;
}
/* check number against all previously-found primes */
prime = 1 ;
for ( q = plist ; prime && ( q < pend ) ; q ++ )
{
/* get the previously-known prime */
t = *q ;
/* if t > r, we know t > sqrt(n), therefore n is prime */
if ( t > r )
{
break ;
}
/* if n is evenly divisible by t, n is not prime */
if ( ! ( n % t ) )
{
prime = 0 ;
break ;
}
}
/* search is done. is n prime? */
if ( prime )
{
/* add it to the end of the list */
*pend++ = n ;
/* increment the counter as well */
pcount ++ ;
/* maybe tell the user we found one */
/* printf ( "(%d) " FORMAT "\n" , pcount , n ) ; */
}
}
/* finish the timer */
gettimeofday ( &te , &tz ) ;
/* show the user what we found */
/*
for ( q = plist ; q < pend ; q++ )
{
printf ( FORMAT "\n" , *q ) ;
}
*/
/* or maybe show them just the last prime */
printf ( "Prime #%d is " FORMAT "\n" , pcount , plist[pcount-1] ) ;
/* and show them how long it took */
unsigned long long stim , etim , rtim ;
stim = ts.tv_sec * 1000000 + ts.tv_usec ;
etim = te.tv_sec * 1000000 + te.tv_usec ;
rtim = etim - stim ;
printf ( "Elapsed: %llu.%06llu sec\n" ,
( rtim / 1000000 ) , ( rtim % 1000000 ) ) ;
return 0 ;
}