Curs 13 Rezolvarea sistemelor de ecuaţii liniare
Rezolvarea sistemelor de ecuaţii liniare
Metoda de eliminare a lui Gauss
Metoda lui Gauss este una din metode tradiţionale directe de rezolvare a sistemelor de
ecuaţii liniare. Ideea de bază a metodei constă în aducerea sistemului de ecuaţii prin
transformări elementare la o formă echivalentă, având matrice superior sau inferior
triunghiulară, urmată de rezolvarea sistemului rezultat prin procedee recurente specifice,
foarte eficiente.
Transformarea sistemului iniţial într-un sistem de formă triunghiulară se realizează cu
ajutorul a trei operaţii elementare:
a) interschimbarea a două ecuaţii între ele;
b) înmulţirea unei ecuaţii cu o constantă nenulă;
c) scăderea unei ecuaţii din alta şi înlocuirea celei de-a doua ecuaţii cu rezultatul
scăderii.
Transformarea sistemului este echivalentă cu eliminarea succesivă a necunoscutelor din
ecuaţii, prin această operaţie ajungându-se la o matrice triunghiulară a sistemului (faza
eliminării). Rezolvarea sistemului cu matrice triunghiulară constă în determinarea
necunoscutelor şi substituţia lor în ecuaţiile sistemului în ordine inversă, fiind denumită din
acest motiv faza substituţiei inverse.
Fie pentru exemplificare un sistem de trei ecuaţii cu trei necunoscute:
(1)
sau sub formă matricială
respectiv:
A · x =b
cu A=[a ij] – matricea sistemului,
x=[xi] – matricea coloană a necunoscutelor
b=[bi] – matricea coloană a termenilor liberi.
Pivotarea parţială
Pentru majoritatea metodelor de rezolvare a sistemelor de ecuaţii liniare - fie ele
directe sau iterative - se ajunge ca la un moment dat sa fie necesară împărţirea la un
element diagonal din matricea A, ajj. Acest element poartă numele generic de pivot. Desigur,
dacă acest element este nul metoda respectivă eşuează, deoarece - din punct de vedere
numeric - operaţia de împărţire la 0 este imposibilă şi, dacă nu se iau masuri speciale, orice
program de calcul se întrerupe din execuţie in urma producerii unei erori de tipul "Run time
error" şi un eventual mesaj de eroare de genul "Floating point overflow".
Situaţia extremă în care într-o matrice A poate să apară un element diagonal nul este
cea în care matricea respectivă este singulară. Pe de alta parte, pivotul se poate anula fară
ca matricea A să fie singulară. Mai mult decât atât, pivotul poate fi nenul, dar cu o valoare
foarte mică, astfel încât împărţirea la el să conducă la producerea unor erori de rotunjire,
care prin acumulare pot denatura rezultatul.
1
Curs 13 Rezolvarea sistemelor de ecuaţii liniare
Evitarea unor asemenea situaţii se poate face prin adoptarea unei masuri care să
permită aducerea pe diagonală a unui element suficient de mare. O asemenea tehnică este
aceea a pivotării care constă în schimbarea între ele a două linii sau a două coloane, astfel
încât noul pivot sa aibă o valoare absolută cât mai mare posibil. Căutarea noului pivot se
face pe coloana curentă j*, pe liniile situate sub linia j*, inclusiv aceasta (pivotarea
partiala) sau pe liniile şi coloanele situate sub linia j* şi la dreapta coloanei j*, inclusiv
acestea (pivotarea completă).
In cazul pivotării parţiale, la un pas j*, se spune ca se executa pivotarea partială pe
coloana j*. Prin această tehnică se caută pe coloana j* elementul subdiagonal maxim în
valoare absolută şi se schimbă linia j* cu linia in care apare acel element. Deoarece
schimbarea a două linii în matricea A înseamnă de fapt schimbarea ordinii a două ecuaţii din
sistem, se impune ca - simultan - să se realizeze şi schimbarea termenilor liberi
corespunzatori din vectorul b şi de asemenea, schimbarea semnului determinantului matricii.
În cazul concret al unui sistem de trei ecuaţii cu trei necunoscute, la primul pas al
etapei eliminării urmărim eliminarea necunoscutei x1 din toate ecuaţiile sistemului, cu
excepţia primei ecuaţii. Pentru aceasta, împărţim mai întâi prima linie la elementul pivot a11,
presupus nenul, adică a11 ≠0 (dacă nu este îndeplinită această condiţie se reordonează
ecuaţiile sistemului):
Scădem apoi prima ecuaţie înmulţită cu primul coeficient al celei de-a doua ecuaţii
din a doua ecuaţie şi, respectiv, prima ecuaţie înmulţită cu primul coeficient al celei de-a
treia ecuaţii din cea de-a treia ecuaţie. Obţinem astfel sistemul:
unde:
(2)
Aşadar, primul pas al metodei eliminării a lui Gauss conduce la ecuaţia matricială:
(3)
La următorul pas, eliminăm necunoscuta x2 din ultima ecuaţie. Pentru aceasta, împărţim mai
întâi a doua ecuaţie la elementul pivot , diferit de zero (dacă nu este aşa, atunci realizăm
interschimbarea ecuaţiilor a doua şi a treia) şi apoi scădem linia obţinută, înmulţită cu ,
din ecuaţia a treia. Acest al doilea pas conduce la ecuaţia matricială:
2
Curs 13 Rezolvarea sistemelor de ecuaţii liniare
(4)
unde:
(5)
Faza eliminării se încheie împărţind cea de a treia ecuaţie la elementul pivot ,
care, pentru un sistem cu matrice nesingulară (matrice pătratică, cu determinant nenul),
trebuie să fie diferit de zero. Rezultă după acest pas sistemul:
(6)
cu
(6)
sau matricial A(3) · x =b(3).
Se observă că matricea A(3) este superior triunghiulară, iar sistemul obţinut este
echivalent cu cel iniţial A · x =b, adică are soluţia (x1, x2, x3).
Faza substituţiei inverse (mersul înapoi) presupune parcurgerea în sens invers a ecuaţiilor
sistemului cu matrice triunghiulară (6), rezultat în faza eliminării şi stabilirea soluţiei
sistemului potrivit unui calcul recursiv:
(7)
Aşadar, determinarea soluţiilor se face recursiv, de la indici mari spre indici mici,
fiecare nouă soluţie depinzând în mod explicit numai de soluţiile determinate anterior.
Determinantul matricii sistemului
Metoda de eliminare Gauss permite şi calcularea determinantului matricii sistemului.
Se observă că, matricea A(3) a sistemului (6) fiind triunghiulară, are determinantul egal cu
produsul elementelor diagonale, adică:
det A(3)=1
Având în vedere, însă, că prin împărţirea liniilor matricii sistemului la elementele pivot rezultă
o matrice a cărui determinant este egal cu determinantul matricii iniţiale împărţit la produsul
elementelor pivot, rezultă:
det A
det A (3) = =1 adica : det A = a11a(221) a(332)
a11a(221) a(332)
3
Curs 13 Rezolvarea sistemelor de ecuaţii liniare
Metoda descrisă mai sus poate fi generalizată şi aplicată la rezolvarea unui sistem cu
un număr oarecare de n ecuaţii liniare cu n necunoscute de forma (8):
(8)
sau sub formă matriceală:
La fiecare pas k al fazei eliminării se efectuează împărţirea liniei pivot k la elementul
diagonal . Deoarece este posibil ca un element pivot să fie egal cu zero (imposibilitate
de continuare a fazei eliminării !) se procedează în practică astfel:
a) se determină elementul maxim al matricii A(k–1) situat pe coloana k şi pe liniile l ≥ k
(pe diagonala principală şi sub ea), adică elementul alk;
b) se aduce elementul maxim alk găsit pe diagonala principală, interschimbând între ele
liniile l şi k.
Relaţiile de calcul a elementelor matricii sistemului şi a matricii termenilor liberi sunt:
akij = akij −1 − mki −1 ⋅ akkj−1
cu
akik−1
mki −1 =
akkk−1
respectiv :
bki = bki −1 − mki −1 ⋅ bkk −1
cu: i,j=k+1, ...,n şi k=1,n-1
Inversa unei matrici
Fie o matrice pătratică A=[aij]. Pentru determinarea inversei acestei matrici se foloseşte
relaţia A·A-1=I unde I este matricea unitate.
Fie matricea inversă
⎡ x11 x12 ... x1n ⎤
⎢ ⎥
x x 22 ... x 2n ⎥
A −1 = ⎢ 21
⎢ ... ... ... ... ⎥
⎢ ⎥
⎣ xn1 xn2 ... xnn ⎦
Cum A·A-1 trebuie să fie I, obţinem:
n
∑a
k =1
ik xkj = δij i, j = 1,...,n
Aşadar, pentru determinarea inversei matricii A trebuie rezolvate n sisteme de ecuaţii
care au toate aceeaşi matrice A, fiind diferite doar matricile termenilor liberi. Cele n2
elemente ale matricii A-1 pot fi determinate folosind metoda Gauss.
4
Curs 13 Rezolvarea sistemelor de ecuaţii liniare
Program Metoda Gauss
#include<stdio.h>
#include<stdlib.h>
#include<math.h>
#include<conio.h>
#define dm 10
// Constanta dm reprezinta dimensiunea maxima a matricii sistemului
double a[dm][dm],ain[dm][dm],inv[dm][dm],prod[dm][dm]; //a este matricea
sistemului
double b[dm],x[dm],det=1.0,dtm; //b este matricea termenilor liberi
int n,opt; //n este dimensiunea actuala a matricii sistemului
FILE *f; //descriptorul fisierului in care se vor scrie toate rezultatele
void stop()
{
printf("\n\n====================================");
printf("\nApasa orice tasta pentru continuare.");
printf("\n====================================");
printf("\n\nRezultatele se gasesc si in fisierul:
d:\\[Link]");
getch();
}
void meniu()
//permite selectarea unei optiuni din meniu si apoi calculul corespunzator
{
clrscr();
printf("\n======= M E T O D A G A U S S =================");
printf("\n\nPuteti alege intre:\n\t1) Rezolvarea unui sistem de
ecuatii liniare");
printf("\n\t2) Determinarea inversei unei matrici\n\t3) Calculul
determinantului unei matrici");
do
{
fflush(stdin);
gotoxy(1,9);
printf("Optiunea: ");
gotoxy(12,9);
scanf("%d",&opt);
}
while(((opt!=1)&&(opt!=2))&&(opt!=3));
}
void verif()
//verifica inversarea matricii prin inmultirea acesteia cu inversa ei
//rezultatul trebuie sa fie matricea unitate
{
int i,j,k;
printf("\nVerificarea inversarii: AxA^-1");
fprintf(f,"\n\nVerificarea inversarii: AxA^-1\n");
for(i=0;i<n;i++)
for(j=0;j<n;j++)
{
prod[i][j]=0;
for(k=0;k<n;k++)
prod[i][j]+=ain[i][k]*inv[k][j];
}
5
Curs 13 Rezolvarea sistemelor de ecuaţii liniare
for(i=0;i<n;i++)
{
printf("\n\t\t");
fprintf(f,"\n\t\t");
for(j=0;j<n;j++)
{
if(fabs(prod[i][j])<1e-8)
prod[i][j]=0;
printf("%g\t",prod[i][j]);
fprintf(f,"%g\t",prod[i][j]);
}
}
}
void solutii()
//calculeaza solutiile sistemului liniar
{
int i,j;
double suma;
x[n-1]=b[n-1]/a[n-1][n-1];
for (j=n-2;j>=0;j--)
{
suma=0;
for(i=j+1;i<n;i++)
suma+=a[j][i]*x[i];
x[j]=(b[j]-suma)/a[j][j];
}
}
void tip_sol()
//tipareste solutiile sistemului
{
int i;
for(i=0;i<n;i++)
{
printf("\nx[%d]=%10.5lf",i,x[i]);
fprintf(f,"\nx[%d]=%10.5lf",i,x[i]);
}
}
void citire_mat()
//citeste matricea de intrare
{
int i,j;
gotoxy(1,10);
printf("Ordinul matricii sistemului: ");
scanf("%d",&n);
gotoxy(1,11);
printf("Matricea sistemului: ");
fprintf(f,"Matricea sistemului: \n");
for(i=0;i<n;i++)
{
fprintf(f,"\n");
for(j=0;j<n;j++)
{
gotoxy(7*j+1,i+12);
scanf("%lf",&a[i][j]);
ain[i][j]=a[i][j];
fprintf(f,"%10.5lf",a[i][j]);
6
Curs 13 Rezolvarea sistemelor de ecuaţii liniare
}
}
}
void citire_tl()
{
//citeste matricea termenilor liberi
int i;
fprintf(f,"\nTermenii liberi:\n");
for(i=0;i<n;i++)
{
gotoxy(7*n+1,i+12);
scanf("%lf",&b[i]);
fprintf(f,"%10.5lf",b[i]);
}
}
double gauss()
{
//realizeaza eliminarea Gaussiana
double t,m;
int i,j,k,l;
for(i=0;i<n;i++)
for(j=0;j<n;j++)
a[i][j]=ain[i][j];
//Se cauta elementul pivot
//de fapt linia care are in coloana k cel mai mare element
for(k=0;k<n;k++)
{
l=k; // l este indexul liniei pivot
for (i=k;i<n;i++)
if (fabs(a[i][k])>fabs(a[k][k]))
l=i;
//Daca elementul pivot este intr-o alta linie decat cea curenta
// atunci se interschimba linia curenta cu cea in care se afla elem. pivot
if(l!=k)
{
for(j=0;j<n;j++)
{
t=a[k][j];a[k][j]=a[l][j];a[l][j]=t;
}
t=b[k];b[k]=b[l];b[l]=t;
det*=-1;
}
//Se calculeaza elementele matricii la eliminarea k
for(i=k+1;i<n;i++)
{
m=a[i][k]/a[k][k];
a[i][k]=0;
for(j=k+1;j<n;j++)
a[i][j]=a[i][j]-m*a[k][j];
b[i]=b[i]-m*b[k];
}
}
//Se calculeaza determinantul matricii
7
Curs 13 Rezolvarea sistemelor de ecuaţii liniare
for(i=0;i<n;i++)
det*=a[i][i];
dtm=det;
return dtm;
}
matinv()
//calculeaza inversa matricii
{
int p,j,i,k;
double suma,t[dm];
for(p=0;p<n;p++)
{
for(j=0;j<n;j++)
{
if(j==p)
t[j]=1;
else
t[j]=0;
b[j]=t[j];
}
gauss();
solutii();
for(k=0;k<n;k++)
inv[k][p]=x[k];
}
printf("\nMatricea inversa este:");
fprintf(f,"\n\nMatricea inversa este:\n");
for(i=0;i<n;i++)
{
printf("\n");
fprintf(f,"\n");
for(j=0;j<n;j++)
{
printf(" %10.5lf",inv[i][j]);
fprintf(f," %10.5lf",inv[i][j]);
}
}
return 0;
}
//========================================================
void main()
{
double x[10],det,suma;
int i,j,k;
f=fopen("d:\\[Link]","w");
meniu();
if(opt==1)
{
citire_mat();
citire_tl();
gauss();
solutii();
tip_sol();
stop();
8
Curs 13 Rezolvarea sistemelor de ecuaţii liniare
}
else
if(opt==2)
{
citire_mat();
matinv();
verif();
stop();
}
else
{
citire_mat();
gauss();
printf("\nDeterminantul matricii este: %lg",dtm);
fprintf(f,"\n\nDeterminantul matricii este: %lg",dtm);
stop();
}
fclose(f);
}