Skip to content

计算几何

2024-05-21

人一定要有梦想。

精度处理

精度设置

精度的最小误差 eps 通常在 108 以内,视题目要求改变。

cpp
using var = double;
const var eps = 1e-8;

符号函数 sign

sign(x)={1,x<eps0,epsxeps1,x>eps
cpp
inline int sign(var x) {
	return x < -eps ? -1 : (x > eps ? 1 : 0);
}

绝对值函数 absolute

absolute(x)=xsign(x)
cpp
inline var absolute(var x) {
	return x * sign(x);
}

符号约定

Point

二维平面点由 x,y 坐标组成。

cpp
struct Point {
	var x, y;
	Point(var x = 0, var y = 0)
		: x(x), y(y) {}
};

在精度范围内判断相等。

cpp
bool operator == (const Point& a, const Point& b) {
	return ! sign(a.x - b.x) && ! sign(a.y - b.y);
}

标准输入输出支持。

cpp
istream& operator >> (istream& in, Point& p) {
	in >> p.x >> p.y;
	return in;
}

ostream& operator << (ostream& out, const Point& p) {
	out << "(" << p.x << ", " << p.y << ")";
	return out;
}

向量 Vector

向量 Vector 的定义沿用 Point 的定义。

cpp
using Vector = Point;

多边形 Poly

多边形 Poly 由一组 Point 构成。

cpp
using Poly = vector<Point>;

标准输入输出支持。

cpp
istream& operator >> (istream& in, Poly& poly) {
	for (auto& p : poly) {
		in >> p;
	}
	return in;
}

ostream& operator << (ostream& out, Poly& poly) {
	for (const auto& p : poly) {
		out << p << ",";
	}
	return out;
}

向量运算

加减 operator +/-

(x1,y1)±(x2,y2)=(x1±x2,y1±y2)
cpp
Vector operator + (const Vector& a, const Vector& b) {
	return { a.x + b.x, a.y + b.y };
}

Vector operator - (const Vector& a, const Vector& b) {
	return { a.x - b.x, a.y - b.y };
}

数乘 operator *

(x,y)k=(xk,yk)
cpp
Vector operator * (const Vector& p, var k) {
	return { p.x * k, p.y * k };
}

点乘 dot

(x1,y1)(x2,y2)=x1x2+y1y2
cpp
var dot(const Vector& a, const Vector& b) {
	return a.x * b.x + a.y * b.y;
}

叉乘 cross

(x1,y1)×(x2,y2)=x1y2x2y1
cpp
var cross(const Vector& a, const Vector& b) {
	return a.x * b.y - a.y * b.x;
}

取模 module

(x,y)=x2+y2p=pp
cpp
var module(const Vector& p) {
	return sqrt(dot(p, p));
}

旋转 rotate

向量 p 逆时针旋转 θ 弧度。

(x,y)(cosθsinθsinθcosθ)
cpp
Vector rotate(const Vector& p, var theta) {
	return {
		p.x * cos(theta) - p.y * sin(theta),
		p.x * sin(theta) + p.y * cos(theta)
	};
}

向量 vec

从点 A 指向点 B 的向量。

cpp
Vector vec(const Point& a, const Point& b) {
	return b - a;
}

点|点

距离 to_point

cpp
var to_point(const Point& a, const Point& b) {
	return module(a - b);
}

旋转 rotate

P 绕点 C 逆时针旋转 θ 弧度。

cpp
Point rotate(const Point& p, const Point& c, var theta) {
	return c + rotate(p - c, theta);
}

点|线段

相交 on_seg

P 是否在线段 AB 上。

cpp
bool on_seg(const Point& p, const Point& a, const Point& b) {
	return ! sign(cross(p - a, b - a)) and sign(dot(p - a, p - b)) <= 0;
}

距离 to_seg

P 到线段 AB 的距离。

cpp
var to_seg(const Point& p, const Point& a, const Point& b) {
	if (a == b)
		return to_point(p, a);
	Vector x = p - a, y = p - b, z = b - a;
	if (sign(dot(x, z)) < 0)
		return module(x);
	if (sign(dot(y, z)) > 0)
		return module(y);
	return absolute(cross(x, z) / module(z));
}

点|直线

相交 on_line

P 是否在直线 AB 上。

cpp
bool on_line(const Point& p, const Point& a, const Point& b) {
	return ! sign(cross(p - a, b - a));
}

距离 to_line

P 到直线 AB 的距离。

cpp
var to_line(const Point& p, const Point& a, const Point& b) {
	Vector x = p - a, z = b - a;
	return absolute(cross(x, z) / module(z));
}

垂足 foot_point

P 到直线 AB 的垂足。

cpp
Point foot_point(const Point& p, const Point& a, const Point& b) {
	Vector x = p - a, y = p - b, z = b - a;
	var prj1 = dot(x, z) / module(z);
	var prj2 = dot(y, z) / module(z);
	return a + z * (prj1 / (prj1 - prj2));
}

对称点 symmetry_point

P 关于直线 AB 的对称点。

cpp
Point symmetry_point(const Point& p, const Point& a, const Point& b) {
	return foot_point(p, a, b) * 2 - p;
}

线|线

直线|直线 交点 cross_point

直线 AB 和直线 CD 的交点。

cpp
Point cross_point(const Point& a, const Point& b, const Point& c, const Point& d) {
	Vector x = b - a, y = d - c, z = a - c;
	return a + x * (cross(y, z) / cross(x, y));
}

直线|线段 相交 inter_line_seg

直线 AB 和线段 CD 是否相交。

cpp
bool inter_line_seg(const Point& a, const Point& b, const Point& c, const Point& d) {
	return on_seg(cross_point(a, b, c, d), c, d);
}

线段|线段 相交 inter_seg_seg

线段 AB 和线段 CD 是否相交。

cpp
bool inter_seg_seg(const Point& a, const Point& b, const Point& c, const Point& d) {
	var c1 = cross(b - a, c - a), c2 = cross(b - a, d - a);
	var d1 = cross(d - c, a - c), d2 = cross(d - c, b - c);
	return sign(c1) * sign(c2) < 0 and sign(d1) * sign(d2) < 0;
}

多边形

有向面积 area

逆时针为正向。

cpp
var area(const Poly& poly) {
	var s = 0;
	int n = poly.size();
	for (int i = 0; i < n; i ++)
		s += cross(poly[i], poly[(i + 1) % n]);
	return s / 2.0;
}

点|多边形

相交 on_poly

P 与多边形 Poly 的关系。

  • 在外部:返回 0
  • 在内部:返回 1
  • 在边界:返回 2
cpp
int on_poly(const Point& p, const Poly& poly) {
	int counter = 0;
	int n = poly.size();
	
	for (int i = 0; i < n; i ++) {
		const Point& a = poly[i];
		const Point& b = poly[(i + 1) % n];
		
		if (on_seg(p, a, b))
			return 2;
		
		if (p.y >= min(a.y, b.y) and p.y < max(a.y, b.y)) {
			var product = a.x + (p.y - a.y) / (b.y - a.y) * (b.x - a.x);
			counter += sign(product - p.x) > 0;
		}
	}
	
	return counter % 2;
}