HDU 4305 Lightning (高斯消元解kirchhoff矩阵+逆元)

　　　题意是：给一些坐标点，如果两点之间的距离小于R，并且两点之间没有其他点，则这两个点保持连通，这样构成了一个图。问这个图中生成树的个数。

因为数据量并不大，O(N^3)的建图没有问题。

建好图以后就可以用kirchhoff矩阵计算生成树的个数，之所以写这道题的解题报告是因为在高斯消元解kirchhoff矩阵时，需要用到逆元。

(a/b)% mod，如果a,b的范围很大，结果会有很大的误差，这里可以转换一下 b*x = 1(% mod)则x为b的逆元 (a/b)%mod = (a*x)%mod

求逆元的过程就是解线性同余方程b*x ≡ 1(% mod)。

详见代码：

/*
高斯消元的整个过程中不能出现负数
*/

//#pragma comment(linker,"/STACK:327680000,327680000")
#include <iostream>
#include <cstdio>
#include <cmath>
#include <vector>
#include <cstring>
#include <algorithm>
#include <string>
#include <set>
#include <functional>
#include <numeric>
#include <sstream>
#include <stack>
#include <map>
#include <queue>

#define CL(arr, val) memset(arr, val, sizeof(arr))
#define REP(i, n)for((i) = 0; (i) < (n); ++(i))
#define FOR(i, l, h)for((i) = (l); (i) <= (h); ++(i))
#define FORD(i, h, l)for((i) = (h); (i) >= (l); --(i))
#define L(x)(x) << 1
#define R(x)(x) << 1 | 1
#define MID(l, r)(l + r) >> 1
#define Min(x, y)x < y ? x : y
#define Max(x, y)x < y ? y : x
#define E(x)(1 << (x))
#define iabs(x) (x) < 0 ? -(x) : (x)
#define OUT(x)printf("%I64d\n", x)
#define lowbit(x)(x)&(-x)
#define Read()freopen("data.in", "r", stdin)
#define Write()freopen("data.out", "w", stdout);

const double eps = 1e-8;
typedef long long LL;
const int inf = ~0u>>2;

using namespace std;

const int N = 330;
const int MOD = 10007;

struct node {
    double x;
    double y;
    double dis(node cmp) {
        return (cmp.x - x)*(cmp.x - x)  + (cmp.y - y)*(cmp.y - y);
    }
}p[N];

int exp_gcd(int a, int b, int &x, int &y) {
    if(b == 0)  {
        x = 1; y = 0;
        return a;
    }
    int p = exp_gcd(b, a%b, x, y);
    int tmp = x;
    x = y; y = tmp - (a/b)*y;
    return p;
}

int det(int a[][N], int n) {
    int i, j, k;
    int ans = 1, x, y;
    bool flag = true;

    for(i = 0; i < n; ++i) {
        if(!a[i][i]) {
            for(j = i + 1; j < n; ++j) {
                if(a[j][i]) {
                    for(k = i; k < n; ++k)  swap(a[i][k], a[j][k]);
                    flag = !flag;
                    break;
                }
            }
            if(j == n)  return -1;
        }
        ans = (ans*a[i][i]%MOD + MOD)%MOD;
        exp_gcd(a[i][i], MOD, x, y);
        x = (x%MOD + MOD)%MOD;

        for(k = i + 1; k < n; ++k)  a[i][k] = (a[i][k]*x%MOD + MOD)%MOD;

        for(j = i + 1; j < n; ++j) {
            if(a[j][i]) {
                for(k = i + 1; k < n; ++k) {
                    a[j][k] = ((a[j][k] - a[j][i]*a[i][k])%MOD + MOD)%MOD;
                }
                a[j][i] = 0;
            }
        }
    }

    if(flag)    return ans;
    return (-ans%MOD + MOD)%MOD;
}

int A[N][N], D[N][N], C[N][N];

int main() {
    //Read();

    int T, n, r, i, j, k;
    cin >> T;
    while(T--) {
        cin >> n >> r;
        REP(i, n)   cin >> p[i].x >> p[i].y;
        CL(A, 0);
        REP(i, n) {
            for(j = i + 1; j < n; ++j) {
                //printf("%d %d %f\n", i, j, p[i].dis(p[j]));
                if(p[i].dis(p[j]) > r*r)    continue;

                for(k = 0; k < n; ++k) {
                    if(k == i || k == j )   continue;
                    if((p[k].y - p[i].y)*(p[j].x - p[k].x) == (p[k].x - p[i].x)*(p[j].y - p[k].y) &&
                       p[i].dis(p[j]) > p[i].dis(p[k]) && p[i].dis(p[j]) > p[k].dis(p[j]))   break;
                }
                if(k == n)  A[i][j] = A[j][i] = 1;
            }
        }

        CL(D, 0);
        bool flag = true;
        REP(i, n) {
            D[i][i] = 0;
            REP(j, n) {
                if(A[i][j]) D[i][i]++;
            }
            if(D[i][i] == 0)    {flag = false; break;}
        }

        if(!flag)   {puts("-1"); continue;}
        REP(i, n) {
            REP(j, n) {
                C[i][j] = D[i][j] - A[i][j];
                C[i][j] = (C[i][j]%MOD + MOD)%MOD;
            }
        }

        printf("%d\n", det(C, n-1));
    }
    return 0;
}