Pagini recente » Cod sursa (job #2568461) | Cod sursa (job #1979180) | Cod sursa (job #2153436) | Cod sursa (job #1624167) | Cod sursa (job #861089)
Cod sursa(job #861089)
#include <cstdio>
#define MAXN 302
#define EPS 0.000000001
using namespace std;
int N, M;
double A[MAXN][MAXN];
double X[MAXN];
int main()
{ freopen("gauss.in", "r", stdin); freopen("gauss.out", "w", stdout);
//Citire
scanf("%d%d", &N, &M);
for(int i = 1; i <= N; ++i)
for(int j = 1; j <= M+1; ++j)
scanf("%lf", &A[i][j]);
int i = 1, j = 1, k;
double aux;
//Algoritmul lui Gauss propriu-zis
while(i <= N && j <= M)
{ //Cautam o linie k pentru care A[k][j] sa fie nenul. Folosim epsilonul pentru a nu lua in considerare micile erori de calcul.
for(k = i; k <= N; ++k)
if(A[k][j]<-EPS || A[k][j]>EPS) break;
//Daca nu gasim linia, necunoscuta j este variabila libera, incrementam pe j si trecem la pasul urmator.
if(k == N+1) {++j; continue;}
//Interschimbam linia k cu linia i, daca k != i.
if(k != i)
for(int l = 1; l <= M+1; ++l)
{aux = A[i][l]; A[i][l] = A[k][l]; A[k][l] = aux;}
//Pentru a usura calculele, impartim toate valorile de pe linia i la A[i][j], A[i][j] devenind 1.
//Observam ca valorile de pe linia i si coloanele 1..j-1 sunt egale cu 0 de la pasii precedenti ai algoritmului,
//deci nu e necesar sa le parcurgem pentru a le imparti.
for(int l = j+1; l <= M+1; ++l) A[i][l] = A[i][l] / A[i][j];
A[i][j] = 1;
//Scadem din ecuatiile i+1...N ecuatia i inmultita cu A[u][j], pentru a egala toti coeficientii de coloana j
//a acestor linii la 0.
for(int u = i+1; u <= N; ++u)
{ for(int l = j+1; l <= M+1; ++l) A[u][l] -= A[u][j] * A[i][l];
A[u][j] = 0;
}
++i; ++j;
}
//Calculul necunoscutelor
for(int i = N; i>0; --i)
for(int j = 1; j <= M+1; ++j)
if(A[i][j]>EPS || A[i][j]<-EPS)
{ //Singura valoare nenegativa de pe linia i este rezultatul => sistemul nu are solutie.
if(j == M+1) {printf("Imposibil\n"); return 0;}
//Calculam pe necunoscuta j = rezultatul ecuatiei i - necunoscutele j+1...M inmultite cu coeficientii lor din aceasta ecuatie.
X[j] = A[i][M+1];
for(int k = j+1; k <= M; ++k) X[j] -= X[k] * A[i][k];
break;
}
//Afisare
for(int i = 1; i <= M; ++i) printf("%.8lf ", X[i]);
printf("\n"); return 0;
}