Problem with computing inverse using LU

Hi! I seem to have a problem getting the exact inverse of a matrix using LU. This is the code I made, I already have the code for formulating the L and U, this is just a the inverse part for testing.

l = [  2            0            0            0
      -1          1.5            0            0
       0           -1          4/3            0
       0            0           -1         1.25];

u = [1         -0.5            0            0
     0            1         -2/3            0
     0            0            1        -0.75
     0            0            0            1];

n = length(a);
x = zeros(n,1);
c = zeros(n,1);
d = zeros(n,1);
inverse = zeros(n);
c(1) = 1;

d(1) = c(1) / l(1,1);
for k=1:n
    for i=2:n
        sum = 0;
        for j=1:i-1
            sum = sum + l(i,j) * d(j);
        end
        d(i) = (c(i) - sum) / l(i,i);
    end
    x(n) = d(n) / u(n,n);
    for i=n-1:-1:1
        sum = 0;
        for j=i+1:n
            sum = sum + u(i,j) * x(j);
        end
        x(i) = [d(i) - sum] / u(i,i);
    end

      c(k)=0;
      c(k+1)=1;
      inverse(:,k) = x; 
  end

This is the result of my code:

inverse =

            0.8          1.4          1.2            1
            0.6          1.8          1.4            1
            0.4          1.2          1.6            1
            0.2          0.6          0.8            1

while the true inverse is

            0.8          0.6          0.4          0.2
            0.6          1.2          0.8          0.4
            0.4          0.8          1.2          0.6
            0.2          0.4          0.6          0.8

I tested it and I think that the problem may be in the outermost for loop. I just don't know specifically. Thanks in advance!

[Jutaporn Artniyom]

What is the value of c represent for, and if it's possible to explain how this script work thanks a lot

[Mariana Pinheiro]

Can you provide the code, please?

Answer 1

[Yucheng Ma]

It is my understanding that you would like to implement a C-style matrix inverse procedure using LU decomposition in MATLAB. The code above has a minor mistake in computing the inverse of the L matrix, i.e. "d(1)" is initialized but never updated. I rewrote part of the code and pointed out the difference in the comments. Please refer to the attached file "invLU.m".

In MATLAB, you can use the "inv" function to calculate the inverse of a matrix. You can also use the "mldivide" operator("\") to solve systems of linear equations. The "\" operator is more efficient than explicitly calculating the inverse of a matrix, and can handle singular matrices and sparse matrices.

[Marc Edwin Montilla]

Thanks for pointing out my error! It is working fine now. We are tasked to solve for the inverse of a matrix by only using the LU decomposition specifically so I guess the inv function is just for checking. I was just wondering when you said that I am implementing a C-style procedure, Is there any other "style"? Anyways, thanks for the suggestions!

[Mariana Pinheiro]

Can you provide the code, please?