Implementazione in Matlab per risolvere $LUx=Pb$

Il codice che risolve il sistema lineare $Ax=b$ tramite fattorizzazione con pivoting parziale:

function x=solvePALU(A,b)
%SOLVEPALU Risolve il sistema lineare Ax=b fattorizzando la matrice
%   A come LU applicando il pivoting ed infine risolvendo i sistemi
%                Ly=b
%                Ux=y
%
%   x=SOLVEPALU(A,b)
%
%   I parametri della funzione sono:
%       A -> la matrice dei coefficienti del sistema lineare
%       b -> il vettore dei termini noti
%
%   I valori di ritorno sono:
%       x -> il vettore soluzione del sistema lineare
%
%   See Also FATTPALU
  [A,p]=fattPALU(A);
  n=length(b);
  x1=b(p);
  for i=2:n
     for j=1:i-1
         x1(i)=x1(i)-A(i,j)*x1(j);
     end
  end
  x=solveUT(A,x1);