计算几何模板

此部分由我队友完成。浏览体验可能会有所不同因为水平比我强


有些地方无法正常浏览,可能是因为我队友有些用的LateX语法,markdown渲染不出来

基础模板及常用小函数

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
#define Vector Point
#define PI acos(-1)
/*
判断double x 的正负 1:正 0:0 -1:负
*/
int sgn(double x){
if (fabs(x) < ZERO) return 0;
else if (x > 0) return 1;
else return -1;
}

/*
* 封装点 重载 + - * / < == 运算符
*/
struct Point{
double x, y;
Point(double x = 0, double y = 0) : x(x), y(y) { }
Point operator + (Point & b) { return Point(x + b.x, y + b.y); }
Point operator - (Point & b) { return Point(x - b.x, y - b.y); }
Point operator * (double k) { return Point(x * k, y * k); }
Point operator / (double k) {return Point(x / k, y / k); }
bool operator < (Point b) { return equals(x, b.x) ? y < b.y : x < b.x; }
bool operator == (Point b) { return equals(x, b.x) && equals(y, b.y); }
};

// 计算A·B(点积)
double dot(Vector A, Vector B){
return A.x * B.x + A.y * B.y;
}
// 计算AxB(叉积)
double Cross(Vector A, Vector B){
return A.x * B.y - A.y * B.x;
}
// 计算|A|
double get_length(Point A){
return sqrt(dot(A, A));
}

// 计算A,B之间夹角, AB可以交换顺序, 结果不变
double get_angle(Point A, Point B){
return acos(dot(A, B) / get_length(A) / get_length(B));
}

// 将向量A绕原点顺时针旋转theta度
Point rotate(Point A, double theta){
return Point(A.x * cos(theta) + A.y * sin(theta), -A.x * sin(theta) + A.y * cos(theta));
}

// 点C和直线AB的关系: 判断AC和AB的叉积
// -1: C在AB左边 0: C在AB上 1: C在AB右边
int relation(Point A, Point B, Point C){
ll c = Cross(C - A, B - A);
return sgn(c);
}

// 判断点C是否在线段AB上
// 0: 点C不在线段AB上; 1: 点C在线段AB上
bool onSegment(Point A, Point B, Point C){
return relation(A, B, C) == 0 && sgn(dot(C - A, C - B)) <= 0;
}

// 返回点C到直线AB的距离
double dis2Line(Point A, Point B, Point C){
Vector v1 = B - A, v2 = C - A;
return fabs(Cross(v1, v2) / get_length(v1));
}

// 点C到线段AB的距离, 必须保证A在B点左边或者与B重合
double dis2Segment(Point A, Point B, Point C){
if (A == B) return get_length(C - A);
Vector v1 = B - A, v2 = C - A, v3 = C - B;
if (sgn(dot(v1, v2)) < 0) return get_length(v2);
if (sgn(dot(v1, v3)) < 0) return get_length(v3);
return dis2Line(A, B, C);
}

常用定理

余弦定理

c2=a2+b22×a×b×cos(C)c^2=a^2+b^2-2\times a \times b\times \cos (C)

海伦公式

p=(a+b+c)2p=\frac{(a+b+c)}2

SABC=p(pa)×(pb)×(pc)S_{\triangle ABC}=\sqrt{p(p-a)\times (p-b)\times (p-c)}

向量

常用运算

内积(点积)

AB=ABcos(θ)\overrightarrow A·\overrightarrow B = |\overrightarrow A||\overrightarrow B|\cos(\theta)

image-20220717194958500

(1) 几何意义:向量A\overrightarrow A在向量B\overrightarrow B上的投影与B\overrightarrow B的长度的乘积。

(2) 公式推导

定义向量c=ab\overrightarrow c=\overrightarrow a-\overrightarrow b

c2=a2+b22×abcosθc^2=a^2+b^2-2\times |\overrightarrow a||\overrightarrow b|\cos \theta

(ab)(ab)=a2+b22ab=a2+b22×abcosθ(\overrightarrow a -\overrightarrow b)·(\overrightarrow a - \overrightarrow b)=a^2+b^2-2\overrightarrow a·\overrightarrow b=a^2+b^2-2\times |\overrightarrow a||\overrightarrow b|\cos \theta

