Yharim
Snow is on the ground
So NVIDIA, Fuck You🖕
发布于 2024-5-21

计算几何

人一定要有梦想。

精度处理 #

精度设置 #

精度的最小误差 eps 通常在 $10^{-8}$ 以内，视题目要求改变。

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

符号函数 sign #

$$\DeclareMathOperator{\sign}{sign} \DeclareMathOperator{\eps }{eps} \sign(x)= \begin{cases} -1,&x<-\eps\\ 0,&-\eps\le x\le \eps\\ 1,&x>\eps \end{cases}$$

inline int sign(var x) {
	return x < -eps ? -1 : (x > eps ? 1 : 0);
}

绝对值函数 absolute #

$$\DeclareMathOperator{\absolute}{absolute} \absolute(x)=x\cdot\sign(x)$$

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

符号约定 #

点 Point #

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

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

在精度范围内判断相等。

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

标准输入输出支持。

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 的定义。

using Vector = Point;

多边形 Poly #

多边形 Poly 由一组 点 Point 构成。

using Poly = vector<Point>;

标准输入输出支持。

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 +/- #

$$(x_1,y_1)\pm(x_2,y_2)=(x_1\pm x_2,y_1\pm y_2)$$

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)\cdot k=(x\cdot k,y\cdot k)$$

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

点乘 dot #

$$(x_1,y_1)\cdot(x_2,y_2)=x_1\cdot x_2 + y_1\cdot y_2$$

var dot(const Vector& a, const Vector& b) {
	return a.x * b.x + a.y * b.y;
}

叉乘 cross #

$$(x_1,y_1)\times(x_2,y_2)=x_1\cdot y_2 - x_2\cdot y_1$$

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

取模 module #

$$\left\Vert(x, y)\right\Vert=\sqrt{x^2+y^2}$$

$$\left\Vert\vec{p}\right\Vert=\sqrt{\vec{p}\cdot\vec{p}}$$

var module(const Vector& p) {
	return sqrt(dot(p, p));
}

旋转 rotate #

向量 $\vec{p}$ 逆时针旋转 $\theta$ 弧度。

$$(x, y)\cdot \begin{pmatrix} \cos\theta & \sin\theta \\ -\sin\theta & \cos\theta \end{pmatrix}$$

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$ 的向量。

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

点|点 #

距离 to_point # 

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

旋转 rotate #

点 $P$ 绕点 $C$ 逆时针旋转 $\theta$ 弧度。

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

点|线段 #

相交 on_seg #

点 $P$ 是否在线段 $AB$ 上。

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$ 的距离。

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$ 上。

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

距离 to_line #

点 $P$ 到直线 $AB$ 的距离。

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$ 的垂足。

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$ 的对称点。

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

线|线 #

直线|直线 交点 cross_point #

直线 $AB$ 和直线 $CD$ 的交点。

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$ 是否相交。

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$ 是否相交。

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 #

逆时针为正向。

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$。
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;
}
int on_poly_binary(const Point& p, const Poly& poly) {
	
	static auto check = [&](const Point& a, const Point& l, const Point& r) {
		return sign(cross(l - a, r - a)) > 0;
	};
	
	int n = poly.size();
	
	if (check(poly[0], p, poly[1]) or check(poly[0], poly[n-1], p))
		return 0;
	if (on_seg(p, poly[0], poly[1]) or on_seg(p, poly[0], poly[n-1]))
		return 2;
	
	int l = -1, r = n;
	while (l + 1 < r) {
		int m = (l + r) / 2;
		if (check(poly[0], poly[m], p))
			l = m;
		else
			r = m;
	}
	
	if (check(poly[l], p, poly[l+1]))
		return 0;
	if (on_seg(p, poly[l], poly[l+1]))
		return 2;
	
	return 1;
}