计算几何模板与技巧：基础构造与关系处理
本文介绍了一套用于计算几何的高效模板，包括点、向量、直线、线段和圆的基本结构，以及如何处理浮点数比较、基本操作、几何对象关系、平移旋转、对称和交点计算等。通过实例和依赖关系说明，读者可以快速上手解决几何问题。
引用：
https://zhuanlan.zhihu.com/p/338057154
要做计算几何题，首先要有一套基础的模板，以方便处理几何对象。在多次实践后，目前的模板采取以下原则编写：
写全局函数，而非类方法，结构体只存数据；
每个函数标注依赖于哪些函数，且尽量减少依赖；
用较简略的名字，同时传值而非传const引用。
几何对象
// 点
struct Point { double x, y; };
// 向量
using Vec = Point;
// 直线（点向式）
struct Line { Point P; Vec v; };
// 线段（存两个端点）
struct Seg { Point A, B; };
// 圆（存圆心和半径）
struct Circle { Point O; double r; };
当然，向量和点并不是完全一样的，但很多时候这两者可以混用，所以为了方便就这样处理了（这可以使很多函数变得简洁），特殊的情况我后面会指出。
常用常量
// 原点
const Point O = {0, 0};
// 坐标轴
const Line Ox = {O, {1, 0}}, Oy = {O, {0, 1}};
const double PI = acos(-1), EPS = 1e-9;
浮点比较
这个其实不是计算几何内容，但计算几何是ACM中最常涉及浮点数的部分。因为浮点数运算会导致浮点误差，使得比较运算符出问题，所以我们要写一套能够“容忍”一定浮点误差的浮点数比较函数。
// ==
bool eq(double a, double b) { return abs(a - b) < EPS; }
// >
bool gt(double a, double b) { return a - b > EPS; }
// <
bool lt(double a, double b) { return a - b < -EPS; }
// >=
bool ge(double a, double b) { return a - b > -EPS; }
// <=
bool le(double a, double b) { return a - b < EPS; }
基础操作
// 逆时针旋转90度的向量
Vec r90a(Vec v) { return {-v.y, v.x}; }
// 顺时针旋转90度的向量
Vec r90c(Vec v) { return {v.y, -v.x}; }
// 向量加向量
Vec operator+(Vec u, Vec v) { return {u.x + v.x, u.y + v.y}; }
// 向量减向量
Vec operator-(Vec u, Vec v) { return {u.x - v.x, u.y - v.y}; }
// 数乘
Vec operator*(double k, Vec v) { return {k * v.x, k * v.y}; }
// 点乘
double operator*(Vec u, Vec v) { return u.x * v.x + u.y * v.y; }
// 叉乘
double operator^(Vec u, Vec v) { return u.x * v.y - u.y * v.x; }
// 向量长度
double len(Vec v) { return sqrt(v.x * v.x + v.y * v.y); }
// 斜率 // NOTE 不要用isinf判断斜率不存在，用后面的paral_y
double slope(Vec v) { return v.y / v.x; }
向量相关
// 两向量的夹角余弦
// DEPENDS len, V*V
double cos_t(Vec u, Vec v) { return u * v / len(u) / len(v); }
// 归一化向量（与原向量方向相同的单位向量）
// DEPENDS len
Vec norm(Vec v) { return {v.x / len(v), v.y / len(v)}; }
// 与原向量平行且横坐标大于等于0的单位向量
// DEPENDS d*V, len
Vec pnorm(Vec v) { return (v.x < 0 ? -1 : 1) / len(v) * v; }
// 线段的方向向量
// DEPENDS V-V
// NOTE 直线的方向向量直接访问属性v
Vec dvec(Seg l) { return l.B - l.A; }
直线相关
// 两点式直线
// DEPENDS V-V
Line line(Point A, Point B)
{ return {A, B - A}; }
// 斜截式直线
Line line(double k, double b)
{ return {{0, b}, {1, k}}; }
// 点斜式直线
Line line(Point P, double k)
{ return {P, {1, k}}; }
// 线段所在直线
// DEPENDS V-V
Line line(Seg l)
{ return {l.A, l.B - l.A}; }
// 给定直线的横坐标求纵坐标
// NOTE 请确保直线不与y轴平行
double at_x(Line l, double x)
{ return l.P.y + (x - l.P.x) * l.v.y / l.v.x; }
// 给定直线的纵坐标求横坐标
// NOTE 请确保直线不与x轴平行
double at_y(Line l, double y)
{ return l.P.x - (y + l.P.y) * l.v.x / l.v.y; }
// 点到直线的垂足
// DEPENDS V-V, V*V, d*V
Point pedal(Point P, Line l)
{ return l.P - (l.P - P) * l.v / (l.v * l.v) * l.v; }
// 过某点作直线的垂线
// DEPENDS r90c
Line perp(Line l, Point P)
{ return {P, r90c(l.v)}; }
// 角平分线
// DEPENDS V+V, len, norm
Line bisec(Point P, Vec u, Vec v)
{ return {P, norm(u) + norm(v)}; }
线段相关
// 线段的方向向量
// DEPENDS V-V
// NOTE 直线的方向向量直接访问属性v
Vec dvec(Seg l)
{ return l.B - l.A; }
// 线段中点
Point midp(Seg l)
{ return {(l.A.x + l.B.x) / 2, (l.A.y + l.B.y) / 2}; }
// 线段中垂线
// DEPENDS r90c, V-V, midp
Line perp(Seg l)
{ return {midp(l), r90c(l.B - l.A)}; }
几何对象之间的关系
// 向量是否互相垂直
// DEPENDS eq, V*V
bool verti(Vec u, Vec v) { return eq(u * v, 0); }
// 向量是否互相平行
// DEPENDS eq, V^V
bool paral(Vec u, Vec v) { return eq(u ^ v, 0); }
// 向量是否与x轴平行
// DEPENDS eq
bool paral_x(Vec v) { return eq(v.y, 0); }
// 向量是否与y轴平行
// DEPENDS eq
bool paral_y(Vec v) { return eq(v.x, 0); }
// 点是否在直线上
// DEPENDS eq
bool on(Point P, Line l)
{ return eq((P.x - l.P.x) * l.v.y, (P.y - l.P.y) * l.v.x); }
// 点是否在线段上
// DEPENDS eq, len, V-V
bool on(Point P, Seg l)
{ return eq(len(P - l.A) + len(P - l.B), len(l.A - l.B)); }
// 两个点是否重合
// DEPENDS eq
bool operator==(Point A, Point B)
{ return eq(A.x, B.x) && eq(A.y, B.y); }
// 两条直线是否重合
// DEPENDS eq, on(L)
bool operator==(Line a, Line b)
{ return on(a.P, b) && on(a.P + a.v, b); }
// 两条线段是否重合
// DEPENDS eq, P==P
bool operator==(Seg a, Seg b)
{ return (a.A == b.A && a.B == b.B)
|| (a.A == b.B && a.B == b.A); }
// 以横坐标为第一关键词、纵坐标为第二关键词比较两个点
// DEPENDS eq, lt
bool operator<(Point A, Point B)
{ return lt(A.x, B.x)
|| (eq(A.x, B.x) && lt(A.y, B.y)); }
// 直线与圆是否相切
// DEPENDS eq, V^V, len
bool tangency(Line l, Circle C)
{ return eq(abs((C.O ^ l.v) - (l.P ^ l.v)), C.r * len(l.v)); }
// 圆与圆是否相切
// DEPENDS eq, V-V, len
bool tangency(Circle C1, Circle C2)
{ return eq(len(C1.O - C2.O), C1.r + C2.r); }
距离
// 两点间的距离
// DEPENDS len, V-V
double dis(Point A, Point B)
{ return len(A - B); }
// 点到直线的距离
// DEPENDS V^V, len
double dis(Point P, Line l)
{ return abs((P ^ l.v) - (l.P ^ l.v)) / len(l.v); }
// 平行直线间的距离
// DEPENDS d*V, V^V, len, pnorm
// NOTE 请确保两直线是平行的
double dis(Line a, Line b)
{ return abs((a.P ^ pnorm(a.v)) - (b.P ^ pnorm(b.v))); }
平移和旋转
// 平移
// DEPENDS V+V
Line operator+(Line l, Vec v) { return {l.P + v, l.v}; }
Seg operator+(Seg l, Vec v) { return {l.A + v, l.B + v}; }
// 旋转
// DEPENDS V+V, V-V
Point rotate(Point P, double rad)
{ return {cos(rad) * P.x - sin(rad) * P.y,
sin(rad) * P.x + cos(rad) * P.y}; }
// DEPENDS ^1
Point rotate(Point P, double rad, Point C)
{ return C + rotate(P - C, rad); }
// DEPENDS ^1, ^2
Line rotate(Line l, double rad, Point C = O)
{ return {rotate(l.P, rad, C), rotate(l.v, rad)}; }
// DEPENDS ^1, ^2
Seg rotate(Seg l, double rad, Point C = O)
{ return {rotate(l.A, rad, C), rotate(l.B, rad, C)}; }
对称
这里是本文中唯一一处点和向量表现不同的地方。
// 对称
// 关于点对称
Point reflect(Point A, Point P)
{ return {P.x * 2 - A.x, P.y * 2 - A.y}; }
// DEPENDS ^1
Line reflect(Line l, Point P)
{ return {reflect(l.P, P), l.v}; }
// DEPENDS ^1
Seg reflect(Seg l, Point P)
{ return {reflect(l.A, P), reflect(l.B, P)}; }
// 关于直线对称
// DEPENDS V-V, V*V, d*V, pedal
// NOTE 向量和点在这里的表现不同，
// 求向量关于某直线的对称向量需要用reflect_v
// DEPENDS ^1
Point reflect(Point A, Line ax)
{ return reflect(A, pedal(A, ax)); }
// DEPENDS ^1, ^4
Vec reflect_v(Vec v, Line ax)
{ return reflect(v, ax) - reflect(O, ax); }
// DEPENDS ^1, ^4, ^5
Line reflect(Line l, Line ax)
{ return {reflect(l.P, ax), reflect_v(l.v, ax)}; }
// DEPENDS ^1, ^4
Seg reflect(Seg l, Line ax)
{ return {reflect(l.A, ax), reflect(l.B, ax)}; }
交点
// 直线与直线交点
// DEPENDS eq, d*V, V*V, V+V, V^V
vector<Point> inter(Line a, Line b)
{
double c = a.v ^ b.v;
if (eq(c, 0)) return {};
Vec v = 1 / c * Vec{a.P ^ (a.P + a.v), b.P ^ (b.P + b.v)};
return {{v * Vec{-b.v.x, a.v.x}, v * Vec{-b.v.y, a.v.y}}};
}
// 直线与圆交点
// DEPENDS eq, gt, V+V, V-V, V*V, d*V, len, pedal
vector<Point> inter(Line l, Circle C)
{
Point P = pedal(C.O, l);
double h = len(P - C.O);
if (gt(h, C.r)) return {};
if (eq(h, C.r)) return {P};
double d = sqrt(C.r * C.r - h * h);
Vec vec = d / len(l.v) * l.v;
return {P + vec, P - vec};
}
// 圆与圆的交点
// DEPENDS eq, gt, V+V, V-V, d*V, len, r90c
vector<Point> inter(Circle C1, Circle C2)
{
Vec v1 = C2.O - C1.O, v2 = r90c(v1);
double d = len(v1);
if (gt(d, C1.r + C2.r) || gt(abs(C1.r - C2.r), d))
return {};
if (eq(d, C1.r + C2.r) || eq(d, abs(C1.r - C2.r)))
return {C1.O + C1.r / d * v1};
double a = ((C1.r * C1.r - C2.r * C2.r) / d + d) / 2;
double h = sqrt(C1.r * C1.r - a * a);
Vec av = a / len(v1) * v1, hv = h / len(v2) * v2;
return {C1.O + av + hv, C1.O + av - hv};
}
三角形的四心
// 三角形的重心
Point barycenter(Point A, Point B, Point C)
{
return {(A.x + B.x + C.x) / 3, (A.y + B.y + C.y) / 3};
}
// 三角形的外心
// DEPENDS r90c, V*V, d*V, V-V, V+V
// NOTE 给定圆上三点求圆，要先判断是否三点共线
Point circumcenter(Point A, Point B, Point C)
{
double a = A * A, b = B * B, c = C * C;
double d = 2 * (A.x * (B.y - C.y)
+ B.x * (C.y - A.y) + C.x * (A.y - B.y));
return 1 / d * r90c(a * (B - C)
+ b * (C - A) + c * (A - B));
}
// 三角形的内心
// DEPENDS len, d*V, V-V, V+V
Point incenter(Point A, Point B, Point C)
{
double a = len(B - C), b = len(A - C), c = len(A - B);
double d = a + b + c;
return 1 / d * (a * A + b * B + c * C);
}
// 三角形的垂心
// DEPENDS V*V, d*V, V-V, V^V, r90c
Point orthocenter(Point A, Point B, Point C)
{
double n = B * (A - C), m = A * (B - C);
double d = (B - C) ^ (A - C);
return 1 / d * r90c(n * (C - B) - m * (C - A));
}
扫一扫
2万+
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
