how to compute orthonormal vectors via Lanczos process?

2 visualizaciones (últimos 30 días)
Omar B.
Omar B. el 27 de Ag. de 2019
Respondida: Maneet Kaur Bagga el 10 de Oct. de 2024
Hello,
How to initialize vector in Matlab? could you please check what I did ( I want to start with V_0=0 , beta (0)=0 , V(:,1) = V/norm(V,2))
Also, how can we creat a matrix collect all V's we have found ? ( i.e Q=[ V_1 V_2 ... V_j+1] )
My code:
function [] = lanczos(A, m)
[n,k] = size(A);
V(:,0)=0;
V(:,1) = ones(k,1);
V(:,1)=V(:,1)/norm(V(:,1),2);
beta(0)=0;
for j=1:m
w = A*V(:,j) - V(:,j-1)*beta(j-1);
alpha(j-1) = V(:,j)'.* w;
w = w - V(:,j)* alpha(j-1);
beta(j) = norm(w,2);
V(:,j+1) = w/beta(j);
end

Respuestas (1)

Maneet Kaur Bagga
Maneet Kaur Bagga el 10 de Oct. de 2024
Hi,
After going through the code given, I understand that the error you are encountering is due to the following reasons:
  • In MATLAB indexing starts from 1, therefore "V(:,0)" or "beta(0)" throws error. You can use "V(:,1)" and "beta(1)". Also, "j" should start from 2 (since your V(:,1) is the initial normalized vector).
[n, k] = size(A); % n = number of rows, k = number of columns (adjust if A is square)
V = zeros(k, m+1); % Allocate space for V vectors
V(:, 1) = ones(k, 1); % Initialize V1 as ones
V(:, 1) = V(:, 1) / norm(V(:, 1), 2); % Normalize V1
beta(1) = 0; % Initialize beta_0 = 0
  • The expression "alpha(j-1) = V(:,j)' .* w" performs the element wise multiplication where as for the Lancoz process you need to perform scalar multiplication of the two vectors to get the projection of "w" onto "v(:,j)".
alpha(j) = V(:, j)' * w; % scalar product for the two
  • To orthogonalize the vector "w" subtract the projection of "w" onto "V(:,j)" to remove the component of "w" that lies along "V". Therefore when you are updating "w" subtract the projection of "w" onto "V(:,j)" which is "w = w - V(:,j) * alpha(j)".
w = w - V(:, j-1) * beta(j); % Subtract previous beta term
Hope this resolves your query!

Categorías

Más información sobre Manage Products en Help Center y File Exchange.

Etiquetas

Community Treasure Hunt

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

Start Hunting!

Translated by