%% Demonstrate Eigenvalue Calculations via Secant Method
%
clear;

K=[12 -4 0; -4 6 -2; 0 -2 2];
M=[2 0 0; 0 2 0; 0 0 2];
A=K*inv(M);

%Let Matlab Calculate
display('Calculate eigenvectors (V) and eigenvalues (D) ...');
[V,D]=eig(A)
display('Calculate |A-lI| ... ');
det(A-D(1,1)*eye(3))

%Set iteration parameters
toler=0.0001;
maxiter=100;
n=size(A,1);

%Calculate by Secant Method Iteration
display('Secant Method Iteration');
lambda0=2;
lambdak=10;
k=0;
slope=0;
while ((max(abs((lambdak-lambda0)/lambda0))>toler) & (k<maxiter))
    fl0=det(A-lambda0*eye(n));
    flk=det(A-lambdak*eye(n));
    display(['k=',num2str(k,'%3.3i'),' | lambdak=',num2str(lambdak,'  %1.16e'),' | det()=',num2str(flk,'  %1.16e'),' | slope=',num2str(slope,'  %1.16e')]);
    k=k+1;
    slope=(flk-fl0)/(lambdak-lambda0);
    lambda0=lambdak;
    lambdak=lambdak-flk/slope;
end;

