HDU 5581 Infinity Point Sets ACM/ICPC 2015 上海区域赛 I 计算几何+组合计数

2015年上海区域赛的题目，这道题还是比较有趣的，反正我是WA哭了。。

题目在HDU上也有，链接：

http://acm.split.hdu.edu.cn/showproblem.php?pid=5581

题目的意思是，给出二维空间里n个点的坐标，求有多少个不同的子点集不是无限点集。无限点集的定义是，将点集中的点两两相连，线段产生的交点加入点集中，继续上面的操作，如果操作能够无限地进行下去，则称之为无限点集。

 

首先我们可以分析，只有4个点或更少的点集一定不是无穷点集，如图：

 


其次，5个点及以上的点集大概有以下三种情况：

 


在第一种情况中，五个蓝色的点交点为五个新的红色的点，这样的操作显然能够无限地进行下去。

在第二种情况中，是有三个以上的点共线，还有两个点分别在线段的两边，这样一来生成的新的黑色的点也在线段上，操作不会无限地进行。

第三种情况中，是四个以上的点共线，还有一个点在外面，这种情况不会产生交点，于是操作也不会无限进行。

当然还有第四种情况，所以的点都共线。

 

但是在第二种情况中，会出现重复计算的状况，比如下图：

 


这种情况下，分别以AOC共线，BD在上下两边 和BOD共线，AC在左右两边，将第二种情况重复计算了一次。

 

综上所述，答案由以下几部分组成：

1.四个点及以下

2.五个及以上的共线的点

3.四个及以上共线的点+线外一个点

4.三个及以上共线的点+线外两侧各一个点

5.减去重复计算的第四部分

 
那么实现方法就是，1直接用排列组合做了。
枚举每一个点i，对所以其他的点做对点i的极角排序，这样就能找出共线的点，这种情况下情况2就能计算出来。在极角排序之后同时统计直线两侧点的个数，就能够计算出3和4。最后统计以点i为中心，每条直线两个方向上各有多少个点，然后枚举两两直线，将重复的情况4给计算出来，然后减掉。
 
这样这个问题就完全解决了。对于更具体的实现，计算直线两侧的点的个数 以及 直线上两边点的个数，直接枚举大于等于0的极角，然后二分查找相反方向的角，就能很快地处理出来这些数据。

最后的结果要对10^9+7取膜，于是还有预处理出1-1000的逆元，在计算排列的时候除要用乘逆元来代替。
 
AC代码：

#include <iostream>
#include <cmath>
#include <cstdio>
#include <cstdlib>
#include <cstring>
#define MAXN 1005
#define MOD 1000000007
#define LL long long
#define PI (acos(-1))
#define EPS 1e-10
using namespace std;

struct Point {
double x, y;
Point(double x = 0, double y = 0):x(x), y(y) {}
};

typedef Point Vector;

Vector operator - (Point A, Point B) { return Vector(A.x - B.x, A.y - B.y); }

LL C[MAXN], n;

inline int sign(double x) {
if (fabs(x) < EPS) return 0;
else return x < 0 ? -1 : 1;
}

inline double Dot(Vector A, Vector B) { return A.x*B.x + A.y*B.y; }//????
inline double Length(Vector A) { return sqrt(Dot(A, A)); }
inline double Cross(Vector A, Vector B) { return A.x*B.y - A.y*B.x; }
inline double Angle(Vector A, Vector B) {
double res;
res = acos(Dot(A, B) / Length(A) / Length(B));
if (sign(PI - res) == 0 || Cross(B, A) >= 0) return res;
else return -res;
}

inline LL Power(LL a, LL b) {
LL ans = 1, now = a;
if (a == 0) return 0;
while (b != 0) {
if (b & 1) {
ans *= now;
ans %= MOD;
}
now *= now;
now %= MOD;
b = b >> 1;
}
return ans;
}

inline LL c(LL a, LL b) {
LL ans = 1;
b = b < a - b ? b : a - b;
for (int i = 0; i < b; i++) {
ans *= a - i;
ans %= MOD;
ans *= C[i + 1];
ans %= MOD;
}
return ans;
}

Point node[MAXN];
double beta[MAXN];
int s1[MAXN], e1[MAXN], s2[MAXN], e2[MAXN], len;
int tot, tot2;
Vector vx;

inline int comp(const void* a, const void *b) {
return sign((*(double*)a) - (*(double*)b));
}

