%Conditionarea sistemelor de ecuatii liniare, sistem rezolvat cu metoda Gauss
clearvars
clc

%Date de intrare
%Matricea coeficientilor A
A= [ 10 2 -1
    -3 -6  2
     1  1  5]

%Vectorul termenilor liberi B
B=[1 2 3]'

%Vectorul perturbare R
R=[0.01 -0.01 0.01]'

B_pert=B+R
perturbare=norm(R)/norm(B);% ?R?/?B?, perturbarea trebuie sa ramana mica

% Determinarea numarului de conditie
disp('Numarul de conditie a matricei A este:')
numar_conditie=cond(A,"inf") % calculul numarului de conditie folosind norma liniilor
%numar_conditie=cond(A',"inf") % calculul numarului de conditie folosind norma coloanelor (inlocuind A cu A transpus)
fprintf('\n')

%Determinarea cu metoda Gauss a solutiilor sistemului neperturbat

%start proces eliminare
i = 1; % initializare
AB = [ A B ]; %matrice extinsa formata din A si B
[ nX mX ] = size( AB);   
while i <= nX 
   if AB(i,i) == 0 % se verifica daca termenii de pe diag principala au valoarea egala cu 0
            disp('Element egal cu zero pe diagonala principala a matricii coeficientilor') % se aloca mesajul in cazul in care avem termen egal cu zero pe diag principala
            return
   end
   % eliminarea elementelor de pe coloana i, de sub diagonala principala
   for k =  (i+1):nX 
    a=AB(k,i)/AB(i,i);
    AB(k,i:mX) = AB(k,i:mX) - a*AB(i,i:mX);  
   end
    i = i +1;
end  %final proces eliminare
    
U=AB(1:nX,1:nX); %matricea superior triunghiulara
B1=AB(:,mX); %termenii liberi dupa procesul de eliminare

%determinare solutie prin substitutie inapoi
x(nX)=AB(nX,mX)/AB(nX,nX);
for i=nX-1:-1:1
    S=B1(i);
    for j=i+1:nX
    S=S-U(i,j)*x(j);
    end
    x(i)=S/U(i,i);
end

%afisare rezultate sistem neperturbat

disp('Solutia sistemului neperturbat este:')
for k=1:nX 
fprintf ( '   x(%1.0f) = %2.10f  \n',k,x(k))
end
fprintf('\n')

C=x; % X

%Determinarea cu metoda Gauss a solutiilor sistemului perturbat

B=B+R; % B'

%start proces eliminare
i = 1; % initializare
AB = [ A B ]; %matrice extinsa formata din A si B
[ nX mX ] = size( AB);   
while i <= nX 
   if AB(i,i) == 0 % se verifica daca termenii de pe diag principala au valoarea egala cu 0
            disp('Element egal cu zero pe diagonala principala a matririi coeficientilor') % se aloca mesajul in cazul in care avem termen egal cu zero pe diag principala
            return
   end
   % eliminarea elementelor de pe coloana i, de sub diagonala principala
   for k =  (i+1):nX 
    a=AB(k,i)/AB(i,i);
    AB(k,i:mX) = AB(k,i:mX) - a*AB(i,i:mX);  
   end
    i = i +1;
end  %final proces eliminare
    
U=AB(1:nX,1:nX); %matricea superior triunghiulara
B1=AB(:,mX); %termenii liberi dupa procesul de eliminare

%determinare solutie prin substitutie inapoi
x(nX)=AB(nX,mX)/AB(nX,nX);
for i=nX-1:-1:1
    S=B1(i);
    for j=i+1:nX
    S=S-U(i,j)*x(j);
    end
    x(i)=S/U(i,i);
end

%afisare rezultate sistem perturbat
disp('Solutia sistemului perturbat este:')
for k=1:nX 
fprintf ( '   x(%1.0f) = %2.10f  \n',k,x(k))
end

%afisare rezultate sistem perturbat
fprintf('\n')
disp('Perturbarea privita in raport cu efectul pe care il provoaca asupra rezultatelor:')
D=x-C; %?E?
fprintf ('   ')
perturbare %?R?/?B?
fprintf ('   ')
efect=norm(D)/norm(C) %?E?/?X?
%daca efect este mai mic sau are cel mult acelasi ordin de marime cu perturbare atunci sistemul este bine conditionat