Why is there an "inifite recursion"?

Hello,

I am currently trying to implement an more efficient algorithm to compute the greatest common divisor, which uses the hgcd. The Problem is that MATLAB goes into an infite recursion, although I don't know why.

function [M] = hgcd(p, q)
%hggT Ergibt den halbrößten gemeinsamen Teiler
%   ||p|| > ||q||

GradP = size(p)-1; %Hier wird Grad(p) und Grad(q) bestimmt
GradQ = size(q)-1;
GradP = GradP(1,2);
GradQ = GradQ(1,2);
m = ceil(GradP / 2); 
if GradQ < m %1. Fall: Der hgcd ist die Einheitsmatrix
  M = [1 0; 0 1];
else
  XhochM = [1, zeros(1, m)]; %Hier wird das Polynom x^m gebildet
  p_0 = deconv(p, XhochM); %Polynomdivision ist der deconv()-Befehl
  q_0 = deconv(q, XhochM);
  R = hgcd(p_0, q_0); %Rekursiver aufruf
  q = [zeros(1, GradP - GradQ), q]; %Fügt führende Nullen ein, damit p und q gleich groß sind
  R_inv = RegPolMatrixInv(R); %Invertiert R
  v = ConvPolMatrix3(R_inv, [p; q]); %Mulitipliziert R_inv und die Matrix [p \\ q]
  Groesse = size(v); %Im Folgenden wird v in p und q aufgeteilt
  Groesse = Groesse(1,2);
  for i = 1:Groesse
      pStrich(i) = v(1,i);
      qStrich(i) = v(2,i);
  end
  while qStrich(1)==0 %Löscht führende Nullen
     qStrich(1) = [];
  end
  GradQStrich = size(qStrich)-1; %Bestimmt den Grad
  GradQStrich = GradQStrich(1,2);
  if GradQStrich < m %2. Fall: Hier sind wir mit M = R zuende
      M = R;
    else
        [Q, s] = deconv(pStrich, qStrich); %Q ist Ergebnis der Polynomdivision, s ist der Rest
        r = pStrich; 
      l = size(r)-1; %l = grad(r)
      l = l(1,2);
      k = 2*m-l;
      XhochK = [1, zeros(k)]; %Erzeugt das Polynom x^k
      r_0 = deconv(r, XhochK);
      s_0 = deconv(s, XhochK);
      S = hgcd(r_0, s_0);
      Groesse2 = size(Q); %Dies ist notwendig, um <Q> zu erzeugen
      Groesse2 = Groesse2(1,2);
      Eins = [zeros(1, Groesse2)-1, 1]; %Das Polynom 1 mit führenden Nullen
      Nulll = zeros(1, Groesse2); %Das Nullpolynom mit führenden Nullen
      Q = [Q, Eins; Eins, Nulll]; %Erzeuge <Q>
      M = R*Q*S;
  end
end
end

Here is ConvPolMatrix3 a program that multiplies two Polynomial matrices, and RegPolMatrixInv a program that inverts any regular polynomial matrix. Their codes are the followings:

function [M] = ConvPolMatrix3(M, N)
%ConvPolMatrix3 N ist nur 2x1 Matrix.
% -----

h_1 = size(M);
h_2 = size(N);
h_1 = h_1(1,2);
h_2 = h_2(1,2);
h_1 = h_1 / 2;
for i = 1:h_1
 p_11(1,i) = M(1,i);
 p_12(1,i) = M(1,i+h_1);
 p_21(1,i) = M(2,i);
 p_22(1,i) = M(2,i+h_1);
end
for i = 1:h_2
 q_11(1,i) = N(1,i);
 q_21(1,i) = N(2,i);
end

