Go Down

Topic: Arduino generating prime numbers (Read 10 times) previous topic - next topic

Nick Gammon

Oct 08, 2013, 01:52 pm Last Edit: Nov 25, 2013, 10:54 pm by Nick Gammon Reason: 1
I got inspired after watching the guys on Numberphile unbox a gadget that displays prime numbers:

http://www.youtube.com/watch?v=8UqCyepX3AI

I made a quick version that did a similar thing with an 8-digit LED display:

http://vimeo.com/80314148

Totally useless of course, but strangely mesmerizing to watch.

Hardware: Arduino Uno, plus one of these 8-digit display boards from eBay for about $10.





References:

http://en.wikipedia.org/wiki/Primality_test

http://primes.utm.edu/

Code:

Code: [Select]

// Prime number generator
//
// Author: Nick Gammon
// Date:   8th October 2013

#include <SPI.h>

const int SHOW_EVERY = 500;  // how often to echo a prime to the serial port
const unsigned long TIME_TO_WAIT = 500;  // mS

// MAX7219 registers

const byte MAX7219_REG_NOOP        = 0x0;
// codes 1 to 8 are digit positions 1 to 8
const byte MAX7219_REG_DECODEMODE  = 0x9;
const byte MAX7219_REG_INTENSITY   = 0xA;
const byte MAX7219_REG_SCANLIMIT   = 0xB;
const byte MAX7219_REG_SHUTDOWN    = 0xC;
const byte MAX7219_REG_DISPLAYTEST = 0xF;

// 7-segments patterns for digits 0 to 9
const byte digits [10] = {
 0b1111110,  // 0
 0b0110000,  // 1
 0b1101101,  // 2
 0b1111001,  // 3
 0b0110011,  // 4
 0b1011011,  // 5
 0b1011111,  // 6
 0b1110000,  // 7
 0b1111111,  // 8
 0b1111011,  // 9
};


// send one byte via  SPI to MAX7219
void sendByte (const byte reg, const byte data)
 {    
 digitalWrite (SS, LOW);
 SPI.transfer (reg);
 SPI.transfer (data);
 digitalWrite (SS, HIGH);
 }  // end of sendByte

// send one character (data) to position (pos) with or without decimal place
// pos is 0 to 7
// data can be '0' to '9' plus various other things as below
void sendChar (const byte pos, const char data, const bool dp = false)
 {
 byte converted;  
 switch (data)
   {
   case '0' ... '9' : converted = digits [data - '0']; break;
   case '-':          converted = 0b0000001;  break;  
   case 'A':          converted = 0b1110111;  break;  
   case 'b':          converted = 0b0011111;  break;  
   case 'c':          converted = 0b0001101;  break;  
   case 'C':          converted = 0b1001111;  break;  
   case 'd':          converted = 0b0111101;  break;  
   case 'E':          converted = 0b1001111;  break;  
   case 'F':          converted = 0b1000111;  break;  
   case 'h':          converted = 0b0010111;  break;  
   case 'H':          converted = 0b0110111;  break;  
   case 'L':          converted = 0b0001110;  break;  
   case 'o':          converted = 0b0011101;  break;  
   case 'P':          converted = 0b1100111;  break;  
   case 'r':          converted = 0b0000101;  break;  
   case 'S':          converted = 0b1011011;  break;  
   case 't':          converted = 0b0000111;  break;  
   case 'u':          converted = 0b0011100;  break;  
   case ' ':          converted = 0b0000000;  break;  
   case 'I':          converted = digits [1]; break;
   
   default:           converted = 0b0000001;  break;  // -
   } // end of switch  
 if (dp)
   converted |= 0b10000000;
 sendByte (8 - pos, converted);
 }  // end of sendChar
 
// write an entire null-terminated string to the LEDs
void sendString (const char * s)
{
 byte pos;
 
 for (pos = 0; pos < 8 && *s; pos++)
   {
   boolean dp = s [1] == '.';
   sendChar (pos, *s++, dp);   // turn decimal place on if next char is a dot
   if (dp)  // skip dot
     s++;
   }
   
 // space out rest
 while (pos < 8)
   sendChar (pos++, ' ');
}  // end of sendString


