Pagini recente » Cod sursa (job #2490469) | Cod sursa (job #731364) | Cod sursa (job #2174246) | Cod sursa (job #1888261) | Cod sursa (job #2497069)
#include <iostream>
#include <fstream>
#include <iomanip>
#include <sys/mman.h>
#include <sys/stat.h>
#include <sys/types.h>
#include <fcntl.h>
#include <unistd.h>
std::ofstream t("gauss.out");
#define NMAX 310
#define EPSILON 1e-8
template <typename T>
using Matrix = T[NMAX][NMAX];
template <typename T>
using Vector = T[NMAX];
class mapping_parser {
public:
inline mapping_parser() {
/// default c-tor
}
inline mapping_parser(const char* file_name) {
int fd = open(file_name, O_RDONLY);
index &= 0;
fstat(fd, &sb);
buffer = (char*)mmap(0, sb.st_size, PROT_READ, MAP_PRIVATE | MAP_NORESERVE, fd, 0);
close(fd);
}
inline mapping_parser& operator >> (double& n) {
aux = 0;
sign &= 0;
for (; buffer[index] < '0' or buffer[index] > '9'; ++index);
sign |= (buffer[index - 1] == '-');
for (; '0' <= buffer[index] and buffer[index] <= '9'; ++index)
aux = (aux << 3) + (aux << 1) + buffer[index] - '0';
aux ^= ((aux ^ -aux ) & -sign);
n = aux;
return *this;
}
inline mapping_parser& operator >> (int& n) {
n = 0;
for (; buffer[index] < '0' or buffer[index] > '9'; ++index);
for (; '0' <= buffer[index] and buffer[index] <= '9'; ++index)
n = (n << 3) + (n << 1) + buffer[index] - '0';
return *this;
}
~mapping_parser() {
munmap(buffer, sb.st_size);
}
private:
int aux;
struct stat sb;
int index, sign;
char* buffer;
};
mapping_parser f("gauss.in");
int main() {
Matrix<double> A;
Vector<double> b;
Vector<double> x = {};
int n, m;
f >> n >> m;
auto abs = [](double target) {
return target < 0 ? -target : target;
};
auto subtract_line_with_scalar = [&A, &b, m](int dst, int src, double scalar, int start = 0) mutable {
for (int i = start; i < m; ++i)
A[dst][i] -= A[src][i] * scalar;
b[dst] -= b[src] * scalar;
};
auto divide_line_with_scalar = [&A, &b, m](int dst, double scalar, int start = 0) mutable {
for (int i = start; i < m; ++i)
A[dst][i] /= scalar;
b[dst] /= scalar;
};
// reading the matrix and results
for (int i = 0; i < n; ++i) {
for (int j = 0; j < m; ++j)
f >> A[i][j];
f >> b[i];
}
for (int i = 0; i < n; ++i) {
for (int j = 0; j < m; ++j)
std::cout << A[i][j];
std::cout << std::endl;
}
for (int i = 0; i < m; ++i) { // for each column
divide_line_with_scalar(i, A[i][i], i); // no pivoting, but scaling
for (int j = i + 1; j < n; ++j) { // rest to 0
subtract_line_with_scalar(j, i, A[j][i], i);
}
}
double sigma;
// now perform back substitution
for (int i = n - 1, aux = n < m ? n - 1 : m - 1; i >= 0; --i) {
sigma = .0;
for (int j = n < m ? n - 1 : m - 1; j > aux; --j) {
sigma += x[j] * A[i][j];
}
if (abs(sigma + A[i][aux]) < abs(EPSILON) && abs(b[i]) > abs(EPSILON)) {
t << "Imposibil\n";
return 0;
}
else if (A[i][aux] > abs(EPSILON)) {
x[aux] = (b[i] - sigma) / A[i][aux];
--aux;
}
}
t << std::setprecision(10) << std::fixed;
for (int i = 0; i < m; ++i)
t << x[i] << " ";
return 0;
}