Math 2325
Project 3 (Chapter 3) - The Euclidean Algorithm

The Four Basic Programs

Euclid1

"Euclid1" computes the GCD of two positive integers. The output shows and counts all intermediate steps.

Euclid1[a_,b_]:=Module[{c=a,d=b,ctr=1},{c,d}=Sort[{c,d}];Print["Step ", ctr, ":\t ",c,"\t",d];While[c≠0,{c,d}={Mod[d,c],c};ctr++;Print["Step ", ctr, ": \t",c,"\t",d]]]

Here is an example: 
Let us compute the GCD of 15 and 24.

We begin in Step 1 with 15 and 24.

We divide the smaller number into the larger number; the remainder is 9. We discard the larger number and keep the remainder: Step 2 thus produces 9 and 15.

We divide the smaller number into the larger number; the remainder is 6. We discard the larger number and keep the remainder: Step 3  produces 6 and 9.

We divide the smaller number into the larger number; the remainder is 3. We discard the larger number and keep the remainder: Step 4  produces 3 and 6.

We divide the smaller number into the larger number; the remainder is 0. A remainder of 0 means the algorithm stops! We discard the larger number and keep the remainder: Step 3 thus produces 0 and 3. The non-zero number 3 is the GCD of 15 and 24. 

The program "Euclid1" tells the whole story:

Euclid1[24,15]

Here is the first example in the book:

Euclid1[6923, 5117]

Euclid2

"Euclid2" picks two random integers and computes their GCD. The rest works the same as "Euclid1". 

"Euclid2" takes one integer input: N. The two numbers whose GCD is being computed are chosen at random between 1 and N.

Euclid2[N_Integer]:=Module[{c=RandomInteger[{1,N}],d=RandomInteger[{1,N}],ctr=1},{c,d}=Sort[{c,d}];Print["Step ",ctr,":\t ",c,"\t",d];While[c≠0,{c,d}={Mod[d,c],c};ctr++;Print["Step ",ctr,": \t",c,"\t",d]]]

Execute the command below several times to check that different pairs of numbers between 1 and 1,000,000,000 are picked each time. Step 1 lists the two random integers chosen.

Euclid2[10^9]

Euclid3

"Euclid3" takes two positive integers as inputs, n and M.

The program picks M pairs of randomly chosen numbers between (10^n)and (10^(n+1))and computes how many steps are necessary in the Euclidean Algorithm to compute their GCD. 

The output lists how many times a certain number of steps was required to compute the GCD. The information is also listed as a bar graph. Make sure you understand the meaning of the output!

Also: While M can be rather large, keep n  reasonably small.

Euclid3[n_Integer,M_Integer]:=(Print["Frequency count for ",M," pairs of random integers chosen between ",(("10")^n)," and ",(("10")^(n+"1"))];store={};Do[Module[{c=RandomInteger[{(10^n),(10^(n+1))}],d=RandomInteger[{(10^n),(10^(n+1))}],ctr=1},{c,d}=Sort[{c,d}];While[c≠0,{c,d}=Sort[{Mod[d,c],c}];ctr++;];z=ctr;store=Append[store,z];];,{M}];
Print[TableForm[Sort[Tally[store]], TableHeadings->{None,{"Number","Frequency"}}]];
poss=Range[Max[store]];
Print[BarChart[(Count[store,#]& /@ poss), ChartLabels->poss,AxesLabel-> {"Number of steps","Frequency"}]]
)

Euclid3[3,50]

Euclid4

"Euclid4" is similar to "Euclid3".

"Euclid4" takes two positive integers as inputs, n and M.

The program picks M pairs of randomly chosen numbers between (10^n)and (10^(n+1))and computes their GCD. 

The output lists how many times a certain GCD was encountered. The information (for the 20 smallest GCDs only) is also listed as a bar graph. Make sure you understand the meaning of the output!

Also: While M can be rather large, keep n  reasonably small.

Euclid4[n_Integer,M_Integer]:=(Print["GCD distribution for ",M," pairs of random integers chosen between ",(("10")^n)," and ",(("10")^(n+"1"))];store={};Do[store=Append[store,GCD[RandomInteger[{(10^n),(10^(n+1))}],RandomInteger[{(10^n),(10^(n+1))}]]],{M}];
Print[TableForm[Sort[ Tally[store]], TableHeadings->{None,{"GCD","Frequency"}}]];
Print[BarChart[(Count[store,#]& /@ Range[20]), ChartLabels->Range[20],AxesLabel-> {"GCD","Frequency"}]])

Euclid4[8,500]

Additional Programs

Linear Diophantine Equations

The command below computes integer solutions to the Linear Diphantine Equation  a  x+b y=GCD(a,b).

Dioph[a_,b_]:=Module[{a2=Max[a,b],b2=Min[a,b],a1=Max[a,b],b1=Min[a,b],xold=1,yold=0,x=0,y=1,c=1},Print["The GCD of ",a2," and ",b2," is " ,GCD[a2,b2],".",
"\nWe want to find integer solutions x and y to the equation \n", a2 " * x + ", b2, " * y = ",GCD[a2,b2]];While[c≠0,q=IntegerPart[a2/b2];c=a2-q b2;If[IntegerQ[c],Print[a1 x+b1 y," = ",a1," * (",x,") + ",b1," * (",y,")"],];xnew=xold-q x;ynew=yold-q y;xold=x;yold=y;x=xnew;y=ynew;a2=b2;b2=c]]

Dioph[6923, 5117]

Fibonacci Numbers

A table of the first 40 Fibonacci numbers:

F[1]=1;F[2]=1;F[n_]:=F[n]=F[n-2]+F[n-1]
TableForm[Table[{n,F[n]},{n,1,40}],TableHeadings->{None,{"n","n-th Fibonacci number"}}](*,TableHeadings->{None,{"n", "n-th Fibonacci number"}}]*)

Euclid1[514229,832040]

Helmut Knaust, hknaust@utep.edu.	
September 29, 2003.
Modified for Mathematica 6
Art Duval, artduval@math.utep.edu.
October 7, 2008.
Slightly modified for Mathematica 7
Art Duval, artduval@math.utep.edu.
October 6, 2010.