long candidate;
long found = 5; // Number we found
int count = found - 1;


void setup() {
 Serial.begin(115200);
 
 while (!Serial) { }
 
 Serial.println ();
 Serial.println ();
 
 SPI.begin ();

 sendByte (MAX7219_REG_SCANLIMIT, 7);      // show 8 digits
 sendByte (MAX7219_REG_DECODEMODE, 0);     // use bit patterns
 sendByte (MAX7219_REG_DISPLAYTEST, 0);    // no display test
 sendByte (MAX7219_REG_INTENSITY, 7);      // character intensity: range: 0 to 15
 sendString ("");                          // clear display
 sendByte (MAX7219_REG_SHUTDOWN, 1);       // not in shutdown mode (ie. start it up)

}

unsigned long start;
unsigned long lastShown;

void loop()
 {

 for (candidate = 3; candidate < 99999999; candidate += 2)
   showIfPrime (candidate);

 }  // end of loop

void showIfPrime (long num)
 {
 
 // we are already skipping odd numbers, now skip if divisible by 3
 if (num % 3 == 0)
   return;

 // Only have to check for divisible for the sqrt(number)
 long upper = sqrt(num) + 6;
 
 // Check if the number is divisible (start at 6 going up)
 for (long cnum = 6; cnum <= upper; cnum += 6)
   {
   if (num % (cnum - 1) == 0)
     return;
   if (num % (cnum + 1) == 0)
     return;
   }  // end of checking cnum-1, cnum+1

// echo to serial port for validating
if (++count >= SHOW_EVERY)
 {
 Serial.print(candidate);
 Serial.print(" is prime. Took ");
 Serial.print(millis () - start);
 Serial.print (" mS. Found ");
 Serial.print (found);
 Serial.println (" primes.");
 start = millis ();
 count = 0;
 }

 found++;

 // delay until interval is up
 // (this absorbs the calculation time)
 while (millis () - lastShown < SHOW_EVERY)
   { }
   
 char buf [10];
 sprintf (buf, "%8ld", candidate);
 sendString (buf);
 lastShown = millis ();
 
}  // end of showIfPrime


The code was written in a bit of a hurry. It doesn't show the first few primes (well, we know what they are). I think it starts at 7, because of the nature of the algorithm. When I timed it at around the 128201 mark (the prime number 128201) it was taking about 10 mS per prime. Of course, this makes the display tick over unreadably fast, so it is throttled back to show a number every 500 mS.

The algorithm skips even numbers (naturally) and also numbers divisible by 3. Then it tests in jumps of 6, plus or minus 1, to see if they are primes, up to the square root of the number we are testing for. It isn't particularly fast, but also doesn't use a lot of memory.

Please post technical questions on the forum, not by personal message. Thanks!

More info:
http://www.gammon.com.au/electronics

fungus

The next step is to start launching them into space towards candidate planets taken from the Kepler telescope's web site.


No, I don't answer questions sent in private messages (but I do accept thank-you notes...)

robtillaart

Hi Nick,
Nice prime machine, as every other programmer did some prime-finding too, see the 30k +i optimization here
- http://forum.arduino.cc/index.php?topic=63071.msg457935#msg457935 -

The nice part of this  trick (of the 30K method) is that you can encode the set of eight possible primes into one byte. imagine an array in which the arrayindex indicates the numberx30 and every bit indicates one of the eight mentioned offsets. a 1 bit means prime, a 0 bit means no prime.
Rob Tillaart

Nederlandse sectie - http://arduino.cc/forum/index.php/board,77.0.html -
(Please do not PM for private consultancy)

Nick Gammon

Just as an update, it's been running all night (I'm guessing around 15 hours) and it is up to 1442723. So still the last digit to go.

Once again though, it is throttled back to update every 500 mS, and just by eye it looks like it is able to compute the next number in that time.

If you want to try it at home (without the LED display) there is also serial output. For debugging it prints every 500th prime (you can change that near the start). And you can change the delay between generating primes.

Program size:

Code: [Select]

Binary sketch size: 6,286 bytes (of a 32,256 byte maximum)


