% script lusolve
eeps = 1; %well behaved
% eeps = 0; % 0 pivot
%eeps = .5*10^(-15); % num. 0 pivot
A = [1 1-eeps 3;2 2 20;3 6 4];
b = [5 - eeps; 6; 13];

[L,U] = GE(A);
y = LTriSol(L,b);
x = UTriSol(U,y);

pause

[L,U,piv] = GEpiv(A);
y = LTriSol(L,b(piv));
x = UTriSol(U,y);