function [ndx,L,U] = plu(A)
% usage: [p,L,U] = plu(A)
% description: this function calculates a PLU
% factorization for the (assumed) invertible A.
% L and U are the unit lower and upper triangular
% matrices of the decomposition and p is an array
% of indices that corresponds to the permuations
% performed.

% programmer: your name!

m = size(A,1);
% create index and scales vectors
ndx = (1:m)';
s = zeros(m,1);
for ii = 1:m
  s(ii) = norm(A(ii,:),inf);
end
% main loop follows..overwrite A with L,U
% but remember to use indirect indexing
% to avoid actual row exchanges.


% finally, return the correct L,U
L = eye(m);
U = zeros(m,m);
% now fill in L,U with data from A