抱歉,您的浏览器无法访问本站

本页面需要浏览器支持(启用)JavaScript


了解详情 >

Description

平面内 N 个坐标均为整数的点,其中任意三点不共线,令所有点 P_i 构成一个全集 U ,子集 S\subseteq U |S|\ge 3 。问凸包面积为整数的子集 S 的个数模 10^9 + 7 的值。

3\le N\le 80

Solution

比较神奇的题目。

考虑 Andrew 算法的执行流程:将所有点按照先 x 轴升序后 y 轴升序排序,然后从左往右弄出下凸包,再从右往左弄出上凸包。本题可以利用类似的思路。

然后我们由叉积计算面积的方法知道,整点三角形的面积乘二之后必然是整数。于是我们考虑依次用三角形弄出凸包,具体地,可以 dp:

  • O(n) 枚举凸包最左端的端点 P_{s}

  • \mathrm{upper}_{i,jk} :凸包满足如下条件的点集 S 的数量:

    • P_s 顺时针考虑每个点,最后两个分别是 P_i P_j
    • 所有的点都在直线 P_sP_j 以上;
    • 面积乘以二之后模二为 k
  • \mathrm{lower}_{i,j,k} :凸包满足如下条件的点集 S 的数量:

    • P_s 逆时针考虑每个点,最后两个分别是 P_i P_j
    • 所有的点都在直线 P_sP_j 以下;
    • 面积乘以二之后模二为 k
  • 确定每个最右端的点 P_j ,然后枚举 i k 统计答案:

    \mathrm{ans}_s = \sum_{j = s + 1}^n\sum_{k \in \{0, 1\}}\left(\sum_{i}\mathrm{upper}_{i,j,k} \times \sum_{i}\mathrm{lower_{i,j,k}} \right)

这样做的原理其实就是合并上下凸壳:根据上面设计的 dp 状态可知我们确定了左端点 P_s 和右端点 P_j 后,满足上凸壳的点集个数即为枚举倒数第二个点 i 然后求和的结果 \sum_{i}\mathrm{upper}_{i,j,k} ,满足下凸壳的点集个数即为 \sum_i\mathrm{lower}_{i,j,k}

怎么计算 \mathrm{upper} \mathrm{lower} 呢?我们可以通过 \mathrm{upper}_{i,j,k} 递推出 \mathrm{upper}_{j,l,k} ,其原理就是判断 P_l 能否被加入凸壳中,能则更新,不能就更新 \mathrm{lower}_{j,l,k} (加入下凸壳)。新加入的部分即为 \triangle P_sP_jP_l ,需要根据这个面积的二倍的奇偶性判断更新哪个 k

这个做法要求我们提前预处理每个三角形的面积的二倍的奇偶性 \mathrm{parity}_{i,j,k} 和三角形内点的个数 \mathrm{inside}_{i,j,k} 。总时间复杂度 O(n^4)

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
#include <cstdio>
#include <cctype>
#include <cmath>
#include <algorithm>
#define il inline
#define FOR(i, a, b) for (int i = a; i <= b; ++i)
#define DEC(i, a, b) for (int i = a; i >= b; --i)

const int maxn = 85;

int read()
{
int s = 0, x = 0;
char c = getchar();
while (!isdigit(c))
x |= (c == '-'), c = getchar();
while (isdigit(c))
s = s * 10 + c - '0', c = getchar();
return x ? -s : s;
}

struct Point
{
int x, y;
Point() {}
Point(int _x, int _y) {x = _x, y = _y;}
il bool operator<(const Point &b) const {return x == b.x ? y < b.y : x < b.x;}
} P[maxn];

typedef Point Vector;

struct modint
{
typedef long long ll;
const static ll mod = 1e9 + 7;
ll val;
modint(int _val = 0) {val = _val;}
modint operator+(const modint &b) const {return modint((val + b.val) % mod);}
modint operator-(const modint &b) const {return modint((val - b.val + mod) % mod);}
modint operator*(const modint &b) const {return modint((val * b.val) % mod);}
modint operator+=(const modint &b) {return *this = *this + b;}
modint operator*=(const modint &b) {return *this = *this * b;}
modint operator-=(const modint &b) {return *this = *this - b;}
modint operator+(const int &b) const {return modint((val + b) % mod);}
modint operator-(const int &b) const {return modint((val - b + mod) % mod);}
modint operator*(const int &b) const {return modint((val * b) % mod);}
modint operator+=(const int &b) {return *this = *this + b;}
modint operator-=(const int &b) {return *this = *this - b;}
modint operator*=(const int &b) {return *this = *this * b;}
};

il Vector operator-(const Point &a, const Point &b) {return Point(b.x - a.x, b.y - a.y);}
il int operator^(const Vector &a, const Vector &b) {return a.x * b.y - a.y * b.x;}//叉积

int area2(Point x, Point y, Point z)//返回面积的二倍
{
return abs((y - x) ^ (z - x));
}

int n, parity[maxn][maxn][maxn], inside[maxn][maxn][maxn];
modint pow2[maxn], upper[maxn][maxn][2], lower[maxn][maxn][2];

int main()
{
n = read();
FOR(i, 1, n) P[i].x = read(), P[i].y = read();
std::sort(P + 1, P + n + 1);//对点进行排序

FOR(i, 1, n)
FOR(j, 1, n)
FOR(k, 1, n)
{
if (i == j || j == k || i == k) continue;
parity[i][j][k] = area2(P[i], P[j], P[k]) % 2;//更新面积奇偶性
FOR(l, 1, n)//枚举每个 l 点判断 l 在不在三角形内
if (l != i && l != j && l != k && area2(P[i], P[j], P[k]) == area2(P[l], P[j], P[k]) + area2(P[i], P[l], P[k]) + area2(P[i], P[j], P[l]))
inside[i][j][k] += 1;
}


pow2[0] = 1;
FOR(i, 1, n) pow2[i] = pow2[i - 1] * 2;//预处理二的次幂

modint ans = 0;

DEC(must, n, 1)//枚举最左边的端点
{
FOR(i, must, n)
FOR(j, must, n)
FOR(k, 0, 1)
upper[i][j][k] = lower[i][j][k] = 0;//清空数组
FOR(i, must + 1, n) upper[must][i][0] = lower[must][i][0] = 1;//边界条件
FOR(i, must, n)
FOR(j, must + 1, n)
FOR(k, 0, 1)
FOR(l, j + 1, n)
if (((P[l] - P[j]) ^ (P[j] - P[i])) > 0)//加入上凸壳
upper[j][l][k ^ parity[must][j][l]] += upper[i][j][k] * pow2[inside[must][j][l]];//乘法原理
else lower[j][l][k ^ parity[must][j][l]] += lower[i][j][k] * pow2[inside[must][j][l]];
FOR(j, must + 1, n)//枚举最右端点
FOR(k, 0, 1)
{
modint up = 0, lo = 0;//上下凸壳
FOR(i, must, j - 1)
up += upper[i][j][k], lo += lower[i][j][k];
ans += up * lo;
}
}

printf("%d\n", (ans - n * (n - 1) / 2).val);//最后要去除线段的贡献
return 0;
}

评论




Blog content follows the [Attribution-NonCommercial-ShareAlike 4.0 International (CC BY-NC-SA 4.0) License](https://creativecommons.org/licenses/by-nc-sa/4.0/deed.en)
本站总访问量为 访客数为
Use Volantis as theme