计算几何

叉积

就是判断两个点 \(P_1(x_1,y_1),P_2(x_2,y_2)\) 极角的大小关系。值为 \(x_1y_2-x_2y_1\)。

该值为正，意味着 \((x_1,y_1)\) 极角更大，换言之极坐标系上线段 \(OP_1\) 在 \(OP_2\) 的逆时针 / 正方向。否则类似。

如果钦定一个点 \(P_0\) 为原点，就可以计算以它为一端的两个线段的相对方向，即根据 \((x_1-x_0)(y_2-y_0)-(x_2-x_0)(y_1-y_0)\) 的正负性。

值得一提的是，叉积在数值上则等于两向量围成的平行四边形的有向面积。除以二就是这三个点为顶点的三角形面积。

判断线段相交

使用除法暴力计算交点的可扩展性和精度看起来都不太好。

线段 \(A,B\) 相交当且仅当以下两者之一：

\(A_1,A_2\) 分属直线 \(b\) 的两侧，且\(B_1,B_2\) 分属直线 \(a\) 的两侧。
\(A_1,A_2\) 中有一个点位于 \(b\) 上，或 \(B_1,B_2\) 中有一个点位于 \(a\) 上。

前者相当于找一个线段的一端，然后用叉积判断线段 \(B_1A_1,B_1A_2\) 是不是分属线段 \(B_1B_2\) 的两侧。后者当且仅当判断时出现叉积相等（不在任何一侧）。

求凸包

选择喜闻乐见单调栈。

钦定横坐标最小的点，这个点一定在凸包上。然后其它点按照以其为原点的极角排序，这时候凸包的点一定按照排序后的顺序逆时针出现。

注意：对于极角相等，按照到原点距离排序。否则可能这些点就不满足逆时针出现的性质。

const double eps = 1e-8;
const int maxn = 200005;
struct Point {
	double x, y;
	const double operator* (Point t) {
		return x * t.y - t.x * y;
	}
	const bool operator< (Point t) const {
		return y == t.y ? x < t.x : y < t.y;
	}
	Point operator- (Point t) {
		return {x - t.x, y - t.y};
	}
	Point operator+ (Point t) {
		return {x + t.x, y + t.y};
	}
} a[maxn];
inline double dis(Point a, Point b) {
	Point c = a - b;
	return c.x * c.x + c.y * c.y;
}
// 我们充分发扬类人智慧
Point o;
inline bool cmp(Point a, Point b) {
	Point u = a - o, v = b - o;
	if (abs(u * v) < eps) return dis(o, a) < dis(o, b);
	else return (u * v) > 0;
}
inline bool check(Point a, Point b, Point c) {
	Point u = b - a, v = c - a;
	return (u * v) > 0;
}
void get(Point p[], int n, int st[], int &top) {
	sort(p + 1, p + n + 1), o = p[1], sort(p + 2, p + n + 1, cmp), st[++top] = 1;
	for (int i = 2; i <= n; ++i) {
		while (top > 1 && !check(p[st[top - 1]], p[st[top]], p[i])) --top;
		st[++top] = i;
	}
}

upd：略微区别的算法：考虑先求上凸壳再拼上下凸壳。直接按照 \(x,y\) 分别为第一，第二关键字排序，正反分别跑一遍即可。

void get(Point p[], int n, int st[], int &top) {
	sort(p + 1, p + 1 + n), st[++top] = 1;
	for (int i = 2; i <= n; i++) {
		while (top > 1 && (p[st[top]] - p[st[top - 1]]) * (p[i] - p[st[top]]) <= 0) top--;
		st[++top] = i;
	}
	int tmp = top;
	for (int i = n - 1; i; i--) {
		// 注意不能把之前的上凸壳除掉，故需要 top >= tmp。
		while (top > tmp && (p[st[top]] - p[st[top - 1]]) * (p[i] - p[st[top]]) <= 0) top--;
		st[++top] = i;
	}
}

旋转卡壳

对于凸包上的一个点引一条不相交于凸包的直线，然后在凸包另一侧作另一条切该凸包的该直线的平行线。

