Algorithms about eigenvalues
- Implement Matlab functions which can compute the corresponding eigenpair(s) via:
(a) Power iteration.
(b) Inverse iteration.
(c) Rayleign Quotient iteration.
(d) Simultaneous iteration.
(e) QR algorithm.
- Compare the differences and performances of these algorithms on a large-scaled problem.
Experimental purpose and requirement
Implement 5 funcions about eigonvalue and eigonvector.
Main equipment:
Hardware: MacBook Pro 2020, Catalina 11.15.7, 2 GHz Quad-Core Intel Core i5, 16 GB 3733 MHz LPDDR4X
MATLAB R2021b
Experimental procedure:
Write 5 function as required, use tic toc to record the time overhead of each calculating method.
Experimental result and analysis:
clc
sz = 50;
epsi = 1e-10;
maxx = 10000;
maxn = 5000;
A = randi([1,10],[sz,sz]);
A = A*A.';
% eig(A)
tic
for i = 1:maxn
PowerIteration(A,sz);
end
toc
tic
for i = 1:maxn
InverseIteration(A,sz);
end
toc
tic
for i = 1:maxn
RQIteration(A,sz);
end
toc
tic
for i = 1:maxn
SimulIteration(A);
end
toc
tic
for i = 1:maxn
QRALgorithm(A);
end
toc
% Power iteration
function [lamb] = PowerIteration(A,sz)
epsi = 1e-10;
maxx = 10000;
v = randi([1,10],[sz,1]);
v = v/norm(v);
for i = 1:maxx
w = A*v;
v = w/norm(w);
if mod(i,5)==0
lamb = v.' * A * v;
if norm(A*v - lamb.*v) < epsi
% disp(['loop ',num2str(i),' time(s)']);
break;
end
end
end
end
% Inverse iteration
function [lamb] = InverseIteration(A,sz)
epsi = 1e-10;
maxx = 10000;
v = randi([1,10],[sz,1]);
v = v/norm(v);
mu = 100;
XX = (A-mu.*eye([sz,sz]))^(-1);
for i = 1:maxx
w = XX*v;
v = w/norm(w);
if mod(i,5)==0
lamb = v.' * A * v;
if norm(A*v - lamb.*v) < epsi
% disp(['loop ',num2str(i),' time(s)']);
break;
end
end
end
end
% Rayleigh Quotient Iteration
function [lamb] = RQIteration(A,sz)
epsi = 1e-10;
maxx = 10000;
v = randi([1,10],[sz,1]);
lamb = v.' * A * v;
for i = 1:maxx
w = (A - lamb.* eye([sz,sz]))\v;
v = w/norm(w);
if mod(i,5)==0
lamb = v.' * A * v;
if norm(A*v - lamb.*v) < epsi
% disp(['loop ',num2str(i),' time(s)']);
break;
end
end
end
end
% Simultaneous Iteration
function [Q] = SimulIteration(A)
maxx = 1000;
epsi = 1e-10;
[Q,R] = qr(A);
for i = 1:maxx
Z = A*Q;
[Q,R] = qr(Z);
if mod(i,5)==0
v = Q(:,1);
lamb = v.' * A * v;
if norm(A*v - lamb.*v) < epsi
% disp(['loop ',num2str(i),' time(s)']);
break;
end
end
end
end
% QR Algorithm
function [A] = QRALgorithm(A)
maxx = 1000;
epsi = 1e-10;
for i = 1:maxx
[Q,R] = qr(A);
A = R*Q;
if mod(i,5)==0
v = Q(:,1);
lamb = v.' * A * v;
if norm(A*v - lamb.*v) < epsi
% disp(['loop ',num2str(i),' time(s)']);
break;
end
end
end
end
Elapsed time is 0.051586 seconds.
Elapsed time is 0.314361 seconds.
Elapsed time is 3.314072 seconds.
Elapsed time is 3.584557 seconds.
Elapsed time is 3.077736 seconds.
>>
本文来自博客园,作者:miyasaka,转载请注明原文链接:https://www.cnblogs.com/kion/p/15975485.html