r_111 = conv(p_11, q_11);
r_112 = conv(p_12, q_21);
r_211 = conv(p_21, q_11);
r_212 = conv(p_22, q_21);
g_111 = size(r_111)-1;
g_112 = size(r_112)-1;
g_211 = size(r_211)-1;
g_212 = size(r_212)-1;
g_111 = g_111(1,2);
g_112 = g_112(1,2);
g_211 = g_211(1,2);
g_212 = g_212(1,2);
Grad = max([g_111, g_112, g_211, g_212]);
r_111 = [zeros(1, Grad - g_111), r_111];
r_112 = [zeros(1, Grad - g_112), r_112];
r_211 = [zeros(1, Grad - g_211), r_211];
r_212 = [zeros(1, Grad - g_212), r_212];
r_11 = r_111 + r_112;
r_21 = r_211 + r_212;
i = 1;
while (r_11(i) == 0) && (r_21(i) == 0)
  r_11(i) = [];
  r_21(i) = [];
end
M = [r_11; r_21];
i = 0;
return
end

and, for RegPolMatrixInv,

function [M_inv] = RegPolMatrixInv(M)
%RegPolMatrixInv Invertiert eine reguläre Polynommatrix
%  Suche M^-1 nach dem Adjunktensatz

h = size(M);
h = h(1,2);
h = h/2;
p_11 = zeros(1, h);
p_12 = p_11;
p_21 = p_11;
p_22 = p_11;
for i = 1:h
 p_11(1,i) = M(1, i);
 p_12(1,i) = M(1, i+h);
 p_21(1,i) = M(2, i);
 p_22(1,i) = M(2, i+h);
end
J = ConvPolMatrix(p_11, p_12, p_21, p_22, p_22, -p_12, -p_21, p_11);
if J ~= [1 0; 0 1]
  det = -1;
else
  det = 1;
end
M_inv = det*[p_22, -p_12; -p_21, p_11];
return
end

Now, the two functions RegPolMatrixInv and ConvPolMatrix3 work just fine, but whenever I try to run the main program hgcd, and the second case is true, MATLAB goes into an infite recursion. Could you please tell me why and how to fix it? Also, sorry for the german comments and bad programming.

Thanks for helping me,

Florian

Answer 1

[Ameer Hamza]

Since you are calling the function hgcd() inside itself, there might be something going wrong in the passing of parameters

hgcd(p_0, q_0);

or

hgcd(r_0, s_0)

Try using MATLAB debugger to see if the recursive function is being called as expected and whether the termination condition will be fulfilled after a certain number of iteration.

[Ameer Hamza]

@Florian Albrecht answer moved here

Thanks for your answer. Unfortunately, I think everything is fine with the recursive call... When I enter p = [1 0 0] and q = [1 0], I get p_0 = [1 0] and q_0 = [0 0]. This should make for an even smaller recursion, but apperently this is not the case. Does one maybe need something like "return" to make the recursion work?

Kind regarts

Florian

[Ameer Hamza]

No, you don't need an explicit return. In the debugging mode, you can use step in button to step inside the recursive call and see what is happening inside the function.

[Florian Albrecht]

Thanks for the tips. I will look tommorow into it. Could you please tell me how to use the debugger correctly, since I never used it beforehand?

[Ameer Hamza]

See the link I shared in the answer. It will tell you how to put a breakpoint on a line. Next time when you call the function, the execution will stop at that line and then you can check the value of the variables which exist inside the function.

Therefore put a breakpoint on the line

R = hgcd(p_0, q_0); %Rekursiver aufruf

and when execution stop, click step in to go one level deeper into the function call stack. In this way, you can check what is going on inside that function i.e. what is the values of the variables. Therefore by going deeper into the recursive function call stack, you can see whether the termination condition is met or not.

[Florian Albrecht]

Thank you for answering! My problem was, that q_0 is, after a certain time, set to zero, and deg(0) should be -Inf. My code thought deg(0) = 0, so deg(p_0) = deg(q_0) instead of deg(p_0) > 2*deg(q_0), which is the abort condition.

Thanks a lot for your help

Florian

[Ameer Hamza]

You are welcome.