高斯消元 hdu3359 Kind of a Blur

最基础的高斯消元(浮点数)版，可以用来测试模板

这个模板直接从大白书上敲下来的，感觉还行

讲解一下这个模板的用法

构造的是A*r=b

行列的下标都是从0开始

A[i][j]是增广矩阵的第i行的第j列

b存放在对应的A[i][n]里面

调用函数后，A[i][n]里面存放的就是对应的r矩阵了

#include<cstdio>
#include<map>
#include<set>
#include<cmath>
#include<cstring>
#include<queue>
#include<stack>
#include<cctype>
#include<vector>
#include<functional>
#include<algorithm>

using namespace std;
typedef long long LL;
typedef pair<int, int> PII;

const int MX = 100 + 5;
const int INF = 0x3f3f3f3f;
typedef double Matrix[MX][MX];

Matrix A, S;

void gauss(Matrix A, int n) {
    int i, j, k, r;
    for(i = 0; i < n; i++) {
        r = i;
        for(int j = i + 1; j < n; j++) {
            if(fabs(A[j][i]) > fabs(A[r][i])) r = j;
        }
        if(r != i) for(j = 0; j <= n; j++) swap(A[r][j], A[i][j]);

        for(k = i + 1; k < n; k++) {
            double f = A[k][i] / A[i][i];
            for(j = i; j <= n; j++) A[k][j] -= f * A[i][j];
        }
    }

    for(i = n - 1; i >= 0; i--) {
        for(j = i + 1; j < n; j++) {
            A[i][n] -= A[j][n] * A[i][j];
        }
        A[i][n] /= A[i][i];
    }
}

int main() {
    //freopen("input.txt", "r", stdin);
    int m, n, d, first = true;
    while(~scanf("%d%d%d", &n, &m, &d), m) {
        memset(A, 0, sizeof(A));

        for(int i = 0; i < m; i++) {
            for(int j = 0; j < n; j++) {
                scanf("%lf", &S[i][j]);
            }
        }

        for(int i = 0; i < m; i++) {
            for(int j = 0; j < n; j++) {
                int cnt = 0, nx, ny;
                for(int dx = -d; dx <= d; dx++) {
                    nx = dx + i;
                    if(nx < 0 || nx >= m) continue;

                    int p = d - abs(dx);
                    for(int dy = -p; dy <= p; dy++) {
                        ny = dy + j;
                        if(ny < 0 || ny >= n) continue;

                        cnt++;
                        A[i * n + j][nx * n + ny] = 1;
                    }
                }
                A[i * n + j][m * n] = cnt * S[i][j];
            }
        }

        gauss(A, m * n);

        if(first) first = false;
        else      printf("\n");

        for(int i = 0; i < m; i++) {
            for(int j = 0; j < n; j++) {
                printf("%8.2lf", A[i * n + j][m * n]);
            }
            printf("\n");
        }
    }
    return 0;
}