00001 /* Compute and print the prime numbers 00002 Copyright (C) 2001, 2002, 2003 Free Software Foundation, Inc. 00003 Written by Stephane Carrez (stcarrez@nerim.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 00054 00055 #include <stdio.h> 00056 #include <stdarg.h> 00057 #include <sys/param.h> 00058 #include <sys/sio.h> 00059 #include <string.h> 00060 #include <imath.h> 00061 00062 /* The DATA_SIZE depends on your board. The more memory, the largest 00063 the bitmap can be. The last number we can record is 8 times the 00064 size of that bitmap (1 bit for each number). */ 00065 00066 /* Remove the following test that reduces the bitmap. 00067 With a 32K memory system, the bitmap can contain up to 00068 224768 bits and the computation takes... a lot of time... */ 00069 #if DATA_SIZE > 8192 + 512 00070 # undef DATA_SIZE 00071 # define DATA_SIZE 8192 + 512 00072 #endif 00073 00074 /* Leave 512 or 128 bytes for stack. */ 00075 #if DATA_SIZE < 512 00076 # define MAX_PRIMES (DATA_SIZE - 128) 00077 #else 00078 # define MAX_PRIMES (DATA_SIZE - 512) 00079 #endif 00080 00081 #define LAST_PRIME (MAX_PRIMES * 8L) 00082 00083 /* Bitmap representing the prime numbers. */ 00084 unsigned char prime_list[MAX_PRIMES]; 00085 00086 /* Returns true if 'n' is a prime number recorded in the table. */ 00087 static inline int 00088 is_prime (unsigned long n) 00089 { 00090 unsigned short bit = (unsigned short) (n) & 0x07; 00091 00092 return prime_list[n >> 3] & (1 << bit); 00093 } 00094 00095 /* Record 'n' as a prime number in the table. */ 00096 static inline void 00097 set_prime (unsigned long n) 00098 { 00099 unsigned short bit = (unsigned short) (n) & 0x07; 00100 00101 prime_list[n >> 3] |= (1 << bit); 00102 } 00103 00104 /* Check whether 'n' is a prime number. 00105 Returns 1 if it's a prime number, 0 otherwise. */ 00106 static int 00107 check_for_prime (unsigned long n) 00108 { 00109 unsigned long i; 00110 unsigned char *p; 00111 unsigned long last_value; 00112 unsigned char small_n; 00113 00114 small_n = (n & 0xffff0000) == 0; 00115 i = 0; 00116 00117 /* We can stop when we have checked all prime numbers below sqrt(n). */ 00118 last_value = lsqrt (n); 00119 00120 /* Scan the bitmap of prime numbers and divide 'n' by the corresponding 00121 prime to see if it's a multiple of it. */ 00122 p = prime_list; 00123 do 00124 { 00125 unsigned char val; 00126 00127 val = *p++; 00128 if (val) 00129 { 00130 unsigned short j; 00131 unsigned long q; 00132 00133 q = i; 00134 for (j = 1; val && j <= 0x80; j <<= 1, q++) 00135 { 00136 if (val & j) 00137 { 00138 val &= ~j; 00139 00140 /* Use 16-bit division if 'n' is small enough. */ 00141 if (small_n) 00142 { 00143 unsigned short r; 00144 00145 /* 'n' is a multiple of prime 'q'. */ 00146 r = (unsigned short) (n) % (unsigned short) (q); 00147 if (r == 0) 00148 return 0; 00149 } 00150 else 00151 { 00152 unsigned long r; 00153 00154 r = n % q; 00155 00156 /* 'n' is a multiple of prime 'q'. */ 00157 if (r == 0) 00158 return 0; 00159 } 00160 } 00161 } 00162 } 00163 i += 8; 00164 } 00165 while (i < last_value); 00166 return 1; 00167 } 00168 00169 #if 0 00170 /* Utility function that can be called from gdb to dump the prime 00171 numbers. Do: `call print_primes()' from gdb. */ 00172 static void 00173 print_primes (void) 00174 { 00175 long i; 00176 00177 for (i = 0; i < LAST_PRIME; i++) 00178 if (is_prime (i)) 00179 printf ("%ld\n", i); 00180 } 00181 #endif 00182 00183 /* Set this variable to 1 if you want to print the prime numbers 00184 while they are found. */ 00185 int verbose = 0; 00186 00187 int 00188 main () 00189 { 00190 long i; 00191 short cnt = 0; 00192 00193 serial_init (); 00194 printf ("Computing prime numbers below %ld\n", 00195 (long) LAST_PRIME); 00196 memset (prime_list, 0, sizeof (prime_list)); 00197 for (i = 2; i < LAST_PRIME; i++) 00198 { 00199 if (check_for_prime (i)) 00200 { 00201 set_prime (i); 00202 cnt ++; 00203 if (verbose) 00204 printf ("%ld\n", i); 00205 } 00206 } 00207 printf ("Found %ld prime numbers below %ld\n", 00208 (long) cnt, (long) LAST_PRIME); 00209 00210 return 0; 00211 }