Arnoldi method to find eigenvalues

21 次查看（过去 30 天）
Haya M
Haya M 2020-4-28
回答： Pavl M. 2024-10-11
I'm trying to implement the Arnoldi method with the inverse power method to find eigenvalues of large pencil matrix.
I have followed the practical implementation in Saad's book and I started with a small matrix to check if the code work well. As I understand H is a square matrix and has size of the number of the iterations but the resulted H is of size 3x2 and V is 4x3. I'm not sure if this is correct and I do'nt know how I can find the eigenvalues of H and the corresponding eigenvectors. Here is my attempt, and I really appreciate any help..
A=[1 0 0 0;0 17 0 0;0 0 -10 0 ;0 0 0 5]
n=length(A)
B=eye(n)
a=0 %interval [a,b]
b=10
m=2; %iterations
x=ones(n,1); %constant vector
shift=0.5;
V = zeros(n,m); %orthonormal basis of Krylov space
V(:,1) = x/norm(x);
H = zeros(n,m); %upper Hessenberg matrix
C=A-shift*B;
[L, U] = lu(C, 'vector');
for j=1:m
w=U\(L\(B*V(:,j)));
for i=1:j %Gram-schmidt
H(i,j) = V(:,i)'*w;
w = w - V(:,i)*H(i,j);
end
H(j+1,j) = norm(w,2);
V(:,j+1) = w/H(j+1,j);
end

请先登录，再回答此问题。

回答（1 个）

Pavl M.
Pavl M. 2024-10-11
Ran in:
Hi, happy to help:
Use next code:
clc
clear all
close all
%non-square matrix:
%A=[1 0 0 0;0 17 0 0;0 0 -10 0 ;0 0 0 5;0 0 0 0]
A=[1 0 0 0;0 17 0 0;0 0 -10 0 ;0 0 0 5]
A = 4×4
1 0 0 0 0 17 0 0 0 0 -10 0 0 0 0 5
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
%DMD code:
%version 1:
n1=size(A,1);
n2=size(A,2);
B=eye(n1,n2);
a=0 %interval [a,b]
a = 0
b=10
b = 10
m=n2; %iterations
x=ones(n2,1); %constant vector
shift=0.5;
V = zeros(n2,m); %orthonormal basis of Krylov space
V(:,1) = x/norm(x);
H = zeros(min(m+1,m),n1); %upper Hessenberg matrix
W = zeros(n1,m);
Z = zeros(n1,m);
C=A-shift*B;
[L, U] = lu(C, 'vector');
for j=1:m
w=U\(L\(B*V(:,j)))
z = A*V(:,j)
W(:,j) = w;
Z(:,j) = z;
for i=1:j %Gram-schmidt
H(i,j) = V(:,i)'*w;
w = w - V(:,i)*H(i,j);
end
if j < n1
H(j+1,j) = norm(w,2);
if H(j+1,j) == 0, return, end
V(:,j+1) = w/H(j+1,j);
end
end
w = 4×1
1.0000 0.0303 -0.0476 0.1111
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
z = 4×1
0.5000 8.5000 -5.0000 2.5000
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
w = 4×1
1.7168 -0.0174 0.0361 -0.0426
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
z = 4×1
0.8584 -4.8835 3.7932 -0.9590
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
w = 4×1
0.2294 -0.0042 -0.0645 -0.1607
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
z = 4×1
0.1147 -1.1711 -6.7746 -3.6164
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
w = 4×1
0.0110 0.0493 0.0365 -0.0969
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
z = 4×1
0.0055 13.8395 3.8362 -2.1799
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
H
H = 4×4
0.5469 0.8464 -0.0000 0.0000 0.8464 1.4731 0.2534 0.0000 0 0.2534 0.0991 0.0927 0 0 0.0927 0.0684
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
V
V = 4×4
0.5000 0.8584 0.1147 0.0055 0.5000 -0.2873 -0.0689 0.8141 0.5000 -0.3793 0.6775 -0.3836 0.5000 -0.1918 -0.7233 -0.4360
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
W
W = 4×4
1.0000 1.7168 0.2294 0.0110 0.0303 -0.0174 -0.0042 0.0493 -0.0476 0.0361 -0.0645 0.0365 0.1111 -0.0426 -0.1607 -0.0969
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
Z
Z = 4×4
0.5000 0.8584 0.1147 0.0055 8.5000 -4.8835 -1.1711 13.8395 -5.0000 3.7932 -6.7746 3.8362 2.5000 -0.9590 -3.6164 -2.1799
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
%find dynamical modes:
%vector of mode amplitudes and ? estimation:
v1 = H*A
v1 = 4×4
0.5469 14.3892 0.0000 0.0000 0.8464 25.0426 -2.5343 0.0000 0 4.3084 -0.9915 0.4634 0 0 -0.9269 0.3422
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
v2 = A*diag(Z)
v2 = 4×1
0.5000 -83.0188 67.7460 -10.8993
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
v23 = diag(z)*A
v23 = 4×4
0.0055 0 0 0 0 235.2707 0 0 0 0 -38.3617 0 0 0 0 -10.8993
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
v11 = W*A
v11 = 4×4
1.0000 29.1848 -2.2943 0.0550 0.0303 -0.2960 0.0417 0.2467 -0.0476 0.6141 0.6452 0.1827 0.1111 -0.7245 1.6073 -0.4844
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
v22 = Z*A
v22 = 4×4
0.5000 14.5924 -1.1471 0.0275 8.5000 -83.0188 11.7107 69.1973 -5.0000 64.4848 67.7460 19.1809 2.5000 -16.3023 36.1643 -10.8993
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
display('Accuracy measure 1:')
Accuracy measure 1:
acc = mean(mean(abs(V-Z)))
acc = 3.4670

类别

在 Help Center 和 File Exchange 中查找有关 Programming Utilities 的更多信息

标签

产品

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by