Page 118 - Compact Numerical Methods For Computers

The algebraic eigenvalue problem 107

Algorithm 10. Inverse iteration via Gauss elimination (cont.)
while (msame<nRow) and (itcount<itlimit) do
begin {STEP 11 -- perform the back-substitution first}
itcount:=itcount+1; {to count the iterations}
m:=nRow; s:=X[nRow];
X[nRow]:=Y[nRow]; {save last trial vector -- zeros on iteration 1}
if abs(A[nRow,nRow])<tol then Y[nRow]:=s/tol
else Y[nRow]:=s/A[nRow,nRow];
t:=abs(Y[nRow]);{ to set the first trial value for vector Y}
for i:=(nRow-1) downto 1 do {STEP 12}
begin {back-substition}
s:=X[i]; X[i]:=Y[i];
for j:=(i+1) to nRow do
begin
s:=s-A[i,j]*Y[j];
end;
if abs(A[i,i])<tol then Y[i]:=s/tol else Y[i]:=s/A[i,i];
if abs(Y[i])>t then
begin
m:=i; t:=abs(Y[i]);
end; {to update new norm and its position}
end; {loop on i}
ev:=shift+X[m]/Y[m]; {current eigenvalue approximation -- STEP 13}
writeln(‘Iteration’ ,itcount,’ approx. ev=’,ev);
{Normalisation and convergence tests -- STEP 14}
t:=Y[m]; msame:=0;
for i:=1 to nRow do
begin
Y[i]:=Y[i]/t;
if reltest+Y[i] = reltest+X[i] then msame:=msame+1;
{This test is designed to be machine independent. Mixed mode
arithmetic is avoided by the use of the constant reltest. The
variable msame is used in place of m to avoid confusion during the
scope of the ‘while’ loop.}
end; {loop on i}
{STEP 15 -- now part of while loop control}
if msame<nRow then
begin {STEP 16 -- multiplication of vector by transformed B matrix}
for i:=1 to nRow do
begin
s:=0.0;
for j:=1 to nRow do s:=s+A[i,j+nRow]*Y[j];
X[i]:=s;
end; {loop on i}
end; {if msame < nRow}
end; {while loop -- STEP 17 is now part of the while loop control}
if itcount>=itlimit then itcount:=-itcount; {set negative to
indicate failure to converge within iteration limit}
shift:=ev; {to pass eigenvalue to calling program}
end; {alg10.pas == gii}