0. Code

>	restart; with(plots):

Warning, the name changecoords has been redefined

>	Primes  :=  proc( N )     local n, m, M;     n := min(14, floor( sqrt(N)*1.5 ));     m := ceil( N/n);     print(cat(`The first `, n*m, ` primes : `)); print(` `);     array( [seq( [  seq( ithprime( (j-1)*n+i), j = 1..n ) ],   i = 1..m) ]);  end proc:

>

Sieve  :=  proc( N )     local n, m, A,M, p,i,j,k, K;     n := min(20, floor( sqrt(N)*1.5 ));     m := ceil( N/n);  M := n*m;     A:= array( [seq( [  seq( (i-1)*n+j+1, j = 1..n ) ],   i = 1..m) ]):     print(cat(`Here are the integers from 2 to `, M+1)); print(` `);     print(A); print(` `);print(` `);print(` `);     K := 0;     for K while ( ithprime(K+1) < evalf(sqrt(M))) do end do;          for k from 1 to K do       p := ithprime(k);       for i from 1 to m do           for j from 1 to n do               if(A[i,j]<>p)               then if( irem(A[i,j],p)=0) then A[i,j]:=` ` fi; fi;       od;od;       print(cat(`Here are the integers from 2 to `, M+1,                 ` after removing multiples of `,p,` (except `,p,` itself)`));                             print(` `);       print(A);print(` `);print(` `);     od;     print(` `);     print(`This is the Sieve of Erasathones. EVERY number left in this sieve is a prime number`);     print(`because we have removed all of the composite numbers!!`);  end proc:

>

LargestPrime := proc(p) local n,k, P, N,NS,product; n := 0; for n while ( ithprime(n) < p) do end do; P:= ithprime(n); print(cat(`If `, P, ` (the `,n,`th prime) is the largest prime, then "all known" primes would be :`)); product := 1; for k from 1 to n do     product := product* ithprime(k); print(ithprime(k)); od; print(` `); print(cat(`The product of all these "known" primes is `, product)); N := product +1 ; print(cat(`If we add one to this product, we get the number `, N)); NS :=  numtheory[factorset](N); if ( nops(NS) = 1)    then print(cat(N,` is a "new" prime number itself!`));    else print(cat(N,` is not prime, but it contains "new" primes.`));         print(N = ifactor(N)); fi; end proc:

>

FindSquare := proc(N)   local n,k,FS, square;   if( numtheory[issqrfree](N) )     then print(cat(N,` is squarefree!`));     else        FS := numtheory[factorset](N);        square := 1;   n:= nops(FS);        for k from 1 to n do              if( irem(N, FS[k]^2) = 0) then square := square*FS[k]^2; fi;        od;        print(cat(square, ` is the greatest square that divides `,N));   fi;    end proc:

1. What is a Prime?

In considering the natural numbers other than one, there are two types of numbers :
             - prime numbers  - those numbers which can only be divided by 1 and the number itself, and
             - compositie numbers  - those numbers which can also be factored by other numbers

>	`These are prime numbers : `; 7: %=ifactor(%);  31: %=ifactor(%);  173: %=ifactor(%);

>	`These are composite numbers : `;   6: %=ifactor(%);  35: %=ifactor(%);  221: %=ifactor(%);

All natural numbers greater than one must be either prime or composite. Lets look at a larger set of integers - the numbers from 2 to 40

>	Integers := { seq(k, k = 2..40) } ;

Here are the prime numbers ² 40. Notice that NONE of them can be factored.

>	Primes := select( isprime, Integers ); for k from 1 to nops(Primes) do Primes[k] = ifactor(Primes[k]); od;

Here are the composite numbers ² 40. Note that each one of them can be factored.

>	Composites := remove( isprime, Integers);

>	for k from 1 to nops(Composites) do Composites[k] = ifactor(Composites[k]); od;

There is an ancient method of finding primes. You start with a table of natural numbers, then cast out multiples of primes. You leave the prime itself but remove all other multiples of that prime. You need only use primes less than or equal to the square root of the largest number in your table.

>	Sieve( 100);

2. The Sieve of Erasathones

There is an ancient method of finding a whole bunch of primes at once. You start with a table of natural numbers, then cast out multiples of primes. You leave the prime itself but remove all other multiples of that prime. You need only use primes less than or equal to the square root of the largest number in your table.

