计算几何基础

基础知识

极坐标

在平面内取一个定点$O$,叫极点,引一条射线$Ox$,叫做极轴,再选定一个长度单位和角度的正方向(通常取逆时针方向)。对于平面内任何一点$M$,用$ρ$表示线段$OM$的长度(有时也用$r$表示),$θ$表示从$Ox$到$OM$的角度,$ρ$叫做点$M$的极径,$θ$叫做点$M$的极角,有序数对$(ρ,θ)$就叫点$M$的极坐标,这样建立的坐标系叫做极坐标系。通常情况下,$M$的极径坐标单位为$1$(长度单位),极角坐标单位为$rad$(或°)。

与直角坐标的变换

$r = \sqrt{y^2 + x^2}$

$\theta = atan2(y,x)$

atan2是已将象限纳入考量的反正切函数,单位是弧度。

三角形

余弦公式

$c = \sqrt{a^2 + b^2 - 2ab\times cos(C)}$

$cos(C) = \frac{a^2 + b^2 - c^2}{2ab}$

正弦公式

$\frac{a}{sin(A)} = \frac{b}{sin(B)} = \frac{c}{sin(C)} = 2R = D$

向量的减法

$\vec{c} = \vec{a} + (-\vec{b})$

数量积(点乘)

$\vec{a} \cdot \vec{b} = |\vec{a}| \cdot |\vec{b}| \cdot cos(\theta) = a_x \cdot b_x + a_y \cdot b_y$

$0 \leq \theta \leq \pi $

作用:

  1. 可以用来求两向量的夹角。
  2. $\vec{a} \cdot \vec{b} < 0$ 夹角在$[90,180)$度内
  3. $\vec{a} \cdot \vec{b} = 0$ 两向量垂直
  4. $\vec{a} \cdot \vec{b} > 0$ 夹角在$[0,90)$度内

向量积(叉积)

简单地说,两个向量$\vec{v}$和$\vec{w}$的叉积等于$\vec{v}$和$\vec{w}$组成的三角形的有向面积的两倍。

如果$\vec{w}$在$\vec{v}$的左边,那么$cross(v,w) > 0$

$|\vec{a} \times \vec{b}| = x_a \times y_b - x_b \times y_a$

$\vec{c} = \vec{a} \times \vec{b}$

$|\vec{c}| = |\vec{a}| \cdot |\vec{b}| \cdot sin(\theta)$

$0 \leq \theta \leq \pi $

其中$\vec{c}$垂直于$\vec{a},\vec{b}$所组成的平面,满足右手定则。

几何意义: 以$\vec{a},\vec{b}$为棱的平行四边形面积。

那么多边形面积就可以通过叉积求出。

一般意义上,竖直向上的符号为正。

$
\vec{a} \times \vec{b} = \left| \begin{array}{ccc}
\vec{i} & \vec{j} & \vec{k}\\
a_x & a_y & a_z\\
b_x & b_y & b_z\\
\end{array} \right|
= (a_y \cdot b_z - a_z \cdot b_y) \cdot \vec{i} + (-a_x \cdot b_z + a_z \cdot b_x) \cdot \vec{j} + (a_x \cdot b_y - a_y \cdot b_x) \cdot \vec{k}
$

混合积

$
\vec{a} \cdot (\vec{b} \times \vec{c}) = \left| \begin{array}{ccc}
a_x & a_y & a_z\\
b_x & b_y & b_z\\
c_x & c_y & c_z\\
\end{array} \right|
$

证明:

几何意义为: 以$\vec{a},\vec{b},\vec{c}$为棱的平行六面体的体积。

Pick定理

设以整数点为顶点的多边形的面积为S,多边形内部的整数点数为N,多边形边界上的整数点为L,
则 S = L/2 + N - 1.

重要定理

三点共线定理

已知O是AB所在直线外一点,若$\vec{OC} = \lambda \cdot \vec{OA} + \mu \cdot \vec{OB}$, 且$\lambda + \mu = 1$,则A,B,C三点共线。

实际判断其实是用: $\vec{AB} \times \vec{AC} = 0 $,则三点共线

重心判断式

在$\Delta ABC$中,若$\vec{GA} + \vec{GB} + \vec{GC} = 0$,则G为$\Delta ABC$的重心。

重心: 三边中线的交点。

垂心判断式

