Why doesn't the L in my partial pivot LU decomposition switch rows?

Question

Jessica Arroyo 2023-2-19

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1915310-why-doesn-t-the-l-in-my-partial-pivot-lu-decomposition-switch-rows

回答： Piyush Patil 2023-3-1

Everything else works well besides my L, please help.

A= [1 2 3 1; 4 5 6 2; 7 8 0 4; 0 1 3 1];
GE = GE_partial(A);
L = 4×4
     1     0     0     0
     0     1     0     0
     0     0     1     0
     0     0     0     1
L = 4×4
     1     0     0     0
     0     1     0     0
     0     0     1     0
     0     0     0     1
L = 4×4
     1     0     0     0
     0     1     0     0
     0     0     1     0
     0     0     0     1
%%
function[L,U,P] = GE_partial(A)
[n,n] = size(A); %nxn matrix A
L = eye(n,n); 
P = eye(n,n);
U = A;
for j =1:n-1 % looping over the columns of U
    %Finding the row width with the largest magnitude 
    [~, r] = max(abs(U(j:n,j)));
    % Calling the index of this row k
    k = r + j - 1;
    % Swap rows of j and k for U and P
    U([j,k],:) = U([k,j],:); %Swaping the first and the third row.
    P([j, k],:) = P([k, j],:); % swaping the first ad third row
    % Swap rows j and k in L for coulmns 1 to j-1
    for c = 1:j-1
        L([j,k],c) = L([k,j],c) 
    end
end
end

3 个评论
显示 1更早的评论隐藏 1更早的评论

Jessica Arroyo 2023-2-19

For U and P my rows are swapped. The output for L is the identity.

Torsten 2023-2-20

编辑：Torsten 2023-2-20

j=1, k=3 means: your loop to modify L is not entered. Do you understand ?

请先登录，再进行评论。

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

Answer 1

Piyush Patil 2023-3-1

1
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1915310-why-doesn-t-the-l-in-my-partial-pivot-lu-decomposition-switch-rows#answer_1182465

在 MATLAB Online 中打开

Hello Jessica,

As @Torsten correctly mentioned, you should check the values of variables j and k in the outer “for” loop.

You will notice that –

For j = 1, value of k = 3 but since value of j is 1, so it won’t enter the for loop that modifies the matrix L

For j = 2, value of k = 2. Therefore, both j and k are pointing to same row and hence swapping of rows won’t happen.

Similarly, for j = 3, value of k = 3 and hence swapping won’t take place.

Now, to add on to this discussion, since you are trying to perform LU Decomposition with partial pivoting, you need to modify matrix L and U by computing multipliers and eliminating the elements below pivot element.

You can do this by adding the following logic to the code –

% Compute multipliers and eliminate elements below pivot
for i=j+1:n
    L(i,j) = U(i,j)/U(j,j);
    U(i,:) = U(i,:) - L(i,j)*U(j,:);
end  

Also, you can replace the for loop that switches rows of matrix L by a single line –

L([j, k], 1:j-1) = L([k, j], 1:j-1);

Both will give the same result.

So, finally your entire function should look like –

function[L,U,P] = GE_partial(A)
    [n,n] = size(A); %nxn matrix A
    L = eye(n,n); 
    P = eye(n,n);
    U = A;
    for j =1:n-1 % looping over the columns of U
        %Finding the row width with the largest magnitude 
        [~, r] = max(abs(U(j:n,j)));
        % Calling the index of this row k
        k = r + j - 1;
        % Swap rows of j and k for U, P, and L
        U([j, k],:) = U([k, j],:);
        P([j, k],:) = P([k, j],:);
        L([j, k], 1:j-1) = L([k, j], 1:j-1);
        % Compute multipliers and eliminate elements below pivot
        for i=j+1:n
            L(i,j) = U(i,j)/U(j,j);
            U(i,:) = U(i,:) - L(i,j)*U(j,:);
        end
    end
end