RAM usage unknown, presumably not much.

Please post technical questions on the forum, not by personal message. Thanks!

More info:
http://www.gammon.com.au/electronics

Nick Gammon

I was right, it seems. Using WolframAlpha to calculate:

Code: [Select]

PrimePi(1442723) = 110117


Therefore 1,442,723 is the 110,117th prime number. So at two primes a second it should have taken:

Code: [Select]

110117 / 2 / 60 / 60 = 15.29 hours


Please post technical questions on the forum, not by personal message. Thanks!

More info:
http://www.gammon.com.au/electronics

Jack Christensen


Totally useless of course...


No such thing! I will admit to a penchant for "useless" projects myself but I'd also argue that it's the rare project that is really useless. Usually there is something to be learned or some fun to be had.

At work we used to say that no project is useless because it can always be used as a bad example.

But I like this one! Can't wait to try it later.
MCP79411/12 RTC ... "One Million Ohms" ATtiny kit ... available at http://www.tindie.com/stores/JChristensen/

Nick Gammon

I thought you might want to know ...



That's at 5:20 pm the day after it started. So, around 24 hours later.

As you can see, it will be around 5 days before it ticks over to the 8th digit.

This sort of thing might look good in a shop window. ;)
Please post technical questions on the forum, not by personal message. Thanks!

More info:
http://www.gammon.com.au/electronics

Jack Christensen


I thought you might want to know ...
That's at 5:20 pm the day after it started. So, around 24 hours later.
As you can see, it will be around 5 days before it ticks over to the 8th digit.
This sort of thing might look good in a shop window. ;)


About how long between numbers at this point?
MCP79411/12 RTC ... "One Million Ohms" ATtiny kit ... available at http://www.tindie.com/stores/JChristensen/

Nick Gammon

Next morning now. I am up to 4020503 and the time between numbers seems to be 500 mS still.

Please post technical questions on the forum, not by personal message. Thanks!

More info:
http://www.gammon.com.au/electronics

Nick Gammon

By my reckoning the number of tests it has to do in the worst case is:
(?4020503) / 4 = 501

Then in most cases it will quickly find the number divisible by 3, 5, 7 and so on. So probably testing all 501 would be something fairly infrequent.

To test that timing, I wrote this:

Code: [Select]

volatile long foo = 4020503;
volatile long bar;

void setup ()
  {
  Serial.begin (115200);
  Serial.println ();
  unsigned long start = micros ();
  for (long i = 0; i < 501; i++)
    {
    bar = foo % i; 
    }
  Serial.println (micros () - start);
  }  // end of setup

void loop () { }


Output:

Code: [Select]

20924


So, 21 mS to do 501 modulus operations. It's no wonder it can keep up, since I am giving it 500 mS.
Please post technical questions on the forum, not by personal message. Thanks!

More info:
http://www.gammon.com.au/electronics

Nick Gammon

Of course, that would only apply for adjacent primes. Still, the gap between primes doesn't seem to often be much greater than 50 or so.
Please post technical questions on the forum, not by personal message. Thanks!

More info:
http://www.gammon.com.au/electronics

Nick Gammon

http://en.wikipedia.org/wiki/Prime_gap

The average gap between primes is apparently ln(p) where p is the prime we are looking at. So:

Code: [Select]

ln (4020503) = 15.2


Thus I would expect an average gap of 15. And of course, not all of those are being tested (like even numbers).
Please post technical questions on the forum, not by personal message. Thanks!

More info:
http://www.gammon.com.au/electronics

Jack Christensen

Sweet. I need to get one of those 8-digit displays. I'm sure I have the parts to put one together but it'd be nice to have it on a compact board like that.
MCP79411/12 RTC ... "One Million Ohms" ATtiny kit ... available at http://www.tindie.com/stores/JChristensen/

Nick Gammon

I'm up to 8 digits!



Still going strong, it's still doing one every 500 mS.
Please post technical questions on the forum, not by personal message. Thanks!

More info:
http://www.gammon.com.au/electronics

odometer

Those display modules aren't chainable, are they?

How is your chip at dividing long longs? Too bad there isn't a 48-bit integer datatype.

Go Up