Pagini recente » Cod sursa (job #758575) | Cod sursa (job #581276) | Monitorul de evaluare | Diferente pentru happy-birthday-infoarena-2014 intre reviziile 3 si 4 | Cod sursa (job #1283170)
#include <fstream>
#include <iostream>
#include <iomanip>
#include <algorithm>
#include <cmath>
#define EPS 0.0000001
using namespace std;
ofstream g("gauss.out");
void debug(double**& A, const int& N, const int& M)
{
for(int i = 1; i <= N; i++, cout << endl)
for(int j = 1; j <= M + 1; j++)
cout << A[i][j] << " ";
cout << endl;
}
void read_input(double**& A, int& N, int& M)
{
ifstream f("gauss.in");
f >> N >> M;
A = new double*[N + 1];
for(int i = 1; i <= N; i++)
{
A[i] = new double[M + 2];
}
for(int i = 1; i <= N; i++)
{
for(int j = 1; j <= M + 1; j++)
{
f >> A[i][j];
}
}
f.close();
}
int argmax(double**& A, const int& N, const int& col)
{
int row_max = col;
for(int j = col + 1; j <= N; j++)
{
if(abs(A[j][col]) > abs(A[row_max][col]))
{
row_max = col;
}
}
return row_max;
}
void swap_rows(double**& A, const int& N, const int& M, const int& k, const int& pivot_index)
{
for(int i = 1; i <= M + 1; i++)
{
swap(A[k][i], A[pivot_index][i]);
}
}
// Compute in-place row echelon form
bool compute_row_echelon_form(double**& A, const int& N, const int& M)
{
//debug(A, N, M);
for(int k = 1; k < min(N, M); k++)
{
int pivot_index = argmax(A, N, k);
//cout << "Pivot index " << pivot_index << endl;
if(A[pivot_index][k] > -EPS && A[pivot_index][k] < EPS)
{
return false;
}
if(k != pivot_index)
{
swap_rows(A, N, M, k, pivot_index);
}
for(int i = k + 1; i <= N; i++)
{
for(int j = k + 1; j <= M + 1; j++)
{
A[i][j] = A[i][j] - A[k][j] * (A[i][k] / A[k][k]);
}
A[i][k] = 0;
}
//debug(A, N, M);
}
return true;
}
bool compute_unknowns(double*& X, double**& A, const int& N, const int& M)
{
//cout << "compute_unknowns " << endl;
//Calculul necunoscutelor.
for(int i = N; i >= 1; i--)
{
for(int j = 1; j <= M + 1; j++)
{
if(A[i][j] < -EPS || A[i][j] > EPS)
{
// Toti coeficientii sunt 0 iar rezultatul este diferit de 0.
if(j == M + 1) // adica ultima coloana, b-ul
{
return false;
}
// 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];
}
X[j] /= A[i][i];
break;
}
}
}
return true;
}
int main()
{
int N, M;
double** A;
double* X;
read_input(A, N, M);
X = new double[M + 2];
bool ret = compute_row_echelon_form(A, N, M);
if(!ret)
{
cout << "Imposibil\n";
}
else
{
ret = compute_unknowns(X, A, N, M);
if(!ret)
{
cout << "Imposibil\n";
}
else
{
for(int i = 1; i <= M; i++)
{
cout << std::fixed << std::setprecision(10) << X[i] << " ";
}
}
}
return 0;
}