|
2007-10-13, 17:12
|
#12 | |
|
A Sunny Moo
Aug 2007
USA (GMT-5)
3·2,083 Posts |
Quote:
|
|
|
|
|
2007-10-17, 17:44
|
#13 | |
Dec 2005
22·23 Posts |
Quote:
I posted my program on channel9.msdn.com in a thread I started to find out why it did not compile in Visual Studio (as I originally developed it using GCC on Unix). I subsequently posted later, more enhanced versions of it for additional criticism, as I wrote it to become more proficient in C and I wanted tips from more experienced programmers on how it could have been written to execute more quickly. I am at the university right now, so I do not have access to my copy of the program at home, but here is the most recent version I posted online: Code:
/* prime.c */
/* Sept 29, 2007. */
/* <my name was here> */
/* Program to find prime numbers */
#include <stdlib.h>
int main (void)
{
unsigned int len, b_num, num, bitVal, *arr = NULL, *arrEnd = NULL;
printf ("\nFind prime numbers smaller than \n");
scanf("%i", &len);
if (len < 1)
{
return 0;
}
b_num = ( len % (16 * sizeof(unsigned int)) == 0 ) ?
(len / (16 * sizeof(unsigned int))) :
(len / (16 * sizeof(unsigned int)) + 1);
arr = (unsigned int*)calloc (b_num, sizeof(unsigned int));
arrEnd = arr + b_num;
bin_arr_prime_fill(arr, b_num);
if (len >= 2)
{
printf("Number 2 is a prime number\n");
}
/* Loop through Prime Array and Output Primes */
for ( num = 1; arr < arrEnd ; ++arr )
{
for ( bitVal = 0x80000000; bitVal > 0 && num <= len; bitVal >>= 1, num += 2 )
{
if ((*arr & bitVal) == bitVal)
{
printf("Number %i is a prime number\n", num);
}
}
}
free(arr - b_num);
main();
return 0;
}
int bin_arr_prime_fill (const unsigned int * const arr, const unsigned int len)
{
unsigned int b, num, *a = arr, *aEnd = arr + len;
/* b has a binary 1 where the number is in the array */
/* num is the number being checked for primiality */
/* a is used to seek to the appropriate position in the array * /
/* aEnd is used to determine when to stop */
for ( num = 3, b = 0x40000000; a < aEnd; ++a, b = 0x80000000 )
{
for ( ; b > 0 ; b >>= 1, num += 2 )
{
/* arr is the prime array */
/* a points to the block currently being built */
/* num is the number being tested */
if (isPrime(arr, num) == 1)
{
*a |= b;
}
}
}
return 0;
}
int isPrime(const unsigned int * arr, const unsigned int num)
{
unsigned int bitVal = 0x40000000, div = 3;
/* arr points to an array; it is used to seek to the correct block */
/* num is the number being checked for primiality */
/* bitVal is the bit position */
/* div is the divisor, starting with 3 */
/* Set bitVal = 0x2000 so it corresponds to m = 3 for the */
/* first inner loop */
/* Set bitVal = 0x8000 so it corresponds to the first number */
/* in the block in every suceeding loop */
/* The outer loop loops through each block of the array */
for ( ; div <= (num / div); ++arr, bitVal = 0x80000000)
{
/* The inner loop checks each number in the block for\ */
/* primality, skipping even numbers */
for ( ; div <= (num / div) ; bitVal >>= 1, div += 2 )
{
/* If the divisor is prime and it divides */
/* n evenly, return composite */
if ((*arr & bitVal) == bitVal && (num % div) == 0)
{
return 0;
}
}
}
/* Return Prime */
return 1;
}
If my memory serves me correctly, it took 12 to 16 hours to do such a computation on my computer, although I did not take precise time measurements as I only wanted a ball park estimate of how long it took to compute. The only potential single threaded method I can think of doing to compute this more quickly would be to take advantage of differences in the computational speed of a processor's ALUs and FPUs by using floats instead of unsigned integers, but I have not tried it as I am working on another project. I am curious, how did you compute all of the primes between 1 and 6 billion in 24 hours on a 400MHz Pentium 3? |
|
|
|
|
|
2007-10-17, 19:26
|
#14 | |
Jun 2003
2·3·7·112 Posts |
Quote:
Please take a look at here for an explanation of SoE. A properly implemented SoE can do the range (2^31) in a few seconds (not hours). |
|
|
|
|
2007-10-17, 19:34
|
#15 | |
|
Dec 2005
22×23 Posts |
Quote:
|
|
|
|
|
|
2007-10-17, 21:15
|
#16 |
|
Dec 2005
22·23 Posts |
I rewrote my program. It required far less code to start with a bitmap of ones and then mark off the zeroes. It also can sieve from 1 to 1073741824 in next to no time at all.
Here is the code: Code:
/* prime.c */
/* Oct 17, 2007. */
/* <my name was here> */
/* Program to find prime numbers */
#include <stdlib.h>
int main (void)
{
unsigned int len, b_num, num, bitVal, *arr = NULL;
printf ("\nFind prime numbers smaller than \n");
scanf("%i", &len);
if (len < 1)
{
return 0;
}
b_num = ( len % (16 * sizeof(unsigned int)) == 0 ) ?
(len / (16 * sizeof(unsigned int))) :
(len / (16 * sizeof(unsigned int)) + 1);
arr = (unsigned int*)calloc (b_num, sizeof(unsigned int));
arr[0] = 0x7FFFFFFF;
for ( num = 1 ; num < b_num ; ++num )
{
arr[num] = 0xFFFFFFFF;
}
bin_arr_prime_fill(arr, b_num);
if (len >= 2)
{
printf("Number 2 is a prime number\n");
}
/* Loop through Prime Array and Output Primes */
for ( num = 1; num <= len ; ++arr )
{
for ( bitVal = 0x80000000; bitVal > 0 && num <= len; bitVal >>= 1,
num += 2 )
{
if ((*arr & bitVal) == bitVal)
{
printf("Number %i is a prime number\n", num);
}
}
}
free(arr - b_num);
main();
return 0;
}
int bin_arr_prime_fill (unsigned int * const arr, const unsigned int b_num)
{
unsigned int num, i, len = b_num * 64;
/* num is the number being checked for primiality */
for ( num = 3; num * num <= len ; num += 2 )
{
if ((arr[num >> 6] & (0x80000000 >> ((num >> 1) % 32))) == 0)
{
continue;
}
for ( i = num * 3; (i >> 6) < b_num ; i += num << 1 )
{
if (arr[i >> 6] & (0x80000000 >> ((i >> 1) % 32)))
{
arr[i >> 6] ^= (0x80000000 >> ((i >> 1) % 32));
}
}
}
return 0;
}
Last fiddled with by ShiningArcanine on 2007-10-17 at 21:16 |
|
|
|
|
2007-10-17, 23:01
|
#17 |
"Mike"
Aug 2002
202A16 Posts |
|
|
|
|
|
2007-10-17, 23:22
|
#18 |
"Ben"
Feb 2007
22×3×293 Posts |
Getting rid of the function call in the innermost loop certainly helped too. That was a fast rework of your code, good job.
|
|
|
|
|
2007-10-18, 01:09
|
#19 | |
|
Dec 2005
22×23 Posts |
Quote:
|
|
|
|
|
|
2007-10-18, 01:34
|
#20 |
|
"Ben"
Feb 2007
22×3×293 Posts |
It looks like you got rid of all the division operators in the inner loops. You can also get rid of the % operators when %'ing by a power of 2. since division by a power of 2 (in this case 32) just shifts out the lowest bits, using & accomplishes the same thing. % 32 becomes & 31.
|
|
|
|
|
2007-10-18, 01:58
|
#21 |
|
"Ben"
Feb 2007
22·3·293 Posts |
I went ahead and tried out the above suggestion on your code because I'm looking for excuses to not do my homework at the moment
 It had negligible impact (Compiling with MSVC 6). To further increase speed (assuming you want to) you'll have to do something about the cache performance of the sieve. Read about cache blocking, or a segmented sieve. This kind of optimization can be addictive. I've spent a fair amount of time tuning mine up, and while it's fast, it's still not the fastest I've seen (but close). There's always something else to try... |
