HDU - 3359 Kind of a Blur - 高斯-约当消元

题目描述

题目有坑，h,w=0,但d可能不等于0就要结束。
h,w输入的顺序也是不合常理，唉，被坑了。

#include<cstdio>
#include<cstring>
#include<cmath>
#include<algorithm>
using namespace std;
#define MAXH 10
#define MAXN 100
const double eps=1e-12;

int h,w,d,id[MAXH+10][MAXH+10];
double a[MAXN+10][MAXN+10],x[MAXN+10];
bool free_x[MAXN+10];

void Init()
{
    memset(a,0,sizeof a);
    memset(x,0,sizeof x);
    memset(id,0,sizeof id);
}
int Dist(int p,int q,int x,int y)
{
    return abs(y-q)+abs(x-p);
}
void read()
{
    double xx;
    for(int i=1;i<=h;i++)
        for(int j=1;j<=w;j++)
            id[i][j]=(i-1)*w+j;
    for(int i=1;i<=h;i++)
        for(int j=1;j<=w;j++){
            scanf("%lf",&xx);
            int cnt=0,t=id[i][j];
            for(int p=1;p<=h;p++)
                for(int q=1;q<=w;q++)
                    if(Dist(p,q,i,j)<=d){
                        a[t][id[p][q]]=1.0;
                        cnt++;
                    }
            a[t][h*w+1]=1.0*cnt*xx;
        }
}
void GJ_elimination(int equ,int var,int &row,int &col)
{
    int mx;
    for(row=col=1;row<=equ&&col<=var;row++,col++){
        mx=row;
        for(int i=row+1;i<=equ;i++){
            if(fabs(a[i][col])>fabs(a[mx][col]))
                mx=i;
        }
        if(mx!=row){
            for(int j=1;j<=var+1;j++)
                swap(a[mx][j],a[row][j]);
        }
        if(fabs(a[row][col])<eps){
            row--;
            continue;
        }
        for(int i=1;i<=equ;i++){
            if(i==row||fabs(a[i][col])<eps)
                continue;
            for(int j=var+1;j>=col;j--)
                a[i][j]-=a[i][col]/a[row][col]*a[row][j];
        }
    }
}
int GJ_Judge(int equ,int var,int row,int col)
{
    //0:唯一解
    for(int i=row-1;i>=1;i--){
        double tmp=0.0;
        for(int j=i+1;j<=var;j++)
            if(fabs(a[i][j])>eps)
                tmp+=a[i][j]*x[j];
        x[i]=(a[i][var+1]-tmp)/a[i][i];
    }
    return 0;
}
void Gauss_Jordan(int equ,int var)
{
    int row,col;
    GJ_elimination(equ,var,row,col);
    GJ_Judge(equ,var,row,col);
}
void print()
{
    for(int i=1;i<=h*w;i++){
        printf("%8.2lf",x[i]);
        if(i%w==0)
            printf("\n");
    }
}
int main()
{
    bool flag=false;
    while(scanf("%d%d%d",&w,&h,&d)&&w&&h){
        if(flag) puts("");
        Init();
        read();
        Gauss_Jordan(h*w,h*w);
        print();
        flag=true;
    }
}