ab=abcosθ\overrightarrow a ·\overrightarrow b=|\overrightarrow a||\overrightarrow b|\cos \theta

(3) 代码实现

1
2
3
4
// 计算A·B(点积)
double dot(Point A, Point B){
return A.x * B.x + A.y * B.y;
}

外积(叉积)

A×B=ABsin(θ)\overrightarrow A\times \overrightarrow B = |\overrightarrow A||\overrightarrow B|\sin(\theta)

image-20220717195009039

(1) 几何意义

向量A\overrightarrow AB\overrightarrow B张成的平行四边形的有向面积。B\overrightarrow BA\overrightarrow A的逆时针方向为正。

(2) 公式推导

$S=|\alpha|\times|\beta|\times \sin(b-a)= |\alpha|\times |\beta|\times (\sin b\times \cos a-\cos b \times \sin a)\=(|\alpha|\cos a)(|\beta|\sin b)-(|\alpha|\sin a)(|\beta| \cos b)=A.x\times B.y - A.y\times B.x $

(3) 代码实现
1
2
3
4
// 计算AxB(叉积)
double Cross(Vector A, Vector B){
return A.x * B.y - A.y * B.x;
}

常用函数

1. 取模
1
2
3
4
// 计算|A|
double get_length(Point A){
return sqrt(dot(A, A));
}
2. 计算向量夹角

θ=arccos(abab)\theta=\arccos(\frac{\overrightarrow a·\overrightarrow b}{|\overrightarrow a||\overrightarrow b|})

余弦公式没啥好说的0.0

1
2
3
4
// 计算A,B之间夹角, AB可以交换顺序
double get_angle(Point A, Point B){
return acos(dot(A, B) / get_length(A) / get_length(B));
}
3. 向量

A\overrightarrow A顺时针旋转θ\theta的角度

image-20220717195028391

x2=OA×cos(αθ)=OA×(cosαcosθ+sinαsinθ)=x1×cosθ+y1×sinθx_2=OA\times \cos(\alpha -\theta)=OA\times(\cos \alpha \cos \theta+\sin\alpha \sin \theta)=x_1\times \cos \theta+y_1\times \sin \theta

y2=OA×sin(αθ)=OA×(sinαcosθcosαsinθ)=x1×sinθ+y1×cosθy_2=OA\times \sin(\alpha -\theta)=OA\times(\sin \alpha \cos \theta-\cos\alpha \sin \theta)=-x_1\times \sin \theta+y_1\times \cos \theta

1
2
3
4
// 向量A顺时针旋转theta度,theta是PI/6之类的, 不是90°,45°这样的
Point rotate(Point A, double theta){
return Point(A.x * cos(theta) + A.y * sin(theta), -A.x * sin(theta) + A.y * cos(theta));
}

点、直线

直线公式

一般式:

ax+by+c=0ax+by+c=0

点向式:

{x=x0+mty=y0+nt\begin{cases} x=x_0+mt \\ y=y_0+nt \end{cases}

设直线LL过点M0(x0,y0)M_0(x_0,y_0)s={m,n}\overrightarrow s=\{m,n\}是直线LL的向量。设MM为直线LL上任意一点,则M0M={xx0,yy0}\overrightarrow {M_0M}=\{x-x_0,y-y_0\},且M0M //  s\overrightarrow {M_0M} ~//~~ \overrightarrow s。由两向量平行的充要条件可知:

t=xx0m=yy0nt=\frac{x-x_0}m =\frac{y-y_0}n,即直线的点向式方程为:

xx0m=yy0n\frac{x-x_0}m =\frac{y-y_0}n(当m,nm,n中有一个为0时,就理解为相应的分子为0)

则方程组{x=x0+mty=y0+nt,tR\begin{cases} x=x_0+mt \\ y=y_0+nt \end{cases},t\in\R

称为直线的参数方程

斜截式

y=kx+by=kx+b

两点式

xx1x2x1=yy1y2y1\frac{x-x_1}{x_2-x_1}=\frac{y-y_1}{y_2-y_1}

常用操作与函数

点与直线

判断一个线段上有多少个整数点
1
2
3
int getPoint(Point A){
return abs(gcd(A.x, A.y)) + 1;
}
判断点和直线的关系

判断点CC和直线ABAB的关系时,可以转为判断x=AC×ABx=\overrightarrow {AC}\times \overrightarrow {AB}的值。

image-20220717195042899

1
2
3
4
5
6
// 点C和直线AB的关系: 判断AC和AB的叉积
// -1: C在AB左边 0: C在AB上 1: C在AB右边
int relation(Point A, Point B, Point C){
ll c = Cross(C - A, B - A);
return sgn(c);
}
点C是否在线段AB上

首先应满足点CC在直线ABAB上,其次满足点CC在线段ABAB间。

即:AC×AB=0,ACBC0\overrightarrow {AC}\times \overrightarrow {AB}=0,\overrightarrow{AC}·\overrightarrow {BC}\leq 0

image-20220717195618777

1
2
3
4
5
// 判断点C是否在线段AB上
// 0: 点C不在线段AB上; 1: 点C在线段AB上
bool onSegment(Point A, Point B, Point C){
return relation(A, B, C) == 0 && sgn(dot(C - A, C - B)) <= 0;
}
点到直线的距离

image-20220717195632580

SABCD=v1×v2=d×v2S_{ABCD}=\overrightarrow {v_1}\times \overrightarrow {v_2}=d\times |\overrightarrow {v_2}|, d=v1×v2v1\therefore d=\frac{\overrightarrow {v_1}\times \overrightarrow {v_2}}{|\overrightarrow {v_1}|}

1
2
3
4
5
// 返回点C到直线AB的距离
double dis2Line(Point A, Point B, Point C){
Vector v1 = B - A, v2 = C - A;
return fabs(Cross(v1, v2) / get_length(v1));
}
CC到线段ABAB的距离

image-20220717195717661

注意:A点必须在B点的左边,否则会计算错误

ABAB为一个点时,返回get_dis(A,C)get\_dis(A,C);当CC在线段ABAB间时,返回dis2Line(C,A,B)dis2Line(C,A,B)

CC在线段的左边,返回v2|\overrightarrow v_2|;当CC在线段的右边,返回v1|\overrightarrow v_1|

1
2
3
4
5
6
7
8
// 点C到线段AB的距离, 必须保证A在B点左边或者与B重合
double dis2Segment(Point A, Point B, Point C){
if (A == B) return get_length(C - A);
Vector v1 = B - A, v2 = C - A, v3 = C - B;
if (sgn(dot(v1, v2)) < 0) return get_length(v2);
if (sgn(dot(v1, v3)) < 0) return get_length(v3);
return dis2Line(A, B, C);
}
点在直线上的投影

image-20220717195732693

如图点DD的坐标就为OA+AD\overrightarrow {OA} +\overrightarrow {AD},因此,求出AD\overrightarrow {AD}即可求出投影点坐标。

AD=AC×cosθ|\overrightarrow {AD}|=|\overrightarrow {AC}|\times \cos \theta, $\overrightarrow {AB}·\overrightarrow {AC}=|\overrightarrow {AC}|\times\cos \theta\times | \overrightarrow {AB}| $

t=ADAB=ABACAB2=ABACABAB\therefore t=\frac{|\overrightarrow {AD}|}{|\overrightarrow {AB}|}=\frac{\overrightarrow {AB}·\overrightarrow {AC}}{|\overrightarrow {AB}|^2}=\frac{\overrightarrow {AB}·\overrightarrow {AC}}{\overrightarrow {AB}·\overrightarrow {AB}}

AD=OA+AB×t\overrightarrow {AD}=\overrightarrow {OA}+\overrightarrow {AB}\times t

1
2
3
4
5
6
// 计算点C到直线AB的投影
Point projection2Line(Point A, Point B, Point C){
Point v = B - A;
v = v * (dot(v, C - A) / dot(v, v)); // 咳咳, 直接写到一行会报错0.0
return A + v;
}

两线段之间的关系

两直线之间的关系:
  • 相交:A×B0\overrightarrow A\times \overrightarrow B \neq 0

  • 平行或共线:A×B=0\overrightarrow A\times \overrightarrow B = 0

快速排斥实验与跨立实验

用来判断两个线段之间的关系。

规定「一条线段的区域」为以这条线段为对角线的,各边均与某一坐标轴平行的矩形所占的区域,那么可以发现,如果两条线段没有公共区域,则这两条线段一定不相交。比如有以下两条线段,它们占用的区域是这样的:

image-20220717195746760

image-20220717195830184

于是可以快速地判断出来这两条线段不相交。

这就是 快速排斥实验。上述情况称作 未通过快速排斥实验

未通过快速排斥实验是两线段无交点的 充分不必要条件,我们还需要进一步判断。

因为两线段 相交, 线段的两个端点一定分布在 线段所在直线两端;同理, 线段的两个端点一定分布在 线段所在直线两端。我们可以直接判断一条线段的两个端点相对于另一线段所在直线的位置关系,如果不同,则两线段相交,反之则不相交。我们可以利用 3.1 中的知识帮助我们判断直线与点的位置关系。

这就是 跨立实验,如果对于两线段 , 线段的两个端点分布在 线段所在直线的两侧, 线段的两个端点分布在 线段所在直线的两侧,我们就说 两线段 通过了跨立实验,即两线段相交。

注意到当两条线段共线但不相交时也可以通过跨立实验,因此想要准确判断还需要与快速排斥实验结合。(y总你的模板又错辣,自己又重写的应该是对的)

1
2
3
4
5
6
7
// 判断线段A1A2, B1B2是否相交, 这个函数不要求A1A2端点的顺序
bool segmentIntersection(Point a1, Point a2, Point b1, Point b2){
int c1 = relation(b1, b2, a1), c2 = relation(b1, b2, a2);
int c3 = relation(a1, a2, b1), c4 = relation(a1, a2, b2);
if(c1 * c2 < 0 && c3 * c4 < 0) return true;
return onSegment(b1, b2, a1) || onSegment(b1, b2, a2);
}
两直线交点
1
2
3
4
5
6
7
8
9
10
// 求直线s1, s2的交点, 这里不要求Line端点的顺序
// 注意, 必须要满足s1,s2不共线,不平行, 存在交点
Point getCrossPoint_Line(Line s1, Line s2){
Vector v1, v2;
v1 = s1.p2 - s1.p1, v2 = s2.p2 - s2.p1;
Vector u = s1.p1 - s2.p1;
double t = Cross(v2, u) / Cross(v1, v2);
Point x = v1 * t;
return s1.p1 + x;
}
两线段交点

直接看全家桶。

多边形

三角形

1. 面积

叉积,放在多边形公式里一起写

海伦公式求面积

p=(a+b+c)2p=\frac{(a+b+c)}2

SABC=p(pa)×(pb)×(pc)S_{\triangle ABC}=\sqrt{p(p-a)\times (p-b)\times (p-c)}

2. 外心, 内心, 垂心, 重心

(1) 外心,外接圆圆心
三边中垂线交点。到三角形三个顶点的距离相等
(2) 内心,内切圆圆心
角平分线交点,到三边距离相等
(3) 垂心
三条垂线交点
(4) 重心
三条中线交点(到三角形三顶点距离的平方和最小的点,三角形内到三边距离之积最大的点)

普通多边形

通常按逆时针存储所有点

1. 相关定义

(1) 多边形:由在同一平面且不再同一直线上的多条线段首尾顺次连接且不相交所组成的图形叫多边形

(2) 简单多边形:简单多边形是除相邻边外其它边不相交的多边形

(3) 凸多边形:过多边形的任意一边做一条直线,如果其他各个顶点都在这条直线的同侧,则把这个多边形叫做凸多边形。
任意凸多边形外角和均为360°360°
任意凸多边形内角和为(n2)×180°(n−2)\times 180°

2. 常用操作与函数

求多边形面积(不一定是凸多边形,凹多边形也可以)

1
2
3
4
5
6
7
8
9
10
// 计算多边形s的面积, s中的点不一定要按照逆时针存储
typedef vector<Point> Polygon;
double get_area(Polygon s){
int n = s.size();
double res = 0;
for (int i=0;i<n;++i){
res += fabs(Cross(s[i], s[(i + 1) % n]) / 2);
}
return res;
}

判断点与多边形(不一定是凸多边形)的关系

射线法:从该点任意做一条和所有边都不平行的射线。交点个数为偶数,则在多边形外,为奇数,则在多边形内。

转角法

判断点在凸多边形内:只需判断点是否在所有边的左边(逆时针存储多边形)

3. 皮克定理

皮克定理指的是一个计算点阵中顶点在格点上的多边形面积公式。该公式可表示为S=a+b21S=a+\frac b 2-1,其中aa表示多边形内部的点数,bb表示多边形落在格点边界上点数,SS表示多边形的面积。

证明:不会证明,略。

例题hdu1706

image-20220717195844971

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
#define ll long long
int a1, b1, a2, b2, a3, b3;

struct Point{
int x, y;
Point(int x = 0, int y = 0) : x(x), y(y) { }
Point operator + (Point & b) { return Point(x + b.x, y + b.y); }
Point operator - (Point & b) { return Point(x - b.x, y - b.y); }
Point operator * (int k) { return Point(x * k, y * k); }
Point operator / (int k) {return Point(x / k, y / k); }
bool operator < (Point b) { return x == b.x ? y < b.y : x < b.x; }
};

int gcd(int a, int b){
return b ? gcd(b, a % b) : a;
}

int Cross(Point A, Point B){
return abs(A.x * B.y - A.y * B.x);
}

int getPoint(Point A){
return abs(gcd(A.x, A.y)) + 1;
}

void solve(){
Point A = Point(a1, b1), B = Point(a2, b2), C = Point(a3, b3);
int area = Cross(B - A, C - A);
int ab = getPoint(B - A), bc = getPoint(B - C), ac = getPoint(C - A);
int b = ab + bc + ac - 3;
int res = (area + 2 - b) / 2;
printf("%d\n", res);
}

int main(void){
while (~scanf("%d%d%d%d%d%d", &a1, &b1, &a2, &b2, &a3, &b3)){
if (a1 == 0 && b1 == 0 && a2 == 0 && b2 == 0) return 0;
solve();
}

return 0;
}

累了,看计算几何全家桶吧orz。

直线与圆、圆与圆的位置关系

1. 直线与圆的位置关系
2. 圆与圆的位置关系
3. 直线到圆的切线

公切线问题

内切圆与外切圆

极角排序

介绍

在平面内取一个定点OO,叫极点,引一条射线OxOx,叫做极轴,再选定一个长度单位和角度的正方向(通常取逆时针方向)。对于平面内任何一点MM,用ρρ表示线段OMOM的长度(有时也用rr表示),θθ表示从OxOxOMOM的角度,ρρ叫做点MM的极径,θθ叫做点MM的极角,有序数对 (ρ,θ)(ρ,θ)就叫点MM的极坐标。

那么给定平面上的一些点,把它们按照一个选定的中心点排成顺(逆)时针。

三种常用方法

1. atan2()函数

image-20220717195902247

1
2
3
4
5
6
7
8
using Points = vector<Point>;
double theta(auto p) { return atan2(p.y, p.x); } // 求极角
void psort(Points &ps, Point c = O) // 极角排序
{
sort(ps.begin(), ps.end(), [&](auto p1, auto p2) {
return lt(theta(p1 - c), theta(p2 - c));
});
}
2.叉积

image-20220717195912242

1
2
3
4
5
6
7
int qua(auto p) { return lt(p.y, 0) << 1 | lt(p.x, 0) ^ lt(p.y, 0); }    // 求象限
void psort(Points &ps, Point c = O) // 极角排序
{
sort(ps.begin(), ps.end(), [&](auto v1, auto v2) {
return qua(v1 - c) < qua(v2 - c) || qua(v1 - c) == qua(v2 - c) && lt(cross(v1 - c, v2 - c), 0);
});
}

我们用0、1、2、3来表示第一、二、三、四象限。这种方法常数可能稍微大一点,但是精度比较好,如果坐标都是整数的话是完全没有精度损失的。另外一定要注意,在角度可能很小时,不要用比较余弦的方法比较角度,那样精度会有严重问题。

应用与例题

二维凸包

这篇博客已经介绍的非常详细了,但是太太太长了,放链接

https://www.luogu.com.cn/blog/ShineEternal/convex-hull

凸包是围住所有点的面积最小,周长也最小的凸多边形。

算法步骤

用一个单调栈来维护凸包形状。

  1. 将点排序,x为第一关键字,y为第二关键字。

  2. 分两步,从左至右维护上半部分(上链);从右至左维护下半部分(下链)。

    维护上链时,首先加入前两个点;维护下链时,首先加入最后两个点。每次看新加进来的点CC是否在向量AB\overrightarrow{AB}的哪一侧。如果是右侧,则将BCBC加入栈中;如果是左侧,则不断弹出点BB,最后加入AC\overrightarrow {AC}如果C在AB\overrightarrow{AB}上,则要根据题目要求确定要不要加入点CC

  3. 注意,在最后一步一定要再用栈底元素更新一遍,否则无法构成凸包。更新这两个边界点非常重要!

image-20220717195922949

image-20220717195946550

代码

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
inline Polygon andrewScan(Polygon p)//求凸包
{
Polygon u, d;
if (p.size() < 3) return p;
sort(p.begin(), p.end());

u.pb(p[0]);
u.pb(p[1]);

d.pb(p[p.size() - 1]);
d.pb(p[p.size() - 2]);

for (int i = 2; i < p.size(); i++)
{
for (int k = u.size(); k >= 2 && ccw(u[k - 2], u[k - 1], p[i]) != C; k--) u.pop_back();
u.pb(p[i]);
}

for (int i = p.size() - 3; ~i; i--)
{
for (int k = d.size(); k >= 2 && ccw(d[k - 2], d[k - 1], p[i]) != C; k--) d.pop_back();
d.pb(p[i]);
}

reverse(d.begin(), d.end());
for (int i = u.size() - 2; i; i--) d.pb(u[i]);
return d;
}

int main(void){
scanf("%d", &n);
Polygon q(n);
for (int i = 0; i < n; i ++ ) scanf("%lf%lf", &q[i].x, &q[i].y);
Polygon d = andrewScan(q);
double res = get_dis(d[0], d.back());
for (int i=1;i<d.size();++i){
res += get_dis(d[i], d[i - 1]);
}

printf("%.2f\n", res);

return 0;
}

最小圆覆盖

image-20220717200003691

image-20220717200011926

image-20220717200026512

image-20220717200048895

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
#include <iostream>
#include <cstdio>
#include <algorithm>
#include <cstring>
#include <cmath>
#define pf(x) (x) * (x)
#define eps 1e-12

using namespace std;

int n;
struct point{
double x;
double y;
}o, a[100000 + 10];

double r;

double dis(point a, point b){
return sqrt(pf(a.x - b.x) + pf(a.y - b.y));
}

void get(point a, point b, point c){ //方程的计算
double a11 = b.x - a.x;
double a12 = b.y - a.y;
double b1 = (pf(b.x) - pf(a.x) + pf(b.y) - pf(a.y)) * 0.5;
double a21 = c.x - a.x;
double a22 = c.y - a.y;
double b2 = (pf(c.x) - pf(a.x) + pf(c.y) - pf(a.y)) * 0.5;

o.x = (b1 * a22 - a12 * b2) / (a11 * a22 - a12 * a21);
o.y = (a11 * b2 - b1 * a21) / (a11 * a22 - a12 * a21);
r = dis(o, a);
}

int main(){
cin >> n;
for(int i = 1; i <= n; i++){
cin >> a[i].x >> a[i].y;
}
random_shuffle(a + 1, a + 1 + n);
o = a[1];
r = 0;
for(int i = 2; i <= n; i++){ //模拟那个过程
if(dis(a[i], o) > r + eps){
o = a[i];
r = 0;
for(int j = 1; j <= i - 1; j++){
if(dis(a[j], o) > r + eps){
o.x = (a[i].x + a[j].x) / 2;
o.y = (a[i].y + a[j].y) / 2;
r = dis(o, a[j]);
for(int k = 1; k <= j - 1; k++){
if(dis(a[k], o) > r + eps){
get(a[i], a[j], a[k]);
}
}
}
}
}
}
printf("%.10lf\n%.10lf %.10lf", r, o.x, o.y);
return 0;
}

半平面交

image-20220717200059604

增量法构造半平面

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
114
115
116
117
118
119
120
121
122
123
#include <bits/stdc++.h>
using namespace std;
const double INF = 1e12;
const double pi = acos(-1.0);
const double eps = 1e-8;
int sgn(double x)
{
if(fabs(x) < eps) return 0;
else return x<0?-1:1;
}
struct Point
{
double x,y;
Point() {}
Point(double x,double y):x(x),y(y) {}
Point operator + (Point B)
{
return Point(x+B.x,y+B.y);
}
Point operator - (Point B)
{
return Point(x-B.x,y-B.y);
}
Point operator * (double k)
{
return Point(x*k,y*k);
}
};
typedef Point Vector;
double Cross(Vector A,Vector B)
{
return A.x*B.y - A.y*B.x; //叉积
}
struct Line
{
Point p;
Vector v;
double ang;
Line() {};
Line(Point p,Vector v):p(p),v(v)
{
ang=atan2(v.y,v.x);
}
bool operator < (Line &L)
{
return ang<L.ang; //用于极角排序
}
};
//点p在线L左边,即点p在线L在外面:
bool OnLeft(Line L,Point p)
{
return sgn(Cross(L.v,p-L.p))>0;
}
Point Cross_point(Line a,Line b) //两直线交点
{
Vector u=a.p-b.p;
double t=Cross(b.v,u)/Cross(a.v,b.v);
return a.p+a.v*t;
}
vector<Point> HPI(vector<Line> L) //求半平面交,返回图多边形
{
int n=L.size();
sort(L.begin(),L.end());//将所有半平面按照极角排序。
int first,last; //指向双端队列的第一个和最后一个元素
vector<Point> p(n); //两个相邻半平面的交点
vector<Line> q(n); //双端队列
vector<Point> ans; //半平面交形成的凸包
q[first=last=0]=L[0];
for(int i=1; i<n; i++)
{
//情况1:删除尾部的半平面
while(first<last && !OnLeft(L[i], p[last-1])) last--;
//情况2:删除首部的半平面:
while(first<last && !OnLeft(L[i], p[first])) first++;
q[++last]=L[i]; //将当前的半平面加入双端队列尾部
//极角相同的两个半平面,保留左边:
if(fabs(Cross(q[last].v,q[last-1].v)) < eps)
{
last--;
if(OnLeft(q[last],L[i].p)) q[last]=L[i];
}
//计算队列尾部半平面交点:
if(first<last) p[last-1]=Cross_point(q[last-1],q[last]);
}
//情况3:删除队列尾部的无用半平面
while(first<last && !OnLeft(q[first],p[last-1])) last--;
if(last-first<=1) return ans; //空集
p[last]=Cross_point(q[last],q[first]); //计算队列首尾部的交点。
for(int i=first; i<=last; i++) ans.push_back(p[i]); //复制。
return ans; //返回凸多边形
}
double Polygon_area(vector<Point> p) //Point *p表示多边形。从原点开始划分三角形
{
int n=p.size();
double area = 0;
for(int i = 0; i < n; i++)
area += Cross(p[i],p[(i+1)%n]);
return area/2; //面积有正负,不能简单地取绝对值
}
int main()
{
int n,m;
vector<Line> L;
cin>>n;
while(n--)
{
cin>>m;
vector<Point> tmp;
for(int i=0;i<m;i++)
{
double a,b;
scanf("%lf%lf",&a,&b);
tmp.push_back(Point(a,b));
}
for(int i=0;i<m;i++)
{
L.push_back(Line(tmp[i],tmp[(i+1)%m]-tmp[i]));
}
}
vector<Point> ans=HPI(L); //得到凸多边形
printf("%.3f",Polygon_area(ans));
return 0;
}

博弈论

image-20220717200109055

image-20220717200118918

image-20220717200128163

image-20220717200136711

image-20220717200149351

image-20220717200154574

image-20220717200200890

image-20220717200206395

image-20220717200211219

image-20220717200215704

博弈论模板

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
#include <bits/stdc++.h>
#define FOR(i,a,b) for(register int i=(a);i<(b);++i)
#define ROF(i,a,b) for(register int i=(a);i>=(b);--i)
#define pi pair<int,int>
#define mk(a,b) make_pair(a,b)
#define mygc(c) (c)=getchar()
#define mypc(c) putchar(c)
#define fi first
#define se second
using namespace std;
typedef long long ll;
typedef double db;
const int maxn = 10000;
const int maxm = 100;
const int inf = 2147483647;
typedef long long ll;
const double eps = 1e-9;
const long long INF = 9223372036854775807ll;
ll qpow(ll a,ll b,ll c){ll ans=1;while(b){if(b&1)ans=ans*a%c;a=a*a%c;b>>=1;}return ans;}
inline void rd(int *x){int k,m=0;*x=0;for(;;){mygc(k);if(k=='-'){m=1;break;}if('0'<=k&&k<='9'){*x=k-'0';break;}}for(;;){mygc(k);if(k<'0'||k>'9')break;*x=(*x)*10+k-'0';}if(m)(*x)=-(*x);}
inline void rd(ll *x){int k,m=0;*x=0;for(;;){mygc(k);if(k=='-'){m=1;break;}if('0'<=k&&k<='9'){*x=k-'0';break;}}for(;;){mygc(k);if(k<'0'||k>'9')break;*x=(*x)*10+k-'0';}if(m)(*x)=-(*x);}
inline void rd(db *x){scanf("%lf",x);}
inline int rd(char c[]){int i,s=0;for(;;){mygc(i);if(i!=' '&&i!='\n'&&i!='\r'&&i!='\t'&&i!=EOF) break;}c[s++]=i;for(;;){mygc(i);if(i==' '||i=='\n'||i=='\r'||i=='\t'||i==EOF) break;c[s++]=i;}c[s]='\0';return s;}
inline void rd(int a[],int n){FOR(i,0,n)rd(&a[i]);}
inline void rd(ll a[],int n){FOR(i,0,n)rd(&a[i]);}
template <class T, class S> inline void rd(T *x, S *y){rd(x);rd(y);}
template <class T, class S, class U> inline void rd(T *x, S *y, U *z){rd(x);rd(y);rd(z);}
template <class T, class S, class U, class V> inline void rd(T *x, S *y, U *z, V *w){rd(x);rd(y);rd(z);rd(w);}
inline void wr(int x){if(x < 10) putchar('0' + x); else wr(x / 10), wr(x % 10);}
inline void wr(int x, char c){int s=0,m=0;char f[10];if(x<0)m=1,x=-x;while(x)f[s++]=x%10,x/=10;if(!s)f[s++]=0;if(m)mypc('-');while(s--)mypc(f[s]+'0');mypc(c);}
inline void wr(ll x, char c){int s=0,m=0;char f[20];if(x<0)m=1,x=-x;while(x)f[s++]=x%10,x/=10;if(!s)f[s++]=0;if(m)mypc('-');while(s--)mypc(f[s]+'0');mypc(c);}
inline void wr(db x, char c){printf("%.15f",x);mypc(c);}
inline void wr(const char c[]){int i;for(i=0;c[i]!='\0';i++)mypc(c[i]);}
inline void wr(const char x[], char c){int i;for(i=0;x[i]!='\0';i++)mypc(x[i]);mypc(c);}
template<class T> inline void wrn(T x){wr(x,'\n');}
template<class T, class S> inline void wrn(T x, S y){wr(x,' ');wr(y,'\n');}
template<class T, class S, class U> inline void wrn(T x, S y, U z){wr(x,' ');wr(y,' ');wr(z,'\n');}
template<class T> inline void wra(T x[], int n){int i;if(!n){mypc('\n');return;}FOR(i,0,n-1)wr(x[i],' ');wr(x[n-1],'\n');}
vector<int>f;
int sg[maxn],vis[maxn];
void getSG(int n,int maxx){//打表sg[0],sg[1]...sg[n],maxx代表sg值的上限
sort(f.begin(),f.end());
memset(sg,0,sizeof(int)*(n+1));
FOR(i,1,n+1){
FOR(j,0,f.size()){
if(f[j]>i)break;
vis[sg[i-f[j]]]=1;
}
FOR(j,0,maxx+1)if(!vis[j]){
sg[i]=j;
break;
}
FOR(j,0,f.size()){
if(f[j]>i)break;
vis[sg[i-f[j]]]=0;
}
}
}
int main(){
//测试
f.push_back(1);
f.push_back(5);
getSG(100,50);
wra(sg,100);
}

image-20220717200227937

image-20220717200233308

image-20220717200237402

image-20220717200242043

image-20220717200247752

image-20220717200252794

image-20220717200302198