|
|
|
|
2007-10-18, 02:23
|
#22 |
|
Dec 2005
22×23 Posts |
I found a slight bug in my code, which lowered upper bound for the sieve fairly substantially. I had only tested very large upper bounds, which caused me to miss the bug. With the fix, for 1 to 1073741824, my program is not quite as fast as I thought it was earlier, but it is still substantially faster than what it was, taking 10 to 11 seconds to find all of the primes between 1 and 1073741824 on my computer versus the approximate 12 to 16 hours it took before the rewrite.
Thanks for the tip regarding using binary AND instead of the arithmetic modulus. I have made that change as well as a few other enhancements to my code. Among them I added a new function to make the code more readable. I have the /Ob2 flag set in Visual Studio instructing it to inline functions whenever it can, so I assume that it is inlined sometime during the compilation process. Here is the code for the new version: Code:
/* prime.c */
/* Oct 17, 2007. */
/* <my name was here> */
/* Program to find prime numbers */
#include <stdlib.h>
int main (void)
{
unsigned int len, num, bitVal, *arr = NULL;
printf ("\nFind prime numbers smaller than \n");
scanf("%i", &len);
if (len < 1)
{
return 0;
}
arr = (unsigned int*)calloc (blocks(len), sizeof(unsigned int));
bin_arr_prime_fill(arr, len);
if (len >= 2)
{
printf("Number 2 is a prime number\n");
}
/* Loop through Prime Array and Output Primes */
for ( num = 1; num <= len ; ++arr )
{
for ( bitVal = 0x80000000; bitVal > 0 && num <= len; bitVal >>= 1,
num += 2 )
{
if ((*arr & bitVal) == bitVal)
{
printf("Number %i is a prime number\n", num);
}
}
}
free(arr - blocks(len));
main();
return 0;
}
unsigned int blocks (const unsigned int len)
{
return ( len & ((sizeof(unsigned int) << 4) - 1) == 0 ) ?
(len / (sizeof(unsigned int) << 4)) :
(len / (sizeof(unsigned int) << 4) + 1);
}
int bin_arr_prime_fill (unsigned int * const arr, const unsigned int len)
{
unsigned int num, i, b_num = blocks(len);
arr[0] = 0x7FFFFFFF;
for ( num = 1 ; num < b_num ; ++num )
{
arr[num] = 0xFFFFFFFF;
}
/* num is the number being checked for primiality */
for ( num = 3; num * num <= len ; num += 2 )
{
if ((arr[num >> 6] & (0x80000000 >> ((num >> 1) & 0x1F))) == 0)
{
continue;
}
for ( i = num * 3; (i >> 6) < b_num ; i += num << 1 )
{
if (arr[i >> 6] & (0x80000000 >> ((i >> 1) & 0x1F)))
{
arr[i >> 6] ^= (0x80000000 >> ((i >> 1) & 0x1F));
}
}
}
return 0;
}
Last fiddled with by ShiningArcanine on 2007-10-18 at 02:25 |
|
|
|
 |
| Thread Tools | |
Similar Threads
|
||||
| Thread | Thread Starter | Forum | Replies | Last Post |
| Finding VERY large primes | c10ck3r | Information & Answers | 34 | 2012-08-29 16:47 |
| Why arent there many softwares for finding Huge Primes | blistervol | Math | 2 | 2012-08-20 17:26 |
| Best Work for Finding Primes | Unregistered | Information & Answers | 9 | 2012-06-24 13:50 |
| Finding primes using modular stacking | goatboy | Math | 1 | 2007-12-07 12:30 |
| Finding primes with a PowerPC | rogue | Lounge | 4 | 2005-07-12 12:31 |
) I had to stop at about the number 6 Billion. It only took about a day to run on a P3 400Mhz.