注意到依次对于各点作直线，其对应的切点（可能为切边）也是单调的。据此可以 \(O(n)\) 算出这些对应的点 / 边。

那么注意到这些对应点 / 边具有极大性，显然可以用来求凸包直径 / 最小矩形覆盖等。

注意：依次扫描各个点，得到每个点的最远点，这并不具有单调性。

具体实现时，我们发现一条边的对应点显然是距离其最大的点，这个是有单调性。这等价于和它围成的三角形面积最大。可以拿叉积利用单调性扫描简化实现。

// 旋转卡壳求直径
// 计算面积统一不除以二
inline double calc(int a, int b, int c) {
	return abs((p[st[a]] - p[st[c]]) ^ (p[st[b]] - p[st[c]]));
}
void rotate_qia_ke(int st[], int &top){
	st[++top] = st[1];
	for (int i = 1, j = 1; i < top; ++i) {
		while (calc(i, i + 1, j) < calc(i, i + 1, j % top + 1))
			j = j % top + 1, ans = max(ans, max(dis(p[st[i]], p[st[j]]), dis(p[st[j]], p[st[i + 1]])));
		ans = max(ans, max(dis(p[st[i]], p[st[j]]), dis(p[st[j]], p[st[i + 1]])));
	}
}

闵可夫斯基和

定义为 \(C=\{a+b\mid a\in A,b\in B\}\) 其中 \(A,B,C\) 为点集。

注意到满足交换律、结合律、可减性。

关于求解，可以观察发现 \(C\) 有 \(|A|\cdot|B|\) 个点，但是其凸包至多由 \(|A|+|B|\) 条来自 \(A,B\) 的凸包的边构成。

人类智慧地，直接对所有 \(A,B\) 中的边归并成极角有序的，然后对相邻两条边，不断平移到上一条的端点，就得到 \(C\) 对应的多边形。

但是显然对于一些平行的 corner case 这样得到的不一定是凸包，再跑一遍求凸包即可。

void solve() {
	get(a, n, sa, ta), get(b, m, sb, tb);
	for (int i = 1; i <= ta; ++i) d1[i] = a[sa[i % ta + 1]] - a[sa[i]];
	for (int i = 1; i <= tb; ++i) d2[i] = b[sb[i % tb + 1]] - b[sb[i]];
	c[1] = a[sa[1]] + b[sb[1]];
	int tot = 1, p1 = 1, p2 = 1;
	while (p1 <= ta && p2 <= tb) ++tot, c[tot] = c[tot - 1] + (d1[p1] * d2[p2] >= 0 ? d1[p1++] : d2[p2++]);
	while (p1 <= ta) ++tot, c[tot] = c[tot - 1] + d1[p1++];
	while (p2 <= tb) ++tot, c[tot] = c[tot - 1] + d2[p2++];
	get(c, ta + tb, sc, tc);
}

求最近点对

最远点对直接旋转卡壳就好了。

关于此，我们充分发扬人类智慧：

如果当前答案为 \(d\)，且不扫描与当前点切比雪夫距离超过 \(d\) 的点（这些答案显然更劣），用这些点到当前点的距离更新答案，那么均摊最多扫描常数个点。

考虑构造一个以当前点为中心，边长为 \(2d\) 的正方形，则可以证明。

但是这个要利用 set，而利用上述思想进行分治可以得到更优的理论复杂度。总共都是大常数 \(O(n\log n)\)。

这个也可以支持一些修改和维护？

2-D Tree 最近点对是假的。

set<Point> st;
bool cmp(Point a, Point b) { return a.y < b.y; }
inline ll sqr(ll x) { return x * x; }
inline ll dis(Point a, Point b) {
	return sqr(a.x - b.x) + sqr(a.y - b.y);
}
void nearest(){
	for (int j = 1, i = 1; i <= n; ++i) {
		while (j < i && sqr(p[i].y - p[j].y) > ans) st.erase(p[j]), ++j;
		auto tt = st.lower_bound(p[i]), it = tt;
		if (it != st.begin()) {
			--it;
			while (sqr(p[i].x - (it -> x)) < ans) {
				ans = min(ans, dis(p[i], *it));
				if (it != st.begin()) --it;
				else break;
			}
		}
		it = tt;
		while (it != st.end() && sqr(p[i].x - (it -> x)) < ans)
			ans = min(ans, dis(p[i], *it)), ++it;
		st.insert(p[i]);
	}
}

