Performing Gauss Elimination with MatLab

K =

   -0.2106    0.4656   -0.4531    0.7106
   -0.6018    0.2421   -0.8383    1.3634
    0.0773   -0.5600    0.4168   -0.2733
    0.7945    1.0603    1.5393    0.0098

I have the above matrix and I'd like to perform Gauss elimination on it with MatLab such that I am left with an upper triangular matrix. Please how can I proceed?

[Raphael Chinchilla]

use the function rref(K)

Answer 1

[József Szabó]

function x = solveGauss(A,b)
s = length(A);
for j = 1:(s-1)
    for i = s:-1:j+1
        m = A(i,j)/A(j,j);
        A(i,:) = A(i,:) - m*A(j,:);
        b(i) = b(i) - m*b(j);
    end
end 

x = zeros(s,1);
x(s) = b(s)/A(s,s);               
for i = s-1:-1:1                    
    sum = 0;
    for j = s:-1:i+1                
        sum = sum + A(i,j)*x(j);    
    end 
    x(i) = (b(i)- sum)/A(i,i);
end 

[Tyvaughn Holness]

Great work, thanks! I found a faster implementation that avoids the double for loop to reduce complexity and time.

https://www.youtube.com/watch?v=juDxSMpzv6c

[Ilyas Nhasse]

what if we got an A(i,i)=0

[Brinzan]

There are no input arguments, what do I do, qnd I dont even know where to put or how much to put like at the Matrix for A do I just write A=[00011110011],[0100001111],[0111110000] or what în the Code cuz it just not working.

[Walter Roberson]

@Brinzan

A = randi([-1 2], 5, 5)

A = 5×5

2 -1 0 0 2 0 -1 -1 2 1 1 2 2 2 0 -1 0 0 -1 2 -1 0 1 2 -1

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

b = randi([-2 2], 5, 1)

b = 5×1

1 1 -2 -1 -2

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

x = solveGauss(A, b)

x = 5×1

0.7059 0.0588 -1.3529 -0.0588 -0.1765

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

function x = solveGauss(A,b)
    s = length(A);
    for j = 1:(s-1)
        for i = s:-1:j+1
            m = A(i,j)/A(j,j);
            A(i,:) = A(i,:) - m*A(j,:);
            b(i) = b(i) - m*b(j);
        end
    end 
    
    x = zeros(s,1);
    x(s) = b(s)/A(s,s);               
    for i = s-1:-1:1                    
        sum = 0;
        for j = s:-1:i+1                
            sum = sum + A(i,j)*x(j);    
        end 
        x(i) = (b(i)- sum)/A(i,i);
    end 
end

Answer 2

[Richard Brown]

The function you want is LU

[L, U] = lu(K);

The upper triangular matrix resulting from Gaussian elimination with partial pivoting is U. L is a permuted lower triangular matrix. If you're using it to solve equations K*x = b, then you can do

x = U \ (L \ b);

or if you only have one right hand side, you can save a bit of effort and let MATLAB do it:

x = K \ b;

[Lukumon Kazeem]

Thank you Richard for your response. I have used this approach a no. of times ([L U]=lu(k)) and the results are always different from that in published literature. I suspect it's because it performs partial and not complete pivoting

[Richard Brown]

I wouldn't expect it would necessarily compare with published literature - what you get depends on the pivoting strategy (as you point out).

Complete pivoting is rarely used - it is pretty universally recognised that there is no practical advantage to using it over partial pivoting, and there is significantly more implementation overhead. So I would question whether results you've found in the literature use complete pivoting, unless it was a paper studying pivoting strategies.

What you might want is the LU factorisation with no pivoting. You can trick lu into providing this by using the sparse version of the algorithm with a pivoting threshold of zero:

[L, U] = lu(sparse(K),0);
% L = full(L); U = full(U); %optionally

Answer 3

[James Tursa]

You could start with this FEX submission:

http://www.mathworks.com/matlabcentral/fileexchange/13451-gauss-elimination-with-complete-pivoting

[Lukumon Kazeem]

James, thank you for your response. I have tried to apply your suggestion to the matrix I posed earlier but it came up with the below prompt. What do you think?

"??? Undefined function or method 'gecp' for input arguments of type 'double'".

[James Tursa]

You need to download the gecp function from the FEX link I posted above, and then put the file gecp.m somewhere on the MATLAB path.