%function [x,a,b,er,o] = ludecompose(a,b)
function ludecompose(a,b)
%Function to use LU decomposition with partial pivoting
%using the Gauss Eliminatin technique
%to solve the linear system of equations [A]{x}={b}
%input [a] square coefficient matrix - size(n x n)
%input {b} right hand side vector - size (n x 1)
%output {x) - solution vector;
%output [a] - decomposed [a] matrix compactly storing [L] and [U]
%Function calls:    
%   from ludecompose: decompose & substitute
%   from decompose: ppivot

%n = length(b);\

a=[1 3 1;2 1 -4;1 4 2];
b=[1 2; 3 -3; -4 1]

[n,num_b_vectors]=size(b)

er = 0;
tol = 0.001;

%Decompose [A]=[L][U]
[a,n,tol,o,s,er] = decompose(a,n,tol,er);

%Solve for x with forward and backward substitution
%[L]{d}={b} and [U]{x}={d}

for i=1:num_b_vectors % loop through all of the b vectors
b_work=b(:,i);
if er~= -1      %if system is well-conditioned
    [a,o,n,b_work,x] = substitute(a,o,n,b_work);
    x
end
end

