% Script lutest.m
% Testing LU factorization based on scaled partial
% pivoting version gaussian elimination

A = [ 1,2,1,4; 2,0,4,3; 4,2,2,1; -3,1,3,2 ]
Acopy = A;

[A,r] = luFactor(A)

b = [13; 28; 20; 6]
x = luSolve(A,r,b)

% check

Acopy*x - b