最小圆覆盖

直接以凸包直径为直径这个想法，随便构造一个筝形就 hack 了。

考虑一个暴力：凸包上三个点确定一个圆，因此朴素枚举之。看起来是 \(O(n^3)\)。然而我们可以证明：枚举前，先随机打乱点集顺序，并且进行剪枝，就可以期望 \(O(n)\)。

具体来说，考虑 \(C_i\) 为前 \(i\) 个点的最小覆盖圆。考虑从 \(C_{i-1}\) 出发递推 \(C_i\)，显然发生改变当且仅当 \(i\) 在 \(C_{i-1}\) 圆外。这时候显然 \(i\) 为 \(C_i\) 上的一个点，继续枚举另外两个点 \(j,k\in[1,i)\) 以确定 \(C_i\)，否则不进入内层循环，直接令 \(C_{i}=C_{i-1}\)

对于 \(j\) 的枚举同理，\(C_{i,j}\) 表示 \(i\cup [1,j]\) 的最小覆盖圆，同样递推求 \(C_{i,i-1}\)。仅当 \(j\) 在 \(C_{i,j-1}\) 圆外时需要继续枚举 \(k\)，否则令 \(C_{i,j}=C_{i,j-1}\)。

对于 \(k\) 的枚举，同样维护 \({i,j}\cup [1,k]\) 的最小覆盖圆，若 \(k\) 在圆外就用 \(i,j,k\) 确定的圆更新之。

考虑证明复杂度，容易发现每次进入内层循环的概率是 \(\dfrac{3}{i}\)。故复杂度为

\[T(n)=O(n)+\sum_{i=1}^n [O(n)+\dfrac{3}{i}\sum_{j=1}^{i-1} (O(n)+\dfrac{3}{j}\sum_{k=1}^{j-1} O(k))]=O(n) \]

杂项

判断点是否在多边形内：从该点向右引一射线，与多边形交奇数次则在其内。

判断点是否在凸包内：稍微简单，直接随便找个原点，在凸包上二分出其极角的前驱和后继点，然后类比求凸包时候拿凸包判定。注意到这个复杂度是对数的。

bool checkin(Point c[], int tc, Point p) {
	if((p * c[1]) > 0 || (c[tc] * p) > 0) return 0;
	int i = lower_bound(c + 1, c + tc + 1, p, cmp) - c - 1;
	return (p - c[i]) * (c[i % tc + 1] - c[i]) <= 0;
}

求凸包面积：直接随便取一个固定顶点加上每个边拆成一堆三角形叉积求和再除以二。

从这里开始是结合数学问题的应用。

半平面交

半平面：一条直线的一侧部分的平面。

其实是一类二维线性规划的求解方法。就是先维护当前的交多边形，然后依次加入边，若相交则把交点加入再弹出不满足的点。注意到这样始终满足交多边形的凸性。

然后对于最小化 / 最大化 \(ax_1+bx_2\) 的目标，在凸包上拿斜率为 \(\dfrac{a}{b}\) 的直线去最大化 / 最小化截距即可。

Quick Hull

是一种 \(O(n^2)\) 的“快速”求凸包方法。但是和 MST 领域同样相对来说难写又难调的 Brouvka 一样，适用于点（边）规模巨大，且有特殊性质的情况。

同样维护一个凸包上的点集，插入横坐标最小的点 \(A\) 与纵坐标最小的点 \(B\)（一定在凸包上）。然后每次寻找一个尚未被加入的点 \(C\) 使得其离线段 \(AB\) 最远。相当于使得 \(|\overrightarrow{AC}\times \overrightarrow{BC}|\) 最大，容易证明 \(C\) 也一定在凸包上。将其加入，然后递归做 \(\overrightarrow{AC},\overrightarrow{CB}\)。

