%-----------------------------------------------------------------------
%             -------------                             --------
%             function file                             LUPivo.m
%             -------------                             --------
%-----------------------------------------------------------------------
%
%               ------------------------------------------
%               Decomposicao  LU  com Pivoteamento Parcial
%---------------------------------------------------------
%
%
%       Chamada do procedimento LUPivo
%-------------------------------------
%       [LU,P] = LUPivo(A)
%
%
%       Variavel de Entrada
%--------------------------
%        A : Matriz de ordem  n
%
%       Variaveis de Saida
%-------------------------
%       LU : Matriz contendo a Decomposicao  LU  da matriz  A.          
%
%        P : Vetor de dimensao  n  contendo a matriz de Permutacao
%
%
%
%
        function [LU,P] = LUPivo(A)
%----------------------------------
%
%
%       Precisao para o teste de Singularidade
        epsilon = 1.0e-10;   
%
%       Ordem da Matriz A
        [n,n] = size(A);    
%
%       Atribuicao da Matriz  A  na matriz  LU
        LU = A;          
%
%       Construcao inicial da Matriz de Permutacao
        P = 1:n;
%
%
%       Inicio da Decomposicao  LU
%---------------------------------
        for j = 1:(n-1)
%
%       Escolha do Pivo
%----------------------
        [maxvalor,i] = max(abs( LU(j:n,j) ));
        i = i + j - 1;
%
%
%       Verifica se a Matriz e' Singular
%---------------------------------------
        if ( maxvalor <= epsilon )
        disp('  '), disp('    ')
        error('      Matriz  Singular')
        end        
%
%
%       Permutacao das Linhas
%----------------------------
        LU([j i],:) = LU([i j],:);
%
        P([j i]) = P([i j]);               
%
%
%       Construcao da Matriz  L
%------------------------------
        LU((j+1):n,j)  =  LU((j+1):n,j) / LU(j,j);
%
%       Construcao da Matriz  U
%------------------------------
        soma = LU((j+1):n,(j+1):n) - LU((j+1):n,j)*LU(j,(j+1):n);
        LU((j+1):n,(j+1):n) = soma;
%
%
        end
%
%-----------------------------------------------------------------------
%-----------------------------------------------------------------------
%--------------------------------------------------
%       Petronio Pulino
%       Departamento de Matemática Aplicada - D.M.A
%       IMECC - UNICAMP
%       e-mail: pulino@ime.unicamp.br
%       http://www.ime.unicamp.br/~pulino/MS411/
%--------------------------------------------------
