Pagini recente » Cod sursa (job #968999) | Cod sursa (job #1049865) | Cod sursa (job #2271402) | Cod sursa (job #2350404) | Cod sursa (job #2301363)
#include <fstream>
#include <iomanip>
#include <vector>
using namespace std;
using System = vector<vector<double>>;
void SwapRows(System &system, int row1, int row2)
{
if (row1 == row2) {
return;
}
for (size_t i = 0; i < system[row1].size(); i += 1) {
swap(system[row1][i], system[row2][i]);
}
}
int FindPivotRow(const System &system, int col)
{
for (size_t i = col; i < system.size(); i += 1) {
if (system[i][col] != 0) {
return i;
}
}
return -1;
}
void TransformRow(System &system, int dest, int src, double coeff)
{
for (size_t i = 0; i < system[src].size(); i += 1) {
system[dest][i] += system[src][i] * coeff;
}
}
void ClearDownwards(System &system, int col)
{
auto coeff = system[col][col];
for (size_t i = 0; i < system[col].size(); i += 1) {
system[col][i] /= coeff;
}
for (size_t i = col + 1; i < system.size(); i += 1) {
auto coeff = system[i][col] / system[col][col];
TransformRow(system, i, col, -coeff);
}
}
bool MakeSolvable(System &system)
{
auto equations = system.size();
auto vars = system[0].size() - 1;
for (size_t i = 0; i < min(equations, vars); i += 1) {
auto row = FindPivotRow(system, i);
if (row < 0) {
return false;
}
SwapRows(system, i, row);
ClearDownwards(system, i);
}
return true;
}
double FindValueAtIndex(const vector<double> &equation, int index,
const vector<double> &vals)
{
double other_sum = 0;
for (size_t i = index + 1; i + 1 < equation.size(); i += 1) {
other_sum += equation[i] * vals[i];
}
return (equation.back() - other_sum) / equation[index];
}
vector<double> ExtractValues(const System &system)
{
auto equations = system.size();
auto vars = system[0].size() - 1;
vector<double> res(vars, 0);
for (int i = min(equations, vars) - 1; i >= 0; i -= 1) {
res[i] = FindValueAtIndex(system[i], i, res);
}
return res;
}
bool IsSolution(const vector<double> &equation, const vector<double> &res)
{
double sum = 0;
for (size_t i = 0; i < res.size(); i += 1) {
sum += equation[i] * res[i];
}
static constexpr auto kEps = 1e-9;
return abs(sum - equation.back()) < kEps;
}
vector<double> Solve(System &system)
{
if (!MakeSolvable(system)) {
return vector<double>();
}
auto vars = system[0].size() - 1;
auto res = ExtractValues(system);
for (auto i = vars; i < system.size(); i += 1) {
if (!IsSolution(system[i], res)) {
return vector<double>();
}
}
return res;
}
int main()
{
ifstream fin("gauss.in");
ofstream fout("gauss.out");
int equations, vars;
fin >> equations >> vars;
System system(equations, vector<double>(vars + 1));
for (auto &row : system) {
for (auto &elem : row) {
fin >> elem;
}
}
auto res = Solve(system);
if (res.empty()) {
fout << "Imposibil\n";
} else {
for (const auto &elem : res) {
fout << setprecision(9) << fixed << elem << " ";
}
fout << "\n";
}
return 0;
}