#include <stdio.h> #include <cmath> #include <cstdlib> #include <iostream> using namespace std; #define N 100 const float eps = 0.0000001; double b[N],x[N]; float r[N]; float a[N][N]; double A[N][N]; float ans[N]; int n; int select(int x) //部分选主元 { int f = x; for (int i = x + 1; i <= n; i++) if (fabs(a[i][x]) > fabs(a[f][x])) f = i; if (f != x) { for (int i = x; i <= n + 1; i++) swap(a[f][i], a[x][i]); } return 0; } int selectFull(int x) //全选主元 { int px = x, py = x; for (int i = x; i <= n; i++) for (int j = x; j <= n; j++) if (fabs(a[i][j]) > fabs(a[px][py])) { px = i; py = j; } if (px != x || py != x) { for (int i = x; i <= n + 1; i++) swap(a[px][i], a[x][i]); for (int i = 1; i <= n; i++) swap(a[i][py], a[i][x]); //swap(a[py][n + 1], a[x][n + 1]); } return 0; } int Gauss() { float t; for (int i = 1; i <= n; i++) { selectFull(i); if (fabs(a[i][i]) <= eps) { printf("No Answers\n"); return 0; } for (int j = i + 1; j <= n; j++) { if (fabs(a[j][i]) > eps) { t = a[i][i] / a[j][i]; for (int k = i; k <= n + 1; k++) a[j][k] = a[j][k] * t - a[i][k]; } } } //解出答案 for (int i = n; i; i--) { for (int j = i + 1; j <= n; j++) { a[i][n + 1] -= a[i][j] * a[j][n + 1]; } a[i][n + 1] /= a[i][i]; } return 0; } int work() //计算残差 { for (int i = 1; i <= n; i++) { b[i] = A[i][n + 1]; x[i] = a[i][n + 1]; } for (int i = 1; i <= n; i++) { double t = 0; for (int j = 1; j <= n; j++) t += A[i][j] * x[j]; r[i] = b[i] - t; } return 0; } int main() { scanf_s("%d", &n); for (int i = 1; i <= n; i++) for (int j = 1; j <= n + 1; j++) { cin >> a[i][j]; A[i][j] = a[i][j]; } Gauss(); swap(a[1][n + 1], a[4][n + 1]); swap(a[3][n + 1], a[2][n + 1]); for (int i = 1; i <= n; i++) printf("%f\n", a[i][n + 1]); work(); printf("残差为:\n"); for (int i = 1; i <= n; i++) printf("%f\n", r[i]); return 0; }
高斯消元和残差计算
发布于 2021-09-21 10 次阅读
Comments NOTHING