inline void Find(double v) {
int l, r, mid;
l = 0; r = tot;
while (l != r) {
mid = (l + r) / 2;
if (sign(beta[mid] - v) >= 0) r = mid;
else l = mid + 1;
}
s2[tot2] = l;
l = 0; r = tot;
while (l != r) {
mid = (l + r) / 2;
if (sign(beta[mid] - v) > 0) r = mid;
else l = mid + 1;
}
e2[tot2] = l;
}

int main() {
int T, Case, i, j, k;
int numu, numd;
LL ans, sum, sum1;

vx.x = 1;
vx.y = 0;

//预处理出1-MAXN的逆元 快速幂
for (i = 1; i < MAXN; i++)
C[i] = Power(i, MOD - 2);

cin >> T;
for (Case = 1; Case <= T; Case++) {
cin >> n;
for (i = 0; i < n; i++)
scanf("%lf%lf", &node[i].x, &node[i].y);
ans = 0;
//1 2 3 4个点都满足条件，故(1,n) (2,n) (3,n) (4,n)必然存在
for (i = 1; i <= 4; i++) {
ans += c(n, i);
ans %= MOD;
}
//枚举每一个点i
for (i = 0; i < n; i++) {
tot = 0; tot2 = 0;
//求出所有点对i点的极角
for (j = 0; j < n; j++) if (j != i)
beta[tot++] = Angle(node[j] - node[i], vx);
//按极角排序
qsort(beta, tot, sizeof(double), comp);

//s1 s2真实存在 e1 e2不存在

for (j = 0; j < tot; j++) {
if (beta[j] < 0 || sign(PI - beta[j]) == 0) continue;
if (sign(beta[j]) >= 0) {

//s1表示 beta[j] 开始的位置；e1表示 beta[j] 结束的位置
//s2表示 PI-beta[j] 即beta[j]的反方向；开始的位置 e2同理

s1[tot2] = j;
while (j + 1 != tot && sign(beta[j + 1] - beta[j]) == 0) j++;
e1[tot2] = j + 1;
len = e1[tot2] - s1[tot2];

//Find用于二分查找s2 e2的位置
if (sign(beta[j]) == 0) Find(PI);
else Find(beta[j] - PI);

if (len >= 2) {
//3 + 2

if (sign(beta[j]) == 0) {
numu = s2[tot2] - e1[tot2];//numu表示直线上方点的个数
numd = s1[tot2];//numd表示直线下方点的个数
} else {
numu = s1[tot2] - e2[tot2];
numd = s2[tot2] + tot - e1[tot2];
}

//上方下方都有点才会出现 3+2
if (numu != 0 && numd != 0) {
sum = c(len, 2);
sum %= MOD;
sum1 = sum;
for (k = 3; k <= len; k++) {
sum *= len - k + 1;
sum %= MOD;
sum *= C[k];
sum %= MOD;
sum1 += sum;
sum1 %= MOD;
}
sum1 *= (LL)numu*(LL)numd;
sum1 %= MOD;
ans += sum1;
ans %= MOD;
}

if (len >= 3) {
//4 + 1
//上方或者下方有点才会出现4+1
if (numu != 0 || numd != 0) {
sum = c(len, 3);
sum %= MOD;
sum1 = sum;
for (k = 4; k <= len; k++) {
sum *= len - k + 1;
sum %= MOD;
sum *= C[k];
sum %= MOD;
sum1 += sum;
sum1 %= MOD;
}
sum1 *= (LL)numu + (LL)numd;
sum1 %= MOD;
ans += sum1;
ans %= MOD;
}

if (len >= 4) {
//5
sum = c(len, 4);
sum %= MOD;
sum1 = sum;
for (k = 5; k <= len; k++) {
sum *= len - k + 1;
sum %= MOD;
sum *= C[k];
sum %= MOD;
sum1 += sum;
sum1 %= MOD;
}
ans += sum1;
ans %= MOD;
}
}
}
if (e2[tot2] - s2[tot2] != 0) tot2++;//如果反方向也有点才记录直线
}
}

//枚举以i为中心两两直线，减去重复计算的次数
for (j = 0; j < tot2; j++)
for (k = j + 1; k < tot2; k++) {
//四个方向的点的数量相乘
sum = 1;
sum *= e2[k] - s2[k];
sum %= MOD;
sum *= e2[j] - s2[j];
sum %= MOD;
sum *= e1[k] - s1[k];
sum %= MOD;
sum *= e1[j] - s1[j];
sum %= MOD;
ans = (ans + MOD - sum) % MOD;
}
}
cout << "Case #" << Case << ": " << ans << endl;
}
return 0;
}