易知复杂度为 \(O(n^2)\)。但是在 timeismoney 或者画框这种题目上就有奇效。

具体来说，可以证明如果把权值和为 \(\sum x,\sum y\) 的方案（对应于被选择元素的集合 \(S\)）看成点 \((\sum x,\sum y)\)，那么最优的方案一定在凸包上，并且可知凸包的点数其实是 \(O(V^2)\) 的，\(V\) 为 \(x_i,y_i\) 的值域。

只需考虑每次计算，可以发现 \(C(\sum x,\sum y)\) 对应的叉积是形如 \(a\sum x+b\sum y+c\) 的式子，然后就可以把原来有两个权值 \(x_i,y_i\) 的点赋予权值 \(w(i)=ax_i+by_i\) 然后直接跑对应的最小化算法（MST / 网络流 / ...）计算。

二维斜率优化

upd：几乎完全打不过李超线段树。

考虑斜率优化的一般形式：\(y=kx+b\)，其中 \(y\) 为仅关于 \(j\) 的项，\(b\) 为仅关于 \(i\) 的项，\(k,x\) 分别为乘积项中关于 \(i,j\) 的部分。

但是这个东西在式子中存在形如 \(f_1(i)g_1(j)+f_2(i)g_2(j)\) 的两个不同乘积项的时候就不适用了。但是其实舍弃仅关于 \(j\) 的项后，上述仍然可以表示为 \((f_1(i),-f_2(i)),(g_2(j),g_1(j))\) 的叉积。而最大化叉积等价于最大化三角形 有向面积，也就是向量 \((f_1(i),-f_2(i))\) 与点 \((g_2(j),g_1(j))\) 的 有向距离。

据此我们仍然直接拿平行于前者的直线去切凸包，找到最远点即可。若整个凸包位于向量同一侧，那么可能找的就是最近点。

可能需要平衡树 / CDQ 来处理询问。

动态凸包

直线型李超树实际上就是动态半平面交弱化版。所以动态凸包是其超集（但不能处理少见的线段型）但复杂得多。

考虑和单调栈求凸包一样，维护一个按极角排序的凸包点集。动态插入时删除不必要的点即可。

Depot

Luogu P6719 [BalkanOI2011] 2circles

给定一个凸包及实数 \(r\)，求能否在平面上放两个半径为 \(r\) 的圆，使其互不相交且完全在凸包内（均允许边界相切）。

Trick：考虑凸包上所有边向凸包所在侧（垂直于边方向）平移 \(r\)，即将凸包内缩 \(r\) 个单位，然后旋转卡壳判断凸包直径是否大于等于 \(2r\)。

Luogu P4293 [WC2010] 能量场

给你 \(n\) 个粒子，每个粒子有两个权值 \(m_i,c_i\)，粒子构成的每个相邻有序对 \((a,b)\) 会产生 \(m_am_b(c_a-c_b)\) 的贡献。要求：1. 找出一个有序对使贡献最大；2. 找出一个序列成环后贡献和最大。

先将贡献转化成叉积：

\[m_am_b(c_a-c_b)=m_ac_am_b-m_bc_bm_a \]

即将每个元素看作 \(x=m_ic_i,y=m_i\) 的向量。

那么第一问就等价于求叉积最大的两个向量，这等价于围成的三角形面积最大。因为所有点都在第一象限，显然这样的向量终点一定是凸包上的点。直接跑旋转卡壳；但是在边界上可能统计不到，考虑正反跑两遍。

第二问可以类似地推出，就是凸包面积的两倍。

Luogu P8101 [USACO22JAN] Multiple Choice Test P

Key Observation：所有最后可能组合出的（指数级别个）点构成的点集中，只有该点集凸包上的（多项式级别个）点可能产生贡献。证明显然。

而求出这个凸包的方法就是套路了。先对每组内点求凸包，将所有凸包做闵可夫斯基和。注意：多凸包做闵可夫斯基和其实和两个凸包并无差别，一起差分排序再归并即可。相当于多数组归并。

posted @ 2023-08-27 21:31 音街ウナ阅读(19) 评论(0) 收藏举报

刷新页面返回顶部

148.