求生成树的个数（矩阵+乘法逆元）
http://acm.hdu.edu.cn/showproblem.php?pid=4305
Lightning
Suddenly the sky turns into gray, and lightning storm comes! Unfortunately, one of the robots is stuck by the lightning!
So it becomes overladen. Once a robot becomes overladen, it will spread lightning to the near one.
The spreading happens when:
Robot A is overladen but robot B not.
The Distance between robot A and robot B is no longer than R.
No other robots stand in a line between them.
In this condition, robot B becomes overladen.
We assume that no two spreading happens at a same time and no two robots stand at a same position.
The problem is: How many kind of lightning shape if all robots is overladen? The answer can be very large so we output the answer modulo 10007. If some of the robots cannot be overladen, just output -1.
The first line is an integer T (T < = 20), indicate the test cases.
For each case, the first line contains integer N ( 1 < = N < = 300 ) and R ( 0 < = R < = 20000 ), indicate there stand N robots; following N lines, each contains two integers ( x, y ) ( -10000 < = x, y < = 10000 ), indicate the position of the robot.
3 3 2 -1 0 0 1 1 0 3 2 -1 0 0 0 1 0 3 1 -1 0 0 1 1 0
3 1 -1
题意：求生成树的个数;
平面上有N个点。每个两个点如果距离小于R且之间没有共线的另一个点，则这两点之间有一条边。求这个图的生成树的个数mod 10007。
Matrix-Tree定理(Kirchhoff矩阵-树定理)。Matrix-Tree定理是解决生成树计数问题最有力的武器之一。它首先于1847年被Kirchhoff证明。在介绍定理之前，我们首先明确几个概念：
1、G的度数矩阵D[G]是一个n*n的矩阵，并且满足：当i≠j时,dij=0；当i=j时，dij等于vi的度数。
2、G的邻接矩阵A[G]也是一个n*n的矩阵， 并且满足：如果vi、vj之间有边直接相连，则aij=1，否则为0。
我们定义G的Kirchhoff矩阵(也称为拉普拉斯算子)C[G]为C[G]=D[G]-A[G]，则Matrix-Tree定理可以描述为：G的所有不同的生成树的个数等于其Kirchhoff矩阵C[G]任何一个n-1阶主子式的行列式的绝对值。所谓n-1阶主子式，就是对于r(1≤r≤n)，将C[G]的第r行、第r列同时去掉后得到的新矩阵，用Cr[G]表示。
矩阵的规则是：
1、在主对角线上的元素为此节点的度数
2、对于其他位置上的元素Matrix(i,j) { i != j }，
(1) 如果节点i和节点j连通，则Matrix(i,j)的值为-k，其中k值为节点i到节点j的平行边个数。如果此图是一个简单图，即任意两点间不存在平行边，那么这个值就为-1.
(2) 但如果节点i和节点j根本不连通，则Matrix(i,j)的值为0。
【注意一下解行列式的过程，由于结果非常大，所以需要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.】
行列式求值 的性质：
（1）对换行列式的任意两行（列），行列式的值变号；
（2）行列式的某一行（列）中所有元素都乘以同一个数k，等于用数k乘以此行列式的值；
（3）把行列式的某一行（列）的各元素都乘以同一数后加到另一行（列）对应元素上，行列式的值不变。
程序：
#include"stdio.h"
#include"string.h"
#include"queue"
#include"stack"
#include"math.h"
#include"iostream"
#define M 555
#define inf 100000000
#define mod 10007
#define eps 1e-10
using namespace std;
struct Node
{
double x,y;
}p[M];
int use[M],G[M][M];
double pow(double x)
{
return x*x;
}
double Len(Node a,Node b)
{
return sqrt(pow(a.x-b.x)+pow(a.y-b.y));
}
void put(int n)
{
printf("------------------------------\n");
for(int i=1;i<=n;i++)
{
for(int j=1;j<=n;j++)
printf("%d ",G[i][j]);
puts("");
}
printf("------------------------------\n");
return ;
}
struct node
{
int v;
node(int vv)
{
v=vv;
}
};
vector<node>edge[M];
void dfs(int u)
{
use[u]=1;
for(int i=0;i<(int)edge[u].size();i++)
{
int v=edge[u][i].v;
if(!use[v])
dfs(v);
}
}
int ok(int n)
{
memset(use,0,sizeof(use));
int flag=0;
for(int i=1;i<=n;i++)
{
if(!use[i])
{
flag++;
dfs(i);
}
}
if(flag==1)
return 1;
else
return 0;
}
int exgcd(int a,int b,int &x,int &y)//乘法逆元返回的d是a,b的公约数，x是a mod b的逆元
{
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;
}
int det(int n)
{
int ans=1;
int flag=1;
int i,j,k;
for(i=1;i<=n;i++)
{
if(G[i][i]==0)
{
for(j=i+1;j<=n;j++)
{
if(G[j][i])break;
}
if(j>n)return 0;
flag=!flag;
for(k=i;k<=n;k++)
swap(G[i][k],G[j][k]);
}
ans=ans*G[i][i]%mod;
int x,y;
int tep=exgcd(G[i][i],mod,x,y);
for(k=i+1;k<=n;k++)
{
G[i][k]=G[i][k]*x%mod;
}
for(j=i+1;j<=n;j++)
{
for(k=i+1;k<=n;k++)
{
G[j][k]=(G[j][k]-G[j][i]*G[i][k])%mod;
if(j==k)
G[j][k]=(G[j][k]+mod)%mod;
}
}
}
ans=(ans%mod+mod)%mod;
if(flag)return ans;
else
return mod-ans;
}
int main()
{
int T;
scanf("%d",&T);
while(T--)
{
int n,i,j;
double R;
scanf("%d%lf",&n,&R);
for(i=1;i<=n;i++)
edge[i].clear();
for(i=1;i<=n;i++)
scanf("%lf%lf",&p[i].x,&p[i].y);
memset(G,0,sizeof(G));
for(i=1;i<=n;i++)
{
for(j=i+1;j<=n;j++)
{
double L=Len(p[i],p[j]);
if(L<=R)
{
int flag=1;
for(int k=1;k<=n;k++)
{
if(i==k||j==k)continue;
double x1=p[k].x-p[i].x;
double y1=p[k].y-p[i].y;
double x2=p[k].x-p[j].x;
double y2=p[k].y-p[j].y;
if(x1*y2-x2*y1==0&&p[k].x>=min(p[i].x,p[j].x)&&
p[k].x<=max(p[i].x,p[j].x)&&
p[k].y>=min(p[i].y,p[j].y)&&
p[k].y<=max(p[i].y,p[j].y))
{
flag=0;
break;
}
}
if(flag)
{
edge[i].push_back(node(j));
edge[j].push_back(node(i));
G[i][j]=G[j][i]=-1;
}
}
}
}
for(i=1;i<=n;i++)
G[i][i]=(int)edge[i].size();
//put(n);
if(!ok(n))
{
printf("-1\n");
continue;
}
printf("%d\n",det(n-1));
}
}
