%Metoda Cholesky pentru rezolvarea sistemelor de ecuatii neliniare
clearvars
clc

%Date  de intrare
%Matricea coeficientilor A. Trebuie sa fie simetrica si pozitiv definita
A= [ 10 -1 2  4 1
     -1  4 2  2 2
      2  2 10 3 3
      4  2 3  4 4
      1  2 3  4 80];
 
 %Vectorul termenilor liberi B
 B=[2 15 10 -8 -12];
 
 %Se verifica daca matricea A este pozitiv definita
 [L, flag]=chol(A);
 if flag~=0 
     disp('Matricea A NU este pozitiv definita, factorizarea NU este posibila');
 return
 end
 
 %Se verifica daca matricea A este pozitiv definita
 test_sim = issymmetric(A);
 if test_sim~=1 
     disp('Matricea A NU este simetrica, factorizarea NU este posibila');
 return
 end
 
 %factorizare: transformare A=L*LT
  L=chol(A,'lower'); %matricea inferior triunghiulara
  LT=L';
  
 %determinare solutii
 n=length(A); %dimensiunea sistemului de ecuatii
       %dupa factorizare sistemul de ecuatii este de forma L*LT*X=B
 
 %determinare solutiilor se face in doua etape
 
 %ETAPA 1:substitutie inainte; practic se rezolva sistemul L*y=B 
            %unde cu y s-a notat LT*X
 y(1)=B(1)/L(1,1);
 for i= 2:n 
    S=B(i);
    for j=1:i-1
      S=S-L(i,j)*y(j);
    end
    y(i)=S/L(i,i);
 end
 
 %ETAPA 2:substitutie inapoi; practic se rezolva sistemul LT*X=Y
LT=L'; 
x(n)=y(n)/LT(n,n);
for i=n-1:-1:1
    S=y(i);
    for j=i+1:n
    S=S-LT(i,j)*x(j);
    end
    x(i)=S/LT(i,i);
end

%afisare solutie
disp('Solutia sistemului de ecuatii este:')
for k=1:n 
fprintf ( '   x(%1.0f) = %2.10f  \n',k,x(k))
end

%proba solutiei; practic se face operatia A*X-B
disp ('   ');
disp('Proba solutiei A*X-B=0:') 
%fprintf ( ' %2.20f  \n', A*x'-B');    

for k1=1:n 
fprintf ( '  Ecuatia %2.0f: %2.20f  \n', k1, A(k1,:)*x'-B(k1));
end