写的好累啊，还没有AC，应该是卡精度了，悲剧，再想想吧。
Wrong Answer
#include <iostream> #include <stdio.h> #include <string.h> #include <math.h> using namespace std; #define zero(x) ((x>0? x:-x)<1e-15) int const MAXN = 1010; double a[MAXN][MAXN]; double b[MAXN][MAXN]; int g[505][505]; int n,m; double R; struct Edge{ double r,k; }edge[505][505]; struct Node{ double x,y; }node[505]; double dis(double x1,double y1,double x2,double y2){ return sqrt((x1-x2)*(x1-x2)+(y1-y2)*(y1-y2)); } double krc(double x1,double y1,double x2,double y2){ if((x1-x2)==0) return 540000.0; else return (y1-y2)/(x1-x2); } bool is_bet(int i,int z,int j){ if(node[i].x<=node[j].x && node[i].y <= node[j].y && node[i].x<=node[z].x && node[z].x<=node[j].x && node[i].y<=node[z].y && node[z].y<=node[j].y) return 1; else if(node[i].x>=node[j].x && node[i].y <= node[j].y && node[i].x>=node[z].x && node[z].x>=node[j].x && node[i].y<=node[z].y && node[z].y<=node[j].y) return 1; else if(node[i].x<=node[j].x && node[i].y >= node[j].y && node[i].x<=node[z].x && node[z].x<=node[j].x && node[i].y>=node[z].y && node[z].y>=node[j].y) return 1; else if(node[i].x>=node[j].x && node[i].y >= node[j].y && node[i].x>=node[z].x && node[z].x>=node[j].x && node[i].y>=node[z].y && node[z].y>=node[j].y) return 1; else return 0; } bool judge_bet(int i,int j){ for(int z=i+1;z<n;z++){ if(edge[i][z].k==edge[i][j].k && edge[i][z].r<edge[i][j].r && is_bet(i,z,j)){ return 1; } } return 0; } double det(double a[MAXN][MAXN],int n) { int i, j, k, sign = 0; double ret = 1, t; for (i = 0; i < n; i++) for (j = 0; j < n; j++) b[i][j] = a[i][j]; for (i = 0; i < n; i++) { if (zero(b[i][i])) { for (j = i + 1; j < n; j++) if (!zero(b[j][i])) break; if (j == n) return 0; for (k = i; k < n; k++) t = b[i][k], b[i][k] = b[j][k], b[j][k] = t; sign++; } ret *= b[i][i]; for (k = i + 1; k < n; k++) b[i][k] /= b[i][i]; for (j = i + 1; j < n; j++) for (k = i + 1; k < n; k++) b[j][k] -= b[j][i] * b[i][k]; } if (sign & 1) ret = -ret; return ret; } void build(){ memset(g,0,sizeof(g)); for(int i=0;i<n;i++){ for(int j=0;j<n;j++){ edge[i][j].r=dis(node[i].x,node[i].y,node[j].x,node[j].y); edge[i][j].k=krc(node[i].x,node[i].y,node[j].x,node[j].y); } edge[i][i].r=0; edge[i][i].k=0; } for(int i=0;i<n;i++){ for(int j=i+1;j<n;j++){ if(edge[i][j].r<=R && !judge_bet(i,j)){ g[i][j]=g[j][i]=1; } } } } void solve(){ memset(a,0,sizeof(a)); for (int i=0;i<n;i++) { int d=0; for (int j=0;j<n;j++) if (g[i][j]) d++; a[i][i]=d; } for (int i=0;i<n;i++) for (int j=0;j<n;j++) if (g[i][j]) a[i][j]=-1; } bool judge(){ for(int i=0;i<n;i++){ int flag=0; for(int j=0;j<n;j++){ if(g[i][j]!=0) flag=1; } if(!flag) return 1; } return 0; } int main() { //freopen("in.txt","r",stdin); int ncas; double x,y; cin >> ncas; while(ncas--){ scanf("%d%lf",&n,&R); for(int i=0;i<n;i++){ scanf("%lf%lf",&x,&y); node[i].x=x; node[i].y=y; } build(); if(judge()){ printf("-1\n"); continue; } solve(); double ans=det(a,n-1); printf("%0.0lf\n",ans+0.001); } return 0; }
终于AC了,原来是忘了MOD的问题了.
注意一下解行列式的过程，由于结果非常大，所以需要mod一个10007.在做除法的时候要改为乘以其逆元。
一定要注意求逆元的那些数必须为正数（即代表度数矩阵的位置）。
b[j][k] =((b[j][k]- b[j][i] * b[i][k]))%MOD;
if(j==k) b[j][k]=(b[j][k]+MOD)%MOD;
如果WA在这里,是因为本以为只要保证最后的结果为正数就可以了，中间高斯消元严格按照过程来就可以了。事实上是不对的。
问题出在扩展欧几里得算法上。
举一个例子。当我们想知道 -1*x=1(mod 7)的时候，如果用exgcd(-1,7,x,y)求出来的x其实是1，而不是我们所期望的6！！ 为什么？扩展欧几里得此时求出来的值并不是最大公约数1，而是-1！计算的是-1*x=-1(mod 7)！！所以我们应当将1转为6才能得到正确答案。
还有一个需要注意的地方是任意交换行列式的两行，行列式计算出来的结果需要乘以-1，这个非常容易证明。而我们需要求的是基尔霍夫行列式的绝对值，因此需要记录一下消元过程中行列式交换的次数，如果是奇数，所得到的ans应该为mod-ans.
RANK-1,实在是个惊喜,也许这就是我最好的生日礼物.
附AC代码:
Accepted
#include <iostream> #include <stdio.h> #include <string.h> #include <math.h> using namespace std; #define MOD 10007 #define zero(x) ((x>0? x:-x)<1e-15) int const MAXN = 1010; int a[MAXN][MAXN]; int b[MAXN][MAXN]; int g[505][505]; int n,m; double R; struct Edge{ double r,k; }edge[505][505]; struct Node{ double x,y; }node[505]; double dis(double x1,double y1,double x2,double y2){ return sqrt((x1-x2)*(x1-x2)+(y1-y2)*(y1-y2)); } double krc(double x1,double y1,double x2,double y2){ if((x1-x2)==0) return 540000.0; else return (y1-y2)/(x1-x2); } int exgcd(int a,int b,int &x,int &y) { if (b==0) { x=1;y=0;return a; } int d=exgcd(b,a%b,x,y); int t=x;x=y;y=t-a/b*y; return d; } bool is_bet(int i,int z,int j){ if(node[i].x<=node[j].x && node[i].y <= node[j].y && node[i].x<=node[z].x && node[z].x<=node[j].x && node[i].y<=node[z].y && node[z].y<=node[j].y) return 1; else if(node[i].x>=node[j].x && node[i].y <= node[j].y && node[i].x>=node[z].x && node[z].x>=node[j].x && node[i].y<=node[z].y && node[z].y<=node[j].y) return 1; else if(node[i].x<=node[j].x && node[i].y >= node[j].y && node[i].x<=node[z].x && node[z].x<=node[j].x && node[i].y>=node[z].y && node[z].y>=node[j].y) return 1; else if(node[i].x>=node[j].x && node[i].y >= node[j].y && node[i].x>=node[z].x && node[z].x>=node[j].x && node[i].y>=node[z].y && node[z].y>=node[j].y) return 1; else return 0; } bool judge_bet(int i,int j){ for(int z=i+1;z<n;z++){ if(edge[i][z].k==edge[i][j].k && edge[i][z].r<edge[i][j].r && is_bet(i,z,j)){ return 1; } } return 0; } int det(int a[MAXN][MAXN],int n) { int i, j, k, sign = 0; int ret = 1, t; for (i = 0; i < n; i++) for (j = 0; j < n; j++) b[i][j] = a[i][j]; for (i = 0; i < n; i++) { if (zero(b[i][i])) { for (j = i + 1; j < n; j++) if (!zero(b[j][i])) break; if (j == n) return 0; for (k = i; k < n; k++) t = b[i][k], b[i][k] = b[j][k], b[j][k] = t; sign++; } ret = ret*b[i][i]%MOD; int x,y; int tp=exgcd(b[i][i],MOD,x,y); for (k = i + 1; k < n; k++) b[i][k] =b[i][k]*x%MOD; for (j = i + 1; j < n; j++) for (k = i + 1; k < n; k++){ b[j][k] =((b[j][k]- b[j][i] * b[i][k]))%MOD; if(j==k) b[j][k]=(b[j][k]+MOD)%MOD; } } ret=(ret%MOD+MOD)%MOD; if (sign & 1) return MOD-ret; return ret; } void build(){ memset(g,0,sizeof(g)); for(int i=0;i<n;i++){ for(int j=0;j<n;j++){ edge[i][j].r=dis(node[i].x,node[i].y,node[j].x,node[j].y); edge[i][j].k=krc(node[i].x,node[i].y,node[j].x,node[j].y); } edge[i][i].r=0; edge[i][i].k=0; } for(int i=0;i<n;i++){ for(int j=i+1;j<n;j++){ if(edge[i][j].r<=R && !judge_bet(i,j)){ g[i][j]=g[j][i]=1; } } } } void solve(){ memset(a,0,sizeof(a)); for (int i=0;i<n;i++) { int d=0; for (int j=0;j<n;j++) if (g[i][j]) d++; a[i][i]=d; } for (int i=0;i<n;i++) for (int j=0;j<n;j++) if (g[i][j]) a[i][j]=-1; } bool judge(){ for(int i=0;i<n;i++){ int flag=0; for(int j=0;j<n;j++){ if(g[i][j]!=0) flag=1; } if(!flag) return 1; } return 0; } int main() { //freopen("in.txt","r",stdin); int ncas; double x,y; cin >> ncas; while(ncas--){ scanf("%d%lf",&n,&R); for(int i=0;i<n;i++){ scanf("%lf%lf",&x,&y); node[i].x=x; node[i].y=y; } build(); if(judge()){ printf("-1\n"); continue; } solve(); int ans=det(a,n-1); printf("%d\n",ans); } return 0; }