在$\Delta ABC$中,若$\vec{HA} \cdot \vec{HB} = \vec{HB} \cdot \vec{HC} = \vec{HC} \cdot \vec{HA}$,则H为$\Delta ABC$的垂心。

将三角形分成了面积相等的三块。

外心判断式

在$\Delta ABC$中,若$\vec{OA} = \vec{OB} = \vec{OC} $,则O为$\Delta ABC$的外心。

外心: 垂直平分线的交点。

模板详解

判断点是否在线段(直线)上

其实就是三点共线的问题,解决办法是判断$\vec{AC} \times \vec{AB} = 0$.

但因为是线段,所以要判断横纵坐标的相应位置。

判断两线段是否相交

这里有两个实验:第一个是快速排斥实验,通过判断相应矩形是否相交来快速得到两线段是否相交。

但可以发现,这不能包含所有情况,因为两个矩形相交不一定能推出两条线段相交。

因此就有了跨立实验。

跨立,顾名思义,就是一条线所在的直线能够切开另一条线段。

如上面两个图,Q1,Q2所在直线可以跨立P1,P2线段; 但P1,P2所在直线无法跨立Q1,Q2线段, 因此两条线段不相交。

所以判断两条线段相交需要同时满足两条线段互相跨立的条件。

判断点在多边形中

射线法:

设要判断的点为P.

以P点为起点,水平向左做一条射线,根据射线与多边形的交点个数来判断P是否在多边形内部。

如果交点个数是奇数的话,说明P点在多边形内部。

但是!

会有一些特殊情况:比如图(b),若射线与多边形的交点是顶点的话,需要判断。

总的来说,步骤如下:

遍历多边形的每一条边:

  1. 如果P在多边形边上,可直接判断。
  2. 如果射线与多边形的边重合,那么可以直接忽略。
  3. 如果射线和多边形的一边有交点,那么交点数加一。(不包括顶点)
  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
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174

#define zero(x) ((fabs(x) < eps ? 1 : 0))
#define sgn(x) (fabs(x) < eps ? 0 : ((x) < 0 ? -1 : 1))

struct point
{
double x, y;
point(double a = 0, double b = 0) { x = a, y = b; }
point operator-(const point& b) const { return point(x - b.x, y - b.y); }
point operator+(const point& b) const { return point(x + b.x, y + b.y); }
// 两点是否重合
bool operator==(point& b) { return zero(x - b.x) && zero(y - b.y); }
// 点积(以原点为基准)
double operator*(const point& b) const { return x * b.x + y * b.y; }
// 叉积(以原点为基准)
double operator^(const point& b) const { return x * b.y - y * b.x; }
// 绕P点逆时针旋转a弧度后的点
point rotate(point b, double a)
{
double dx, dy;
(*this - b).split(dx, dy);
double tx = dx * cos(a) - dy * sin(a);
double ty = dx * sin(a) + dy * cos(a);
return point(tx, ty) + b;
}
// 点坐标分别赋值到a和b
void split(double& a, double& b) { a = x, b = y; }
};
struct line
{
point s, e;
line() {}
line(point ss, point ee) { s = ss, e = ee; }
};


double dist(point a, point b) { return sqrt((a - b) * (a - b)); }


// <0, *> 表示重合; <1, *> 表示平行; <2, P> 表示交点是P;
pair<int, point> spoint(line l1, line l2)
{
point res = l1.s;
if (sgn((l1.s - l1.e) ^ (l2.s - l2.e)) == 0)
return make_pair(sgn((l1.s - l2.e) ^ (l2.s - l2.e)) != 0, res);
double t = ((l1.s - l2.s) ^ (l2.s - l2.e)) / ((l1.s - l1.e) ^ (l2.s - l2.e));
res.x += (l1.e.x - l1.s.x) * t;
res.y += (l1.e.y - l1.s.y) * t;
return make_pair(2, res);
}

// 快速排斥实验 + 跨立实验
bool segxseg(line l1, line l2)
{
return
max(l1.s.x, l1.e.x) >= min(l2.s.x, l2.e.x) &&
max(l2.s.x, l2.e.x) >= min(l1.s.x, l1.e.x) &&
max(l1.s.y, l1.e.y) >= min(l2.s.y, l2.e.y) &&
max(l2.s.y, l2.e.y) >= min(l1.s.y, l1.e.y) &&
sgn((l2.s - l1.e) ^ (l1.s - l1.e)) * sgn((l2.e-l1.e) ^ (l1.s - l1.e)) <= 0 &&
sgn((l1.s - l2.e) ^ (l2.s - l2.e)) * sgn((l1.e-l2.e) ^ (l2.s - l2.e)) <= 0;
}

//l1是直线,l2是线段
bool segxline(line l1, line l2)
{
return sgn((l2.s - l1.e) ^ (l1.s - l1.e)) * sgn((l2.e - l1.e) ^ (l1.s - l1.e)) <= 0;
}

point pointtoline(point P, line L)
{
point res;
double t = ((P - L.s) * (L.e - L.s)) / ((L.e - L.s) * (L.e - L.s));
res.x = L.s.x + (L.e.x - L.s.x) * t, res.y = L.s.y + (L.e.y - L.s.y) * t;
return dist(P, res);
}

point pointtosegment(point p, line l)
{
point res;
double t = ((p - l.s) * (l.e - l.s)) / ((l.e - l.s) * (l.e - l.s));
if (t >= 0 && t <= 1)
res.x = l.s.x + (l.e.x - l.s.x) * t, res.y = l.s.y + (l.e.y - l.s.y) * t;
else
res = dist(p, l.s) < dist(p, l.e) ? l.s : l.e;
return res;
}

bool PointOnSeg(point p, line l)
{
return
sgn((l.s - p) ^ (l.e-p)) == 0 &&
sgn((p.x - l.s.x) * (p.x - l.e.x)) <= 0 &&
sgn((p.y - l.s.y) * (p.y - l.e.y)) <= 0;
}

double area(point p[], int n)
{
double res = 0;
for (int i = 0; i < n; i++) res += (p[i] ^ p[(i + 1) % n]) / 2;
return fabs(res);
}

// 点形成一个凸包, 而且按逆时针排序(如果是顺时针把里面的<0改为>0)
// 点的编号 : [0,n)
// -1 : 点在凸多边形外
// 0 : 点在凸多边形边界上
// 1 : 点在凸多边形内
int PointInConvex(point a, point p[], int n)
{
for (int i = 0; i < n; i++)
if (sgn((p[i] - a) ^ (p[(i + 1) % n] - a)) < 0)
return -1;
else if (PointOnSeg(a, line(p[i], p[(i + 1) % n])))
return 0;
return 1;
}

// 射线法,poly[]的顶点数要大于等于3,点的编号0~n-1
// -1 : 点在多边形外
// 0 : 点在多边形边界上
// 1 : 点在多边形内
int PointInPoly(point p, point poly[], int n)
{
int cnt;
line ray, side;
cnt = 0;
ray.s = p;
ray.e.y = p.y;
ray.e.x = -100000000000.0; // -INF,注意取值防止越界
for (int i = 0; i < n; i++)
{
side.s = poly[i], side.e = poly[(i + 1) % n];
if (PointOnSeg(p, side)) return 0;
//如果平行轴则不考虑
if (sgn(side.s.y - side.e.y) == 0)
continue;
if (PointOnSeg(sid e.s, ray))
cnt += (sgn(side.s.y - side.e.y) > 0);
else if (PointOnSeg(side.e, ray))
cnt += (sgn(side.e.y - side.s.y) > 0);
else if (segxseg(ray, side))
cnt++;
}
return cnt % 2 == 1 ? 1 : -1;
}
// 线段上整数点个数
int OnSegment(line l) { return __gcd(fabs(l.s.x - l.e.x), fabs(l.s.y - l.e.y)) + 1; }

//多边形边界上整数点个数
int OnEdge(point p[], int n)
{
int i, ret = 0;
for (i = 0; i < n; i++)
ret += __gcd(fabs(p[i].x - p[(i + 1) % n].x), fabs(p[i].y - p[(i + 1) % n].y));
return ret;
}

//多边形内部整数点个数
int InSide(point p[], int n)
{
int i, area = 0;
for (i = 0; i < n; i++)
area += p[(i + 1) % n].y * (p[i].x - p[(i + 2) % n].x);
return (fabs(area) - OnEdge(n, p)) / 2 + 1;
}

point waixin(point a, point b, point c)
{
double a1 = b.x - a.x, b1 = b.y - a.y, c1 = (a1 * a1 + b1 * b1) / 2;
double a2 = c.x - a.x, b2 = c.y - a.y, c2 = (a2 * a2 + b2 * b2) / 2;
double d = a1 * b2 - a2 * b1;
return point(a.x + (c1 * b2 - c2 * b1) / d, a.y + (a1 * c2 - a2 * c1) / d);
}
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
// 刘老师的板子
const double eps = 1e-10;
struct Point {
double x, y;
Point(double x = 0, double y = 0) : x(x), y(y) {}
};
typedef Point Vec;

Vec operator+(Vec A, Vec B) { return Vec(A.x + B.x, A.y + B.y); }
Vec operator-(Vec A, Vec B) { return Vec(A.x - B.x, A.y - B.y); }
Vec operator*(Vec A, double p) { return Vec(A.x * p, A.y * p); }
Vec operator/(Vec A, double p) { return Vec(A.x / p, A.y / p); }
bool operator<(const Point& a, const Point& b) {
return a.x < b.x || (a.x == b.x && a.y < b.y);
}
int dcmp(double x) {
if (fabs(x) < eps)
return 0;
else
return x < 0 ? -1 : 1;
}
bool operator==(const Point& a, const Point& b) {
return dcmp(a.x - b.x) == 0 && dcmp(a.y - b.y) == 0;
}
double PolarAngle(Point A) { return atan2(A.y, A.x); } // 极角

double Dot(Vec A, Vec B) { return A.x * B.x + A.y * B.y; } // 点积
double Length(Vec A) { return sqrt(Dot(A, A)); } // 向量长度
double Angle(Vec A, Vec B) { // 两向量夹角
return acos(Dot(A, B) / Length(A) / Length(B));
}
double Cross(Vec A, Vec B) { return A.x * B.y - A.y * B.x; }
double Area2(Point A, Point B, Point C) { // 平行四边形面积
return Cross(B - A, C - A);
}

Vec Rotate(Vec A, double rad) { // 向量绕起点逆时针旋转rad
return Vec(A.x * cos(rad) - A.y * sin(rad),
A.x * sin(rad) + A.y * cos(rad));
}

// 向量的单位法线,即左转90度后把长度归一化(要确保A不是零向量)
Vec Normal(Vec A) {
double L = Length(A);
return Vec(-A.y / L, A.x / L);
}
// 两直线交点,调用前要确保直线P + tv和Q + tw有唯一交点,当且仅当Cross(v,w)!=0
Point GetLineItersection(Point P, Vec v, Point Q, Vec w) {
Vec u = P - Q;
double t = Cross(w, u) / Cross(v, w);
return P + v * t;
}

// P 到直线AB的距离
double DistanceToLine(Point P, Point A, Point B) {
Vec v1 = B - A, v2 = P - A;
return fabs(Cross(v1, v2) /
Length(v1)); // 如果不取绝对值,得到的是有向距离
}

double DistanceToSegment(Point P, Point A, Point B) {
if (A == B) return Length(P - A);
Vec v1 = B - A, v2 = P - A, v3 = P - B;
if (dcmp(Dot(v1, v2)) < 0)
return Length(v2);
else if (dcmp(Dot(v1, v3)) > 0)
return Length(v3);
else
return fabs(Cross(v1, v2)) / Length(v1);
}

// 求点与线段距离的交点
Point GetLineProjection(Point P, Point A, Point B) {
Vec v = B - A;
return A + v * (Dot(v, P - A) / Dot(v, v));
}

// 规范相交: 两线段只有一个交点,且不在任何一条线段的断点
bool SegmentProperIntersection(Point a1, Point a2, Point b1, Point b2) {
double c1 = Cross(a2 - a1, b1 - a1), c2 = Cross(a2 - a1, b2 - a1),
c3 = Cross(b2 - b1, a1 - b1), c4 = Cross(b2 - b1, a2 - b1);
return dcmp(c1) * dcmp(c2) < 0 && dcmp(c3) * dcmp(c4) < 0;
}
// 点是否在线段上(不包括端点)
bool OnSegment(Point p, Point a1, Point a2) {
return dcmp(Cross(a1 - p, a2 - p)) == 0 && dcmp(Dot(a1 - p, a2 - p)) < 0;
}

// 多边形的有向面积
double PolygonArea(Point* p, int n) {
double area = 0;
for (int i = 1; i < n - 1; i++) area += Cross(p[i] - p[0], p[i + 1] - p[0]);
return area / 2;
}
Thank you for your support!