Lets look for all the primes up to about 40 .....
(Note : After you execute this next line, you'll need to manually scroll back up to look at each step).

>	Sieve(40);

You can see the distribution of primes from this diagram also. Notice the twin primes" : (3,5), (5,7),(11,13), (17,19), (29,31),(41,43).

How many primes are there less than 46? Look at the diagram and count how many numbers you see. You can also execute the next line.

>	numtheory[pi](46);

If you're ready for a little more excitement, hit the next command ... and hold on to your hat!

>	Sieve(199);

You can also see the distribution of primes less than 200 here! Notice that they seem to get more sparse as the numbers go up - in general. However, there are still "twin primes" such as 191 and 193 or 197 and 199 even close to the end!

Can you count how many primes are in this chart? You should counted this number :

>	numtheory[pi](200);

We can look at the distribution of primes another way too. If we count how many primes are less than or equal to a given number, then make a graph of that count we get this very interesting graph :

>	plot( numtheory[pi](floor(x)), x = 2..120);

There seems to be a pattern...which was Gauss's conjecture.

>	plot( {numtheory[pi](floor(x)), x/ln(x)}, x = 2..120);

3. Prime Factorization

 In this sense, composite numbers are composed  of smaller numbers, while primes are composed only of themselves.

>	2100 = `(30)(70)`;

However, each of the factors is composite, since they can be factored also.

>	30 = `(6)(5)`;  70 = `(7)(10)`;

And even some of these factors are composite and can be factored further.

>	6 = `(3)(2)`;  10 = `(5)(2)`;

Since each composite number can be factored, and possibly those factors factored further, eventually this process must end . If a factor is not prime, it can be factored further. The process ends when all of the factors of a number are prime.

>	2100 = `(30)(70) = (6)(5)(7)(10) = (3)(2)(5)(7)(5)(2)`;

If we collect all the like primes together, and sort them in ascending order, we get a unique prime factorizaton for every composite :

>	ifactor(2100);

We can easily get the prime factorizations for any composite using maple's integer factorization command, ifactor .

>	32: % = ifactor(%); 105: % = ifactor(%); 192: % = ifactor(%);

>	432: % = ifactor(%); 4332: % = ifactor(%); 4320: % = ifactor(%);

In practice, you can make a factor tree to complete factor a number. Then collect all of the 'leaves' of the tree which are the same prime number, to express in exponent form.

4. How many primes are there? and What is the Largest Prime Number?

There appear to be no shortage of prime numbers :

Here are the first 20 prime numbers ....

>	Primes(20);

and, here are the first 200 primes ...

>	Primes(100);

Here are the first 250 prime numbers ....

>	Primes(250);

We can find the 100th prime :

>	ithprime(100);

or the 200th, 300th, 400th ..., 2000th prime....

>	for k from 1 to 20 do print( cat(100*k, `th prime is `, ithprime(100*k))); od;

Is there a largest prime number? Suppose there was. Let say that p is the largest prime number. Then if we multiply all of the prime numbers together, and add one :
       N = 1 + p1*p2*p3*.....*pn

This new number is not divisible any prime number - which makes a contradiction. So there must be more primes than this supposedly "largest prime." This means there are actually an infinite number of primes.

Lets try this out on some prime. Suppose that 7 were the largest prime.

>	LargestPrime(7);

What if we thought that 17 might be the largest prime? or 37, or 53? (Warning : the bigger the number, the longer it will take. If it takes too long, hit the red stop-sign button on the Maple tool bar at top.)

>	LargestPrime(17);

>	LargestPrime(37);

>	LargestPrime(53);

Thus there are an infinite number of prime numbers, because no matter how large of a prime number we find, there is always a larger one! Mathematicians have found prime numbers with hundreds of thousands of digits. How many pages would it take to display a 300,000 digit number?

Lets do an experiment. Using my computer (a Powerbook with 1024x768 monitor), I find that this number is the largest power of ten I can display on my screen on a single line. Lets assume 52 lines per page.

>

1*10^105;

>	`number of digits` = 300000; `number of lines` = %/105; `number of pages` = %/50;

>	evalf((300000)/(105*52));

It would take most of 55 pages just to display one of these large prime numbers. Why don't we not try to find it and display, but go away with the comforting feeling that we now know such numbers exist..

5. Squarefree Numbers [optional]

 
Later when we want to simplify square roots, we will want to leave only square-free numbers inside the radical.

A square-free number is a number which contains no perfect squares. Look at the following numbers for a minute.

>	42 : % = ifactor(%); 330 : % = ifactor(%); 2431 : % = ifactor(%);

>	FindSquare(42); FindSquare(330); FindSquare(2431);

Each of these composite numbers is factorable - but each factor has an exponent of 1.

Obviously all primes are square free.

>	17 : % = ifactor(%); 41 : % = ifactor(%); 97 : % = ifactor(%);

>	FindSquare(17); FindSquare(41); FindSquare(97);

Now lets look at these numbers, which are NOT square free. Can you identify a perfect square "hiding" in each number?

>	12 : % = ifactor(%); 375 : % = ifactor(%); 72 : % = ifactor(%); 686 : % = ifactor(%);

>	FindSquare(12); FindSquare(375); FindSquare(72); FindSquare(686);

Notice that squarefree numbers are always products of distint primes - each with an exponent of 1! If a prime has a square, cube, 4th power, ...etc... then it contains at least a perfect square of that prime squared.

>	6 : % = ifactor(%);    FindSquare(%%);

>	12 : % = ifactor(%);    FindSquare(%%);

>	24 : % = ifactor(%);    FindSquare(%%);

>	48 : % = ifactor(%);    FindSquare(%%);

>	96: % = ifactor(%);    FindSquare(%%);

>	192 : % = ifactor(%);    FindSquare(%%);