Re: Bairstow's Method

[superklutz]

Here is the code in case the link does not show up:

function [xre, xim, j, ear, eas] = Bairstow(a, deg, es, ri, si, maxit) %a =
matrix for coefficients
printf("Root evaluation using Bairstow's method.\n"); %Title
r = ri;
s = si;
b = zeros(1, deg+1);
c = zeros(1, deg);
xre = xim = zeros(deg, 1);
n = deg;
ier = 0;
ear = eas = 1;
j = 0;
ier = 0;
i = x = n-2;
do
  do
    j = j+1;
    b(1) = a(1);
    b(2) = a(2)+r.*b(n+1);
    c(1) = b(2);
    c(2) = b(n)+r.*c(n);
    do 
      b(x) = a(x) + r.*b(x-1) + s.*b(x-2);
      c(i) = b(x) + r.*c(i-1) + s.*c(i-2);
      i=i+1;
      x=i+1;
    until(i == n && x == n+1)
    det =  c(n-1).^2 - c(n-2).*c(n);
    if (det != 0)
      dr = (-b(n).*c(n-1) + b(n+1).*c(n-2))./det;
      ds = (-b(n+1).*c(n-1) + b(n).*c(n))./det;
      r = r + dr;
      s = s + ds;
    else
      r = r+1;
      s = s+1;
      j = 0;
    endif
  until (ear<es && eas<es || j >= maxit) %Stopping criterion
[r1, i1, r2, i2] = quadform (r, s)
xre(n) = r1;
xim(n) = i1;
xre(n-1) = r2;
xim(n-1) = i2;
n = n-2;
  do
    a(x) = b(x-2);
    x = x-1;
  until (x == 1)
until (n < 3 || j >= maxit)
if (j < maxit)
  if (n == 2)
  r = -a(deg)./a(deg-1);
  s = -a(deg+1)./a(deg-1);
  [r1, i1, r2, i2] = quadform (r, s)
    xre(n) = r1;
    xim(n) = i1;
    xre(n-1) = r2;
    xim(n-1) = i2;
  else
  xre(n) = -a(deg+1)./a(deg);
  xim(n) = 0;
  endif
else
  ier = 1;
endif
printf("After %d iterations, the computed roots are as follows: \n", j);
%Prints the root and the error
xre
xim
clear
endfunction
function [r1, i1, r2, i2] = quadform (r, s)
dis = r^2 + 4*s;
if (dis>0)
  r1 = (r + sqrt(dis))/2;
  r2 = (r - sqrt(dis))/2;
  i1 = i2 = 0;
else
  r2 = r1 = r/2;
  i1 = sqrt(abs(dis))/2;
  i2 = -i1;
endif
endfunction

This is the error:

>> Bairstow([1, -3.5, 2.75, 2.125, -3.875, 1.25], 5, 0.00005, -1, -1, 100)
Root evaluation using Bairstow's method.
error: Bairstow: A(I): index out of bounds; value 7 out of bound 6
error: called from
    Bairstow at line 22 column 12
>>



--
View this message in context: 
http://octave.1599824.n4.nabble.com/Bairstow-s-Method-tp4672471p4672476.html
Sent from the Octave - General mailing list archive at Nabble.com.