primes.c

Description | Download | Table of Contents | Modules | Compound List | File List | Functions


Overview
Compiler
Documentation
Examples
Misc
Help
IDE & Tools

Download
Install

Links
Projects






00001 /* Compute and print the prime numbers
00002    Copyright (C) 2001 Free Software Foundation, Inc.
00003    Written by Stephane Carrez (stcarrez@worldnet.fr)    
00004 
00005 This file is free software; you can redistribute it and/or modify it
00006 under the terms of the GNU General Public License as published by the
00007 Free Software Foundation; either version 2, or (at your option) any
00008 later version.
00009 
00010 In addition to the permissions in the GNU General Public License, the
00011 Free Software Foundation gives you unlimited permission to link the
00012 compiled version of this file with other programs, and to distribute
00013 those programs without any restriction coming from the use of this
00014 file.  (The General Public License restrictions do apply in other
00015 respects; for example, they cover modification of the file, and
00016 distribution when not linked into another program.)
00017 
00018 This file is distributed in the hope that it will be useful, but
00019 WITHOUT ANY WARRANTY; without even the implied warranty of
00020 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00021 General Public License for more details.
00022 
00023 You should have received a copy of the GNU General Public License
00024 along with this program; see the file COPYING.  If not, write to
00025 the Free Software Foundation, 59 Temple Place - Suite 330,
00026 Boston, MA 02111-1307, USA.  */
00027 
00053 
00054 #include <stdio.h>
00055 #include <stdarg.h>
00056 #include <sys/param.h>
00057 #include <imath.h>
00058 
00059 /* The DATA_SIZE depends on your board.  The more memory, the largest
00060    the bitmap can be.  The last number we can record is 8 times the
00061    size of that bitmap (1 bit for each number).  */
00062 
00063 /* Remove the following test that reduces the bitmap.
00064    With a 32K memory system, the bitmap can contain up to
00065    224768 bits and the computation takes... a lot of time... */
00066 #if DATA_SIZE > 8192 + 512
00067 # undef DATA_SIZE
00068 # define DATA_SIZE 8192 + 512
00069 #endif
00070 
00071 /* Leave 512 or 128 bytes for stack.  */
00072 #if DATA_SIZE < 512
00073 # define MAX_PRIMES (DATA_SIZE - 128)
00074 #else
00075 # define MAX_PRIMES (DATA_SIZE - 512)
00076 #endif
00077 
00078 #define LAST_PRIME (MAX_PRIMES * 8L)
00079 
00080 /* Bitmap representing the prime numbers.  */
00081 unsigned char prime_list[MAX_PRIMES];
00082 
00083 /* Returns true if 'n' is a prime number recorded in the table.  */
00084 static inline int
00085 is_prime (unsigned long n)
00086 {
00087   unsigned short bit = (unsigned short) (n) & 0x07;
00088   
00089   return prime_list[n >> 3] & (1 << bit);
00090 }
00091 
00092 /* Record 'n' as a prime number in the table.  */
00093 static inline void
00094 set_prime (unsigned long n)
00095 {
00096   unsigned short bit = (unsigned short) (n) & 0x07;
00097 
00098   prime_list[n >> 3] |= (1 << bit);
00099 }
00100 
00101 /* Check whether 'n' is a prime number.
00102    Returns 1 if it's a prime number, 0 otherwise.  */
00103 int
00104 check_for_prime (unsigned long n)
00105 {
00106   unsigned long i;
00107   unsigned char *p;
00108   unsigned long last_value;
00109   unsigned char small_n;
00110 
00111   small_n = (n & 0xffff0000) == 0;
00112   i = 0;
00113 
00114   /* We can stop when we have checked all prime numbers below sqrt(n).  */
00115   last_value = lsqrt (n);
00116 
00117   /* Scan the bitmap of prime numbers and divide 'n' by the corresponding
00118      prime to see if it's a multiple of it.  */
00119   p = prime_list;
00120   do
00121     {
00122       unsigned char val;
00123       
00124       val = *p++;
00125       if (val)
00126         {
00127           unsigned short j;
00128           unsigned long q;
00129 
00130           q = i;
00131           for (j = 1; val && j <= 0x80; j <<= 1, q++)
00132             {
00133               if (val & j)
00134                 {
00135                   val &= ~j;
00136 
00137                   /* Use 16-bit division if 'n' is small enough.  */
00138                   if (small_n)
00139                     {
00140                       unsigned short r;
00141 
00142                       /* 'n' is a multiple of prime 'q'.  */
00143                       r = (unsigned short) (n) % (unsigned short) (q);
00144                       if (r == 0)
00145                         return 0;
00146                     }
00147                   else
00148                     {
00149                       unsigned long r;
00150                       
00151                       r = n % q;
00152 
00153                       /* 'n' is a multiple of prime 'q'.  */
00154                       if (r == 0)
00155                         return 0;
00156                     }
00157                 }
00158             }
00159         }
00160       i += 8;
00161     }
00162   while (i < last_value);
00163   return 1;
00164 }
00165 
00166 /* Utility function that can be called from gdb to dump the prime
00167    numbers.  Do:  `call print_primes()' from gdb.  */
00168 void
00169 print_primes ()
00170 {
00171   long i;
00172   
00173   for (i = 0; i < LAST_PRIME; i++)
00174     if (is_prime (i))
00175       printf ("%ld\n", i);
00176 }
00177 
00178 /* Set this variable to 1 if you want to print the prime numbers
00179    while they are found.  */
00180 int verbose = 0;
00181 
00182 int
00183 main ()
00184 {
00185   long i;
00186   short cnt = 0;
00187   
00188   printf ("Computing prime numbers below %ld\n",
00189           (long) LAST_PRIME);
00190   memset (prime_list, 0, sizeof (prime_list));
00191   for (i = 2; i < LAST_PRIME; i++)
00192     {
00193       if (check_for_prime (i))
00194         {
00195           set_prime (i);
00196           cnt ++;
00197           if (verbose)
00198             printf ("%ld\n", i);
00199         }
00200     }
00201   printf ("Found %ld prime numbers below %ld\n",
00202           (long) cnt, (long) LAST_PRIME);
00203 
00204   return 0;
00205 }

Description | Download | Table of Contents | Modules | Compound List | File List | Functions

    Last modified,
    Apr 16, 2001
[ Copying ]     [ Feedback to Stephane.Carrez@worldnet.fr ]