Recursive Euclidean algorithm in Mathematica

5

2

Can anyone explain to me how do I use a recursion, if I don't know the limit? For example, I need the remainder $r$ of the Euclidean algorithm for $\gcd(a,b)$ which equals $0$. I figured out that the recursive formula I need is

r[n]=r[n-2]-Floor[r[n-2]/r[n-1] * r[n-1] 

r[1] = a - Floor[a/b]*  b;
r[2] = b - Floor[b/r1] r1;

zermelovac

Posted 2015-05-22T18:47:35.427

Reputation: 199

1

You may want to look into NestWhile. It will allow you to apply a function recursively until a user-specified condition is met. Here is its documentation page. Also, unless you are doing this for practice, Mathematica already has a built-in function to calculate the greatest common divisor, called GCD (docs here). You could use that one as a step towards solving your problem.

– MarcoB – 2015-05-22T18:58:48.880

Thank you for the advice, but I am still stuck, this is what I've tried:

`GreatestCommonDivisor[list_] := Module[{r1, r2, r3, r, a, b, n, k, p},
a = list[[1]]; 
b = list[[2]];
r[1] = a - Floor[a/b]*  b; 
r[2] = b - Floor[b/r1] r1;
r[3] = r1 - Floor[r1/r2] r2;
r[n_] := r[n - 2] - Floor[r[n - 2]/r[n - 1]]*r[n - 1]; 
Return[NestWhile[r, 4, # != 0 &, 1, \[Infinity], -1]]
]`
 – zermelovac  – 2015-05-22T19:15:39.340

If you are trying to implement the standard Euclid algorithm for gcd you should take a look at the recursive implementation here. Your version seems overly complicated.

– C. E. – 2015-05-22T19:53:33.210

Answers

5

As I mentioned in my comment, of course Mathematica has a built-in function to calculate the GCD, called GCD (docs).

My understanding, however, is that you are using the GCD as an example to learn how to apply a function recursively for a number of times that is not decided a priori, but that depends on the inputs and the path of the calculation.

The following is a sample implementation of Euclid's algorithm that strives to use Mathematica's functional style, rather than running through loops explicitly. Hopefully it might help you familiarize yourself with the style.

gcdlist[a_, b_] :=
 NestWhileList[
  {Last[#], (Mod[#1, #2] &) @@ #} &,
  {a, b},
  Last[#] != 0 &, 1
 ]

This function uses Mod to calculate the remainders, and NestWhileList (docs here) to apply a function recursively, as long as the condition is true. This function accomplishes the same purpose of NestWhile, but in addition shows all the intermediate results as well, so you can extract any remainders you need from the generated list:

gcdlist[25, 45]

(* Out: {{25, 45}, {45, 25}, {25, 20}, {20, 5}, {5, 0}} *)

The GCD value itself is the first element of the last element of this list:

gcd[a_, b_] := First@Last@gcdlist[a, b]

(*Or Using Part, gcd[a_, b_] := gcdlist[a, b][[-1, 1]]*)

gcd[4, 3]

(* Out: 1 *)

This is by no means the only way to implement the algorithm in Mathematica; it's just the first one that came to my mind, so feel free to modify it / improve it as you play around. For instance I suspect that the function FixedPoint might be helpful in this context as well.

MarcoB

Posted 2015-05-22T18:47:35.427

Reputation: 53 573

3Simply gcd[a_, b_] := If[b == 0, a, gcd[b, Mod[a, b]]] works too. – C. E. – 2015-05-22T19:50:15.573

@Pickett Oh, it definitely would, and it would be cleaner as well. But then again GCD[a,b] would be even better ;-) I suppose that I am trying to use this toy example to showcase how one might use NestWhile in a more complex problem that might not allow such a simple solution, for instance in the case that the termination condition might depend on more than the last function value, etc. Hopefully the OP might be able to use this as a springboard to dive into the documentation and examples therein. – MarcoB – 2015-05-22T19:55:56.833

@Pickett Thank you, I appreciate it! – MarcoB – 2015-05-22T20:30:24.713

1

@Pickett Related: (676) -- The "original" is gone; it was my reply to a do-my-homework-for-me question so I gave my answer in an obfuscated form: If[#2==0,#,#2~#0~Mod@##]&

– Mr.Wizard – 2015-05-26T13:22:34.770

@Mr.Wizard Nice! +1 – C. E. – 2015-05-26T14:05:23.157

@Mr.Wizard that's very neat! It should be included in the "Neat examples" section of the documentation of Slot to showcase its uses :-) – MarcoB – 2015-05-26T14:06:05.667

2

Method 1

gcdList[x_, y_] :=
  FixedPointList[
    # /. {a_, b_} /; b != 0 :> {b, Mod[a, b]} &, {x, y}]

gcdList[25,45]
{{25, 45}, {45, 25}, {25, 20}, {20, 5}, {5, 0}, {5, 0}}

Another solution

gcd[a_, b_] :=
 Module[{x, y},
  {x, y} = {a, b};
  While[y != 0,
   {x, y} = {y, Mod[x, y]};
  ];
  x
 ]

gcd[15,10]
 5

Method 2

gcdList2[a_, b_] :=
 Module[{x, y, res = {}},
  {x, y} = {a, b};
  While[y != 0,
   AppendTo[res, {x, y}];
   {x, y} = {y, Mod[x, y]}
  ];
  AppendTo[res, {x, 0}]
]

gcdList2[56, 21]
{{56, 21}, {21, 14}, {14, 7}, {7, 0}}

Method 3

With the help of Sow and Reap

gcdList3[a_, b_] :=
 Module[{x, y, res},
  {x, y} = {a, b};
  res =
   Reap[
    While[y != 0,
     Sow[{x, y}];
     {x, y} = {y, Mod[x, y]}
    ]];
  AppendTo[res[[2, 1]], {x, 0}]
 ]

gcdlist3[120, 75]
{{120, 75}, {75, 45}, {45, 30}, {30, 15}, {15, 0}}

xyz

Posted 2015-05-22T18:47:35.427

Reputation: 285