Andrew Que Sites list Photos
Projects Contact
Main

March 20, 2010

Space Port - Part 2 - One big ship

Outer space is big—really big—and almost entirely empty. One thing done in a lot of sci-fi movies is having the characters fly through an asteroid field, dodging asteroids like road hazards. This is silly. You would be lucky to find one asteroid in the Asteroid Belt flying straight through it, yet alone having to worry about dodging such things. Everything in our solar system is separated by distances much larger then anything we can normally work with.

In my story, ships are regularly required to travel to the outer edge of the solar system. Pluto's outer edge of it's orbit is about 50 AU, or 7.3 billion km, or 4.5 billion miles from the sun. So how long would it take to get there?

If you could get a space craft to travel 500,000 mph (not too crazy for space travel), it would take about 1 3/4 years to travel the distance. It you would like to do that in one month, you would have to travel about 10 million miles an hour.

Rocket engines are a source of force. Force in physics is usually measured in Newtons, which is kg·m/s2.

Assuming the force exerted by a rocket engine is constant, we can determine it's acceleration (m/s2) on some given mass (kg). Our mass is the space ship, so knowing it's mass, we can determine acceleration.

What's nice about acceleration is that it keeps increasing velocity. The standard physics equation x = x0 + v0 t + 1/2 a t2 will determine a position (x) after some time (t) with an initial position (x0) and initial velocity (v0). We we start of position 0 (x0 = 0) at a complete stop (v0=0), the equation simplifies to x = 1/2 a t2. And we can solve for t, so t = sqrt( 2 * x / a ). We know x—that is the distance we want to travel. So what is the optimal acceleration? Acceleration due to gravity naturally. This way, we have linear acceleration and perfect gravity. So acceleration is 9.81 m/s2. Thus, if you had the ability to provide a constant force with this rocket engine, it would take about 14 days to reach the outer orbit of Pluto from the Sun. That's not bad considering the distance.

But in working this out, I was faced with one obvious problem. The ship I designed is rotating. It already has gravity. The gravity introduced by linear acceleration will be in the wrong direction. If the rotation provides 1 g in one direction, and linear acceleration provides 1 g in an other direction, people are going to be standing sideways just to stay upright. This isn't going to work.

It was in the car ride back to Wisconsin talking this over with Mikala that a solution was presented. Mikala pointed out what should have been obvious: why not just have the spheres rotate? As the ship accelerates, turn the spheres such that the center of gravity at a 90 degree angle to the floor. In fact, one wouldn't even need to turn the spheres. The bottoms are heavier then the tops, because the tops are empty. If allowed to rotate, the spheres would naturally correct for the center of gravity.

So as the ship accelerates, the rotation is slowed. Since the ship is likely to be large, we will have to gradually accelerate to avoid jarring everything. A change in acceleration is known as jerk. So some constant jerk factor would have to be obtained—likely from that rate at which the rotation could be slowed to a stop.

March 15, 2010

Space Port - Part 1

A couple of weeks ago, I started contemplating an idea of a large transport space ship. On occasion I like to write a little sci-fi, but I'm usually far more interested in the technology then the story. I started to wonder about how a large transport vessel might be laid out.

Not much of a fiction reader, most of my exposure to sci-fi is from movies. And most movies don't try all that hard to create images of accurate space traveling vessels. Most somehow assume there is artificial gravity so they never have to address the issue. But how do you get gravity?

There are two methods that are rather easy: rotation and linear acceleration. Rotation works well if the object is stationary, such as with a space station. Linear acceleration can work if you have the ability to accelerate at 9.8 m/s2—the same acceleration as gravity. That's an acceleration of 79,000 mph2—meaning if you start for a dead stop, after 1 hour would be traveling 79,000 mph—so a fairly hefty acceleration. So I started with rotation.

Rotation provides centripetal force. Most of us as a kid have either tried or have seen someone spin a bucket of water over their head. The centripetal force from the rotation causes the water not to spill out even as the bucket is inverted. One can apply the same principals in space.

Again, most movies get a rotation as a means of artificial gravity wrong. Gravity is a function of the radius of the rotating path, and the speed of rotation. g = R( pi r / 30 )2/9.81 where g is gravity with respect to Earth (so 1.0 is Earth's gravity), R is the radius of rotation and r is rounds/minute of rotation. Two things should stand out about this function: a faster rotation means much higher gravity because the RPM is squared, and gravity depends on the radius. If you picture a circle, the one would be standing on the inside edge. One's feet would be touching the inner edge, and one's head would extend toward the center of the circle. The head, being closer to the center of the circle, would experience less gravity then their feet, because the radius has been reduced by their body length. This effect would be significant for small circles. For example, we could use a radius of 10 meters (so 20 meters, or about 66 feet in diameter). To get gravity at the circles edge would require a rotational speed of about 9.5 RPM. A person who is 5'8" (1.7 meters) tall would only have 82% gravity at the top of their head. Someone lifting 10 pounds would to their head would find it weighting only 8 pounds when level with their head. This may not be terrible. However, there is an other problem. Humans experience the Coriolis effect (i.e. knowing that we are spinning) at speeds above 2 RPM. So someone in such a station would get dizzy.

To get the speed down, it is necessary to increase the radius. One concept studies back in 1975 referred to as the "Stanford torus" was a torus (i.e. the doughnut shape) 900 meters in radius (about 1 mile in diameter). This would require about 1 RPM in order to hold gravity. The "O'Neill cylinder - Island Three" had an outer radius of 16 km (about 20 miles) and would rotate at one round ever 4 1/4 minutes.

After a lot of reading (mostly on Wikipedia), I starting thinking about how to create my own living area. The "Island Three" concept was really neat because it is an open cylinder. You could look up and see the other side—a city scape above your head but not falling on top of you. But this concept was designed to be stationary—orbiting Earth for example. They used large mirrors to reflect in sunlight, which was used for lighting and plant life. My vessel had to move, and the further away from the sun, the less light you get.

An other problem I had with this design was the possibility for catastrophic disasters. The entire atmosphere was contained in a single section. One hull breach, and you risk destroying the entire station. I wanted something more compartmentalized—something that worked in sections that were able to be isolated for other sections.

In other stories I've worked on, I thought up a concept using domes to grow crops. A center light source would shine on a reflective ceiling, and fields would be cultivated by a rotating gantry crane. I wondered if this concept couldn't be expanded to the space port. Rather then just use a dome, I tried using a sphere. It was easier to draw, but after some thought, made more sense. A dome is just a sphere cut in half. So the upper half of this sphere would be the dome, and the lower half could be support equipment. The flat disk in the middle is the plane for which people could live. This idea seemed to work pretty well in my head, so I started to come up with numbers—which is really where I have my fun!

My ring would be 30 km in radius (or 37 miles in diameter). At this size, it would rotate once every 5 3/4 minutes (0.17 RPM). Each dome would be 1 km in radius (1 1/4 miles in diameter), providing 776 acres of land (31 km2). The whole ring could fit 90 such domes for a total of almost 70,000 acres—109 sq. miles (2,800 km2). Because of the large size, even a building 100 meters tall (328 feet, or around 32 stories) will still have 99.67% gravity—so a 150 pound person would weight 149.5 pounds on the top floor.

So far, everything was going well for my space port. But then I remember this thing had to move...

March 09, 2010

Modeling my room

   Yesterday I completed my basic bench layout, and for the most part, everything fits.  Part of the reason is that I used Google Sketchup to do a drawing for placement first.
   The sleeping box has been redesigned.  At first I considered cutting the box down to bench height.  Our new place is quite small.  I instead to use the top of the box as a bunk for our roommate James.  So this leaves the box at a height of 4 feet.  The biggest change is the door.  Previously, the box opened by hinges in the front, done mainly to keep from having to make additional cuts when the box was first assembled.  This was not a good design.  It leaked a lot of light, and was quite heavy.  The sliding door of my box at the Garage works much better.  Since Jai sleeps in the box with me, I decided to recycle the sliding door, but add a second one.  We both go to bed at different times, and this allows her to get into the box whenever she likes without having to crawl over me.

Here is a sketch of the box with one side removed, and the closest door transparent.

   In order to make the new sliding door system, I had to redesign the frame.  My box at the Garage uses a double wall setup, with the door in between.  This is functional, but I used a 2x4 for the gap, resulting in 1.5 inches of space for the 5/8" door.  This means the door flops around and isn't very snug.  I decided this new setup should be a little better.  I omitted the inner wall, mainly because I didn't want all the extra wood it would require.
   So, I build a 2x4 skeleton, edged with 1x3" on the top, left and right, and 1x2" on the bottom.  This frame provides the gab for the door to slide long.  The door is just 1/2" plywood like the sides and top.  It turned out quite well actually.  The doors move very smooth in the space allotted, and close up tight.  After some caulking, the system should be pretty light tight.
   I finished the main construction on Thursday, left for the weekend starting Friday, and caulked on Sunday night.  Yesterday I installed the bench setup, and finished the ventilation.  Still haven't had a full test of light though.
   Tonight I decided to finish up my model.  I measured and drew the bedroom.  I added a the monitors and chair, but most of the items in the bedroom I couldn't find good models for.

Pictured is the sleeping box and my work bench as seen from the north east corner of the room.


Floor plan.

4 comments have been made.

From Erica

March 13, 2010 at 9:05 PM

Sweet! That looks really nice.

From Noah

Cyberspaces

March 15, 2010 at 5:20 PM

It's a nice layout - the only comment I might make is that I prefer to keep a wall to my back when working on the computer. Although perhaps you don't have quite such a "nosy retard" problem as I do. Also, I might add this bit: you've mentioned problems with your doors sticking. It never occurred to me at the time, but I'd guess that a hard wax like those clear paraffin blocks would be ideal for a nonperishable lubricant!

From Erica

March 25, 2010 at 9:41 PM

lmfao @ "nosy retard"

From Liz

Janesville

March 30, 2010 at 5:23 PM

Are there many jobs in Cedar Rapids? After I get my associates I may move somewhere out of state, I'm looking for a town that has work and decent prices for apartments. The kansas city plan that my friends started went out the window.

February 27, 2010

Internet Charlatans

   Having owned and operated a domain for over 5 years I'm use to people trying to sell me junk because the idiots shysters at Network Solutions won't keep my e-mail address private unless I pay for that "privilege".  I've received e-mails in the past from people wanting to sell me a domain name similar to mine.  I ignore them by changing my e-mail address.  So when I received an e-mail the other day, I followed standard procedure.  To my surprise, I received a second e-mail hour latter to the updated e-mail address.  Seems the domain squatters have become more aggressive.  To humor myself, I clicked on the "unsubscribe" link.  Naturally, I received an e-mail today, this time from a different company.
   Let's have a look at the e-mails.  The first one is from 2/20/2010 9:19 AM
Priority Domain Availability Notice for drque.com:

In the next few days, drque.com will be listed for sale.  Since you have a similar domain name, we thought you might be interested in acquiring drque.com.

You can confirm your interest in the domain drque.com by filling out the form here:  http://initrustdomain.net/store/prioritynotice/f5e6fffd-fb88-456a-8d65-c3568bc2522c

After we receive a confirmation that you are interested in the domain drque.com, our sales staff will be in touch with you promptly to make arrangements.

We look forward to hearing back from you.

Kind regards,
InTrust Domains
11605 Meridian Mkt V #124-134
Falcon, CO 80831

 

If you are not interested in future priority notices, please simply reply to this message

No more please: http://initrustdomain.net/store/unsub/f5e6fffd-fb88-456a-8d65-c3568bc2522c
 
The next from 2/22/2010 7:13 AM
Domain Sale Notice:

drque.com is coming available for sale in a few days.

Since you own the domain drque.net, we thought you'd be interested in drque.com.

If you do have interest in acquiring drque.com, please fill up priority notice form availble

here: http://initrustdomainnames.com/store/prioritynotice/f5e6fffd-fb88-456a-8d65-c3568bc2522c

and we will contact you as soon as the domain is available for purchase.

 

We look forward to hearing back from you.

 

Kind regards,

John Timmers

InTrust Domains

11590 Black Forest Road

Colorado Springs, CO 80908-6000

 =========

 

This e-mail is for the sole use of the intended recipient(s) and contains information belonging to InTrust Domains. and/or John Timmers, which is confidential and/or privileged. If you are not the intended recipient, you are hereby notified that any disclosure, copying, distribution or taking of any action in reliance on the contents of this e-mail information is strictly prohibited.  If you have received this e-mail in error, please immediately notify the sender by reply e-mail and destroy all copies of the original message.

 

No more please: http://initrustdomainnames.com/store/unsub/f5e6fffd-fb88-456a-8d65-c3568bc2522c
I received this e-mail at 2/23/2010 7:06 AM
Domain Sale Notice:

drque.com is coming available for sale in a few days.

Since you own the domain drque.net, we thought you'd be interested in drque.com.

If you do have interest in acquiring drque.com, please fill up priority notice form availble

here: http://bestinitrustdomains.biz/store/prioritynotice/f5e6fffd-fb88-456a-8d65-c3568bc2522c

and we will contact you as soon as the domain is available for purchase.

 

We look forward to hearing back from you.

 

Kind regards,

John Timmers

InTrust Domains

11590 Black Forest Road

Colorado Springs, CO 80908-6000

 =========

 

This e-mail is for the sole use of the intended recipient(s) and contains information belonging to InTrust Domains. and/or John Timmers, which is confidential and/or privileged. If you are not the intended recipient, you are hereby notified that any disclosure, copying, distribution or taking of any action in reliance on the contents of this e-mail information is strictly prohibited.  If you have received this e-mail in error, please immediately notify the sender by reply e-mail and destroy all copies of the original message.

 

No more please: http://bestinitrustdomains.biz/store/unsub/f5e6fffd-fb88-456a-8d65-c3568bc2522c
Now I change my e-mail address.  The next e-mail is to the new address at 2/23/2010 6:01 PM
Hey

My name is Eric Taylor from Flex Media.
drque.com will be available for a limited time purchase in a few days.
Since you own drque.net, I thought you might be interested in drque.com.

If you do have an interest in purchasing drque.com, please visit:
http://www.flexmediadomains.com/prioritysales/4b846b3f090cf12526ef82f704f38b09fdbc9deccaad84b84605bd89513c5d12580d3cd5dba07ea060e039ae39

If the link above doesn't work, please copy the address below and paste it into your web browser:
http://www.flexmediadomains.com/prioritysales/4b846b3f090cf12526ef82f704f38b09fdbc9deccaad84b84605bd89513c5d12580d3cd5dba07ea060e039ae39

Once drque.com is available for purchase, one of my account specialists will contact you.

My goal is to make this the easiest transaction possible for both of us.
I can guide you step-by-step over the phone, instant messenger or email on how to secure drque.com if you wish.
If you have any concerns or questions, please do not hesitate to contact me.
No question is a dumb question.

Thank you and have a great day.

Eric Taylor
Flex Media Inc
7777 Davie Road Ext, Suite 106B
Hollywood, FL 33024
954-584-3750
sales@flexmediadomains.com
Skype: flexmediadomains
ICQ: 576313997
AIM: sales@flexmediadomains.com
MSN: flexmediadomains@live.com
www.flexmediadomains.com

This e-mail is for the sole use of the intended recipient(s) and contains information belonging to Flex Media Inc and/or Eric Taylor, which is confidential and/or privileged. If you are not the intended recipient, you are hereby notified that any disclosure, copying, distribution or taking of any action in reliance on the contents of this e-mail information is strictly prohibited. If you have received this e-mail in error, please immediately notify the sender by reply e-mail and destroy all copies of the original message.

No more please: http://www.flexmediadomains.com/opt_out.php?prerelease=&id=4b846b3f090cf12526ef82f704f38b09fdbc9deccaad8
And the last e-mail from 2/26/2010 11:26 AM
Hi,

We thought that you definitely need to know that right now the .com version of your domain is secured by Easy Domain Recovery and is available for purchase at our page.

What we do is take different measures to help you promote your brand on the Internet. Little things' effectiveness can compare to securing .com to a domain with any other extension. Everyone knows .com, a lot of users will intuitively seek you at .com extension so don't give a chance to anyone else to use that advantage of yours.

There are third parties who may be interested in your domain. Many people out there reserve .coms to existing sites and domains in order to make profit of 'em. We may sell the domain to any of these persons and businesses but this is you who decides whether or not we do so. Please choose one of the following:

Secure drque.com right now!
Let it go (unsubscribe).

Best regards,
Tom Harris,
Easy Domain Recovery,
mailto: tom@easydomainrecovery.com 
    Some interesting things can be learned.  First, the use of various domain names to push their products.  They use initrustdomain.net, bestinitrustdomains.biz, flexmediadomains.com and easydomainrecovery.com.  Both initrustdomain.net and bestinitrustdomains.biz are owned by "Moniker Online Services LLC" who according to "SPAM Trackers" wiki have a long history of criminal activity.  The other two domains are "owned" by "DOMAINBULLIES, LLC" and "ENOM, INC."  There isn't much information about them.
   What gets me is that the Internet is a major source of commerce.  Isn't the FTC suppose to be doing something about this?  Why then do places like SPAM Tracker talk about groups like Moniker Online Services and how they have operated for years?  The more I deal with law and law enforcement, the more less I believe it functions.  It is both incapibal of protecting people, and incapibal keeping law enforcement from harassing the people it aims to protect.

3 comments have been made.

From Noah

March 01, 2010 at 4:30 AM

How do they know your current/active email address? Are they simply spamming emails to *@drque.net or is Network Solutions selling that information to these drips?

From Erica

March 13, 2010 at 9:07 PM

That is so pathetic, it's funny.

From Andrew (http://www.project-insomnia.com/)

Palo Alto, CA

November 10, 2010 at 6:07 PM

These scammers are trying to sell me an alternate spelling of my domain name for almost $400. I honestly don't know what to do; I'm not about to pay the danegeld, but I don't like the idea of these criminals holding what should be mine.
    Everyone once and awhile I want to quick put a 3D drawing together, and I'll use Google Sketchup. Being the math nut that I am, I love it when I need to do a bit of geometry.
   Today, I happened to need to place an octagon (8=sides parallelogram) such that northern, and eastern edges were also the edges of a square, and the start of the north west corner (X), and south east corner (Z) of the octagon began at the these corners on the square. When you draw a polygon in sketchup, you need to have a center point (I) and an end point at one of the angles (such as A). So if I wanted an 8-sided polygon (octagon) to calculate where the center points, and some end point. If you happen to have Sketchup, try and see if you can eyeball this--you'll find it isn't easy. A geometric solution is needed. Luckily, it's easy.
  Every angle in an octagon is 45 degrees, which is 360 / 8. To find the center of an octagon, pick any edge, and measure 67.5 degrees (45 / 1.5) toward the center. An example in this drawing is line AB or BC. Where the angles meet is the center, point I. Note also that when drawing this polygon, it will be done from point I to one of the point A through H.
   In my requirement description, I want this octagon (call it octagon I) to meet up with square WXYZ such that vertex A meets vertex X, and D meets Z. Sides AB to XY, and DE to ZY and parallel, and vertex A and the same as X, D the same as Z. Now that we know the 2 verities, and two edges, we can find the center. Measure an angle starting at X along XY that measure 67.5 degrees clockwise. And 67.5 counter clockwise starting Z lone ZY. This is the center point.
   Draw your octagon starting at this location I, going to either vertex X or Z. That's it! :)
  

February 03, 2010

Prime number - Part 2

So my last article on prime numbers led me to ask the question "How many prime numbers are there between 2 and 232". I know the quick answer to this question is "a lot". 232 is over four billion—and that's a fair bit of space to check. It would be easy to make a for-loop and brute force the entire number range, but this is the age of multi-core CPUs—and single threaded is so 2000. Besides, my computer hasn't broken a sweat in awhile.

I decided to do a straight C implementation. I've been doing C and C++ at work, and C++ has been irritating me lately because I tried to do a "quick and dirty" implementation without too much thought into the classes I'd need. That's a recipe for disaster, and how you end up with things like microsoft's windows ME. So vanilla C it was.

Threads and semaphores I used the POSIX standards, and I've tried to demonstrate how to make a dispatcher foreground task with one or more worker threads. The key to this is the dispatch semaphore. Normally, semaphores are used to for mutually exclusive lock, or to send signals between threads. However, they can also be used as a count. Think of them like a library that has a set number (say 5) of books. That means 5 people can check that book out, but if a 6th person wants to check that book out, they will have to wait for one of the first 5 to return their copy. The semaphore does this by keeping a count. When that count reaches zero, the next task that tries to pend the semaphore will have to wait. So the dispatcher is a loop that starts by pending the semaphore. When it gets the semaphore, it will create a task. The task will do some work, then post the semaphore. The semaphore count determines how many tasks can run at the same time. If we limit the number of tasks that can run at once to the number of CPU cores, then we can utilize the CPU to 100%.

First, the dispatcher skeleton code:

  pthread_t Threads[ NUMBER_OF_THREADS ];

  sem_init( &Semaphore, NUMBER_OF_THREADS, NUMBER_OF_THREADS );

  unsigned ThreadIndex = 0;

  while ( NumbersLeft )
  {
    // Wait for a free worker thread.
    sem_wait( &Semaphore );

    // Create a worker thread to check this number set.
    pthread_create
    (
      &Threads[ ThreadIndex ],
      NULL,
      PrimeThread,
      (void *)&Data[ ThreadIndex ]
    );

    ++ThreadIndex;
    if ( ThreadIndex >= NUMBER_OF_THREADS )
      ThreadIndex = 0;
  }

Then, the skeleton of the worker thread:

static void * PrimeThread( void * ArgumentPointer )
{
  //
  // Do prime check.
  //

  // This thread is now complete.  Release one count from the dispatch
  // semaphore.
  sem_post( &Semaphore );

  // End this thread.
  pthread_exit( 0 );

  return 0;

} 

Now just launching a thread to check one number is a bit of a waste—it takes some overhead to start and stop the task. So we give each thread some range of numbers to check. A pthread can be passed parameters, and we pass the starting number, along with the amount of numbers to check. The parameter block is also used for storage of the results. So here is the full worker thread:

static void * PrimeThread( void * ArgumentPointer )
{
  // Get the work data passed to the thread.
  WORK_TYPE * Data = (WORK_TYPE *)ArgumentPointer;
  uint32_t Number = Data->StartNumber;
  uint32_t Count = 0;
  unsigned Index;

  // For all the numbers to check...
  for ( Index = 0; Index < Data->NumberToCheck; ++Index )
  {
    // Is this number a prime?
    if ( IsPrime( Number ) )
      // Then count it.
      ++Count;

    // Next number.
    ++Number;
  }

  // Save results.
  Data->NumberFound = Count;

  // This thread is now complete.  Release one count from the dispatch
  // semaphore.
  sem_post( &Semaphore );

  // End this thread.
  pthread_exit( 0 );

  // Never reached--here for language consistency.
  return 0;

} // PrimeThread

The dispatcher has to accumulate the results. Our could try using a global variable, but you can run into issues with non-atomic access. Adding 1 might seem a basic operation, but (at least in a RISC system) it takes 3 instructions to add one to a variable in memory: load the data to a register, add the register by 1, and store the register. You can not (without a semaphore) guarantee that those instructions will not be interrupted by an other thread. So, the accumulation happens in just one place—the dispatch loop.

The dispatcher will assign all the work, but it is finished once all the work has been assigned. We still need to wait for the work to finish. So the last part of the process looks like this:

  for ( ThreadIndex = 0; ThreadIndex < NUMBER_OF_THREADS; ++ThreadIndex )
  {
    if ( Data[ ThreadIndex ].NumberToCheck )
    {
      // Wait for thread to finish.
      pthread_join( Threads[ ThreadIndex ], NULL );

      // Accumulate the number of primes found in this thread.
      NumberOfPrimes += Data[ ThreadIndex ].NumberFound;
    }
  }

  // Let go of dispatch semaphore.
  sem_destroy( &Semaphore );

Now we can print the results to the screen, and we're done.

I knew this code would tax the CPU rather hard, so I wanted to measure how long the process would take. "time.h" is a C standard unit. Unfortunately, "time_spec", which has a high-resolution method for checking the time, it's part of the C99 standard. But we can get time in the resolution of seconds, and that's good enough. To do this, we simply mark the time when the process begins, and mark it again at the end. The difference (for which there is a function to compute) is the total number of seconds that has elapsed.

The answer: There are 203,280,221 primes between 2 and 4,294,967,295, and it took 6,054 seconds (1 hour, 41 minutes) on a quad-core clocked at 2.5 GHz.

The full source code is here.

1 comment has been made.

From Antony (http://rapid4me.com)

LA

February 17, 2010 at 3:10 AM

Excellent post! I always enjoy a solid technical post (and code) It saved me a good week - Keep up the good work!

February 02, 2010

Prime numbers - Part 1

So I mentioned yesterday that I wrote a simple program to check to see if some numbers in a list were prime or not. I didn't post it because I wanted to spend a little more time with it—giving it an article on it's own.

Finding prime numbers goes back a long time. The easiest way to see if a number is prime is to simply check and see if the number is divisible by all the prime numbers up to the square root of the number in question.

A number (call it x) is prime if there are no two number (call them a and b) such that a * b = x. All non-prime numbers can be expressed as the product of two or more prime numbers. For example, 125 can be made from 5 * 25, but 25 can be made from 5 * 5. So 5 * 5 * 5 = 125, and represents the most factored version of 125. This is true of any number. Since it takes at least two prime numbers to create a factor, we only need to check the primes up to √x (or the square root) of the number. This is because the if x = a * a, then √x = a. If x = a * b, and b is greater a, then a must be less then √x. Thus, we only need to check primes up to √x.

This makes for a nice brute-force test method, and it's pretty simple to implement. If we're going to test a number, we need all the primes up to the square root of that number. For what I needed, the number was capped at 32-bits—or 232. √232 = (232)1/2 = 216 = 65536. So I need a list of all the prime numbers up to 216. We can make the list fairly quick. Have a look at this function:

enum { NUMBER_OF_LOOKUP_PRIMES = 6542 };
unsigned PrimeNumbers[ NUMBER_OF_LOOKUP_PRIMES ];

void GeneratePrimeNumberLookup()
{
  unsigned Index;
  unsigned PrimeNumberCount = 0;

  PrimeNumbers[ PrimeNumberCount++ ] = 2;

  for ( Index = 3; Index < 0x10000ul++Index )
  {
    bool IsPrime = true;

    unsigned SubIndex = 0;
    while ( ( SubIndex < PrimeNumberCount )
         && ( IsPrime ) )
    {
      // Does it divide evenly by this prime number?
      if ( 0 == ( Index % PrimeNumbers[ SubIndex ] ) )
      {
        IsPrime = false;
        break// <- We can stop checking.
      }

      ++SubIndex;
    }
    

    if ( IsPrime )
      PrimeNumbers[ PrimeNumberCount++ ] = Index;
  }

}

We have an array to hold our prime numbers. We set the first entry in the array to two—the first prime number (1 is prime, but who cares). From this, we can test to check the rest of the numbers. We loop from 3 to 216. For each number, we see if any of the primes already in the list will divide evenly into the number under question. If any of them do, they are not prime.

It happens that there are 6542 prime numbers between 2 and 216. Knowing that, we can cap the array size of our lookup table at this value—the number of primes in this range will never change.

We now have enough prime numbers to check any number up to 232. Unlike the function that generates the prime lookup table, we can throw in one more speed improvement: we only need to check up to the square root of the number in question. So here are the guts of the prime test function:

bool IsPrime( uint32_t Number )
{
  bool IsPrime   = true;
  uint16_t Root  = SquareRoot( Number );
  unsigned Index = 0;

  while ( ( Index < NUMBER_OF_LOOKUP_PRIMES )
       && ( PrimeNumbers[ Index ] <= Root )
       && ( IsPrime ) )
  {
    if ( 0 == ( Number % PrimeNumbers[ Index ] ) )
      IsPrime = false;

    ++Index;
  }

  return IsPrime;
}

Note that can add two short cuts. First, checking up to the root of the number like we discussed. And second, we test to see if the number is even—2 is the only even number that is prime.

That's it. Now we have a simple to implement prime number test that works on numbers up to 232. Here is the full source code to a command line version of the test. Compile it, pass it a number (or several), and the program will tell you if the number is prime or not.

 

(thanks Erica for the correction)

February 01, 2010

A good 32-bit LFSR polynomial

At work on Friday I came across an interesting math problem. I needed a pseudo random bit stream, and I needed the generator function to be fast. The typical random number function is a Linear congruential generator, most often f(x) = (1103515245 * x + 12345) mod 232-1. The modulus remainder is implied—the arithmetic is simply allowed to overflow. Multiplies are usually expensive in terms of CPU cycles. A faster generator is a linear feedback shift register (LFSR), which are comprised of shifts and XORes. I've used LFSR a number of times when I needed random data, without looking too much at the "random" stream of data it produced.

Using a 32-bit Galois LFSR with the polynomial x32 + x31 + x29 + x + 1, and a seed of 0x12345678 produced this output:

2B3C091A 159E048D 8ACF0246 4566D123
A2B36891 D158E448 68AC7224 34563912
1A2B1C89 8D14DE44 468A6F22 23453791
91A3CBC8 48D0B5E4 24685AF2 12342D79
891A16BC 448D0B5E 2247D5AF 9122BAD7
C8915D6B E449FEB5 F225AF5A 791387AD
BC8893D6 5E4449EB AF2224F5 D791127A
6BC9D93D B5E5BC9E 5AF38E4F AD789727
D6BC4B93 EB5E25C9 F5AF12E4 7AD6D972
3D6B6CB9 9EB4E65C 4F5A732E 27AD3997
93D7CCCB C9EAB665 E4F55B32 727BFD99
B93CAECC 5C9E5766 2E4F2BB3 9726C5D9
CB9362EC 65C8E176 32E470BB 9972385D
CCB91C2E 665DDE17 B32FBF0B D9968F85
ECCB47C2 7664F3E1 BB3279F0 5D993CF8
2ECDCE7C 1767B73E 0BB28B9F 85D945CF
C2EDF2E7 E177A973 F0BA84B9 F85D425C
7C2FF12E 3E16A897 9F0B544B CF84FA25
E7C27D12 73E13E89 B9F1CF44 5CF9B7A2
2E7D8BD1 973F95E8 4B9E9AF4 25CF4D7A
12E6F6BD 89737B5E 44B8EDAF A25C76D7
D12E3B6B E8971DB5 F44ADEDA 7A256F6D
BD13E7B6 5E88A3DB AF4451ED D7A228F6
6BD1147B B5E9DA3D DAF5BD1E 6D7B8E8F
B6BC9747 DB5E4BA3 EDAF25D1 F6D6C2E8

If you haven't noticed, the pattern "scrolls" digits to the right. This is what one would expect—this is a "shift" register after all. Each number is different, and the period between any two number (except for zero) is 232. However, the pattern between each successive number is fairly obvious. This is because of the polynomial. The shift register only introduces bits 31, 29, 1 and 0 when clocked. Otherwise, output is simply shifted to the right by one. I've arranged the data in 4 columns because each digit (nibble) is 4 bits, so after 4 bits, the nibble is completely shifted right by one.

While this pattern probably would work fine for what I was trying to do, I thought there must be a better solution. A polynomial with more taps should produce a better random stream. But how exactly does one go about generating a polynomial for a shift register?

What is need is a irreducible primitive polynomial. One can put any polynomial (i.e. any number as the polynomial term) into a LFSR they like, but it may not (probably won't) produce a maximum length shift register. That is, a 32-bit word has 232 possibilities. With the correct polynomial, the output of a 32-bit LFSR will produce 232 unique outputs before repeating. Without the correct polynomial, one of two things will happen. The output will happen after less then 232 iterations. Or the output will never repeat and some smaller sub-sequence will repeat.

The polynomial of a LFSR is XORed in with the shifted data. My thought was, if I wanted a better random pattern, we'd better XOR more bits. And not just more bits, but spread out more or less evenly. This should change the stream of data quite a bit.

So how does one create an irreducible primitive polynomial? Good question. I still don't have a complete solution, but I do have part of the solution. I found a program written by a Scott Duplichan that generates primitive polynomials—although not necessarily irreducible. It's designed to generate huge primitive polynomials—on the order of thousands of bits. However, it also would do 32-bit. What was nice about this program was that I could tell it the "weight" I wanted. That is, how many bits I wanted set. The number has to be odd, so I chose to generate polynomials with 15 bits.

There are a lot of primitive polynomials with a weight of 15—millions. And most of them did not have bit patterns that were uniformly distributed. For example, one of the first polynomials produced was 0x00041FFF. While there are 15-bits set, there are pretty much all on one side. The program didn't have a way to tell it to generate polynomials that had uniform distribution, but it did have the ability to find a set number of polynomials. That was good enough for me. I set the program off to find ten million primitive polynomials, and went off to do something else. Just over an hour latter (4152 seconds actually), I had 10 million primitive polynomials with 15-bits set.

The next task was to find polynomials with a more uniform bit distribution. To do this, I wiped up a filter program to only save polynomials that had, at most, 3-bits of the same bits in a row. This significantly reduced the number of polynomials. But I still did not have a list of irreducible primitive polynomials. I started reading about how to test for this, but it was taking too much time. It would be easier to brute force this problem. I wanted a polynomial that was maximum length—that is, produced 232 unique values before repeating. Any polynomial that did this was what I looking for.

It only takes a few seconds to brute-force check all 232 possibilities for a LFSR. With the bit distribution restrictions, there were only 9,441 of the original ten million polynomials to check. So, I launched a brute-force attack, and went home for the weekend.

When I came back this morning, I had a list of 40 irreducible primitive polynomials. I'm not sure how long it took to complete this search, but that was irrelevant—I had my data.

Forty was still more polynomials then I needed. So I decided to reduce the list a little further. I just wanted prime polynomials. This is an other brute-force test. Any number is prime if all none of the prime number up to the square root of the number divide evenly into it. My list was reduced from 40 to 10. I picked the first number from the list, but any would have worked.

Again, my 32-bit Galois LFSR, this time with the polynomial 0x19253292B, and a seed of 0x12345678 produced this output:

2B3C091A 159E048D 8ACF0246 6C4C9370
362649B8 1B1324DC 24A2803D 9251401E
6003B25C 192ACB7D 8C9565BE 6F61A08C
1E9BC215 A666F359 D33379AC 40B2AE85
A0595742 7907B9F2 15A8CEAA 0AD46755
856A33AA 42B519D5 88719EB9 ED13DD0F
DFA2FCD4 6FD17E6A 1EC3AD66 264AC4E0
13256270 20B9A36B B977C3E6 7590F3A0
3AC879D0 1D643CE8 0EB21E74 07590F3A
2A8795CE 3C68D8B4 1E346C5A 0F1A362D
878D1B16 6AED9FD8 1C5DDDBF A705FC8C
7AA9EC15 947FE459 E314E07F F18A703F
F8C5381F D5498E5C 438FD57D 88ECF8ED
C4767C76 623B3E3B 98368D4E 4C1B46A7
8F26B100 47935880 0AE2BE13 85715F09
EB93BDD7 DCE2CCB8 6E71665C 1E13A17D
A622C2ED D3116176 40A3A2E8 097AC327
84BD6193 EB75A29A 5C91C31E 0763F3DC
2A9AEBBD 954D75DE 638DA8BC 18EDC60D
A55DF155 FB85EAF9 D4E9E72F C35FE1C4
4884E2B1 A4427158 522138AC 003B8E05
A936D551 D49B6AA8 4366A707 A1B35383
F9F2BB92 7CF95DC9 9757BCB7 E280CC08
71406604 38A03302 1C501981 8E280CC0
47140660 238A0330 11C50198 21C9929F

And this looks much better. But just because you can't see a pattern doesn't mean one doesn't exist. In fact, one obviously exists—the data is from an equation. But without the equation, how good does the random data look? The quickest way to find out if data is random is to try and compress it. If it compresses, it isn't very random. For example, with 256 MBs of data from the original polynomial, we gets a stream that 7-zip could compress 54% of it's original size. With the new polynomial, 7-zip was only able to compress it to 99.964% of it's original size.

Compression isn't the best judge of a file's randomness. I knew there were better algorithms out for generating an entropy number. Doing a quick search, I found this program originally written by John Walker—one of the co-founders of AutoCAD. It dates back to 1985, but the nice thing about math (unlike computer hardware) is it still works fine after 25 years.

The program produces the chi-square distribution, which is given as percentage. When this percentage is close (or near) to 0% or %100, the data stream is said to not be very random—even 95% or 5% are suspect. The chi-square distribution of the first data set is 0.01%, or not even close to random. The chi-square distribution for the new polynomial is 38.58%. For comparison, the chi-square distribution of some output of the XTEA encryption algorithm in a feedback loop had chi-square distribution of 26.5%. Encryption algorithms desire to produce highly random output.

So our new polynomial of 0x9253292B produces a good psudo-random data stream when used in a 32-bit LFSR.

   It might be a while before I post pictures.  This evening, my camera started acting really strange.  As soon as I put the battery in it, the shuttle started snapping--even with the camera turned off.  I tried several things, but nothing seemed to work.  Once and awhile, I got the message "Err 99".  After some reading on the web, I found this is a fairly common problem with the 20D when the shutter system dies.  Seems the average cost of repair is around $200.  So the question is, repair or upgrade?

4 comments have been made.

From software development in Surrey (http://www.geeks.ltd.uk/Services/software-product-development.html)

http://www.geeks.ltd.uk/Services/software-product-development.html

January 19, 2010 at 4:03 AM

That was an inspiring post, Cool piercings, Thanks for writing, most people don't bother.

From Talon

USS Enterprise The Mobile Chernoble

January 27, 2010 at 6:30 AM

UPGRADE! (if it isn't too late for my opinion)

From Erica Too Lazy Too Log in

January 31, 2010 at 5:27 PM

Moar liek "Dies Irae," amirite?

From Application developers in Sri Lanka (http://www.oganro.com/software-development-sri-lanka)

Sri Lanka

July 10, 2016 at 4:10 AM

Not sure how we can join as well?