半平面交(Half-plane Intersection, HPI)是计算几何中的经典问题。给定一组由有向直线定义的半平面,求这些半平面的交集区域。该算法在可见性计算、机器人路径规划、线性规划等领域有广泛应用。本文用Java实现基于排序的O(n log n)半平面交算法,讲解极角排序、双端队列维护与交点计算。
问题定义与几何基础
半平面:由直线 ax + by + c = 0 将平面分为两部分,有向直线左侧(或右侧)的区域即为一个半平面。
半平面交:多个半平面的交集可能是一个凸多边形、凸无界区域、直线、线段、点或空集。
关键几何工具
/**
* 二维点/向量
*/
class Point {
double x, y;
Point(double x, double y) { this.x = x; this.y = y; }
}
/**
* 有向直线:点p指向点q,左侧为保留的半平面
*/
class Line {
Point p, q; // 直线上两点
double angle; // 极角,用于排序
Line(Point p, Point q) {
this.p = p;
this.q = q;
this.angle = Math.atan2(q.y - p.y, q.x - p.x);
}
}
算法核心思想
基于排序增量法(S&I算法):
- 极角排序:将所有半平面的边界直线按极角排序
- 双端队列维护:用双端队列维护当前半平面交的边界直线
- 冗余判断:新加入的直线可能导致队列首尾直线不再必要,及时弹出
- 交点计算:相邻直线的交点构成多边形的顶点
核心观察:半平面的交集是凸区域,其边界由若干直线段组成,这些直线在原直线集合中按极角有序排列。
Java完整实现
import java.util.*;
/**
* 半平面交(Half-plane Intersection)
* 基于排序增量法,时间复杂度 O(n log n)
*/
public class HalfplaneIntersection {
static final double EPS = 1e-9;
/**
* 二维点/向量,支持基本向量运算
*/
static class Point {
double x, y;
Point(double x, double y) { this.x = x; this.y = y; }
Point() { this(0, 0); }
Point add(Point p) { return new Point(x + p.x, y + p.y); }
Point sub(Point p) { return new Point(x - p.x, y - p.y); }
Point mul(double k) { return new Point(x * k, y * k); }
double cross(Point p) { return x * p.y - y * p.x; } // 叉积
double dot(Point p) { return x * p.x + y * p.y; } // 点积
double len() { return Math.sqrt(x * x + y * y); }
}
/**
* 有向直线:从点p到点q,保留的是直线的左侧半平面
*/
static class Line implements Comparable<Line> {
Point p, q; // 直线上两点
double angle; // 极角
Line(Point p, Point q) {
this.p = p;
this.q = q;
this.angle = Math.atan2(q.y - p.y, q.x - p.x);
}
@Override
public int compareTo(Line other) {
// 按极角排序,极角相同的保留最靠左的
if (Math.abs(angle - other.angle) < EPS) {
// 极角相同,保留左侧更严格的(点在对方左侧)
return onLeft(other.p) ? -1 : 1;
}
return Double.compare(angle, other.angle);
}
/**
* 判断点是否在直线左侧(包含在直线上)
*/
boolean onLeft(Point r) {
return q.sub(p).cross(r.sub(p)) > -EPS;
}
}
/**
* 计算两条直线的交点
* 假设两条直线不平行
*/
static Point intersection(Line a, Line b) {
Point u = a.p.sub(b.p);
Point v = a.q.sub(a.p);
Point w = b.q.sub(b.p);
double t = w.cross(u) / v.cross(w);
return a.p.add(v.mul(t));
}
/**
* 平行判断
*/
static boolean parallel(Line a, Line b) {
return Math.abs(a.q.sub(a.p).cross(b.q.sub(b.p))) < EPS;
}
/**
* 半平面交主算法
* @param lines 输入的有向直线数组(左侧为保留区域)
* @return 交集多边形的顶点列表(逆时针顺序),空集返回空列表
*/
public static List<Point> halfplaneIntersection(Line[] lines) {
int n = lines.length;
// 1. 按极角排序,去除极角相同的直线(保留最严格的)
Arrays.sort(lines);
Line[] filtered = new Line[n];
int m = 0;
for (int i = 0; i < n; i++) {
// 跳过与上一条极角相同的直线
if (i > 0 && Math.abs(lines[i].angle - lines[i - 1].angle) < EPS) {
continue;
}
filtered[m++] = lines[i];
}
// 2. 双端队列维护半平面交
Line[] dq = new Line[m + 5];
Point[] points = new Point[m + 5]; // 相邻直线交点
int head = 0, tail = 0;
for (int i = 0; i < m; i++) {
Line cur = filtered[i];
// 队尾冗余判断:新直线使队尾两条直线交点不在新半平面内
while (tail - head > 1 && !cur.onLeft(points[tail - 1])) {
tail--;
}
// 队首冗余判断
while (tail - head > 1 && !cur.onLeft(points[head + 1])) {
head++;
}
// 如果新直线与队尾直线平行且方向相同,则两者取更严格的
if (tail > head && parallel(dq[tail - 1], cur)) {
if (!cur.onLeft(dq[tail - 1].p)) {
// 当前直线更严格,替换队尾
tail--;
} else {
// 队尾更严格,跳过当前
continue;
}
}
dq[tail++] = cur;
// 计算新加入直线与队尾前一条直线的交点
if (tail - head > 1) {
points[tail - 1] = intersection(dq[tail - 2], dq[tail - 1]);
}
}
// 3. 最终处理:队首和队尾可能仍有冗余
while (tail - head > 1 && !dq[head].onLeft(points[tail - 1])) {
tail--;
}
while (tail - head > 1 && !dq[tail - 1].onLeft(points[head + 1])) {
head++;
}
// 4. 收集多边形顶点
if (tail - head <= 2) {
return new ArrayList<>(); // 空集或面积为零
}
List<Point> polygon = new ArrayList<>();
for (int i = head + 1; i < tail; i++) {
polygon.add(points[i]);
}
// 队首和队尾的交点
polygon.add(intersection(dq[tail - 1], dq[head]));
return polygon;
}
/**
* 计算多边形面积(鞋带公式)
*/
public static double polygonArea(List<Point> polygon) {
if (polygon.size() < 3) return 0;
double area = 0;
int n = polygon.size();
for (int i = 0; i < n; i++) {
int j = (i + 1) % n;
area += polygon.get(i).x * polygon.get(j).y;
area -= polygon.get(j).x * polygon.get(i).y;
}
return Math.abs(area) / 2.0;
}
/**
* 将一般不等式 ax + by <= c 转换为有向直线
* 保留的是不等式表示的半平面
*/
static Line fromInequality(double a, double b, double c) {
// ax + by <= c
// 找直线上两点
Point p, q;
if (Math.abs(b) < EPS) {
// 垂直线: x = c/a
double x = c / a;
p = new Point(x, 0);
q = new Point(x, 1);
} else {
// y = (c - ax) / b
p = new Point(0, c / b);
q = new Point(1, (c - a) / b);
}
// 确保保留的是 ax + by <= c 的区域
// 测试点(0,0):若 a*0 + b*0 <= c,即 c >= 0,则(0,0)在区域内
// 需要确保有向直线的左侧是满足不等式的区域
Line line = new Line(p, q);
Point test = new Point(0, 0);
boolean zeroInRegion = a * 0 + b * 0 <= c + EPS;
boolean zeroOnLeft = line.onLeft(test);
if (zeroInRegion != zeroOnLeft) {
// 方向反了,交换p和q
return new Line(q, p);
}
return line;
}
public static void main(String[] args) {
// 示例1:求一个凸多边形(正方形)
System.out.println("=== 示例1:正方形半平面交 ===");
Line[] square = {
new Line(new Point(0, 0), new Point(1, 0)), // y >= 0(上方)
new Line(new Point(1, 0), new Point(1, 1)), // x <= 1(左侧)
new Line(new Point(1, 1), new Point(0, 1)), // y <= 1(下方)
new Line(new Point(0, 1), new Point(0, 0)) // x >= 0(右侧)
};
// 注意:Line(p,q)保留的是p->q的左侧
// 需要调整方向使内部在左侧
Line[] squareFixed = {
new Line(new Point(1, 0), new Point(0, 0)), // y >= 0
new Line(new Point(1, 1), new Point(1, 0)), // x <= 1
new Line(new Point(0, 1), new Point(1, 1)), // y <= 1
new Line(new Point(0, 0), new Point(0, 1)) // x >= 0
};
List<Point> result1 = halfplaneIntersection(squareFixed);
System.out.println("交点数量: " + result1.size());
for (Point p : result1) {
System.out.printf("(%.2f, %.2f)%n", p.x, p.y);
}
System.out.printf("面积: %.4f%n", polygonArea(result1));
// 示例2:使用不等式形式定义半平面
System.out.println("\n=== 示例2:不等式形式 ===");
// 定义三角形区域:
// x >= 0, y >= 0, x + y <= 1
Line[] triangle = {
fromInequality(-1, 0, 0), // -x <= 0 即 x >= 0
fromInequality(0, -1, 0), // -y <= 0 即 y >= 0
fromInequality(1, 1, 1) // x + y <= 1
};
List<Point> result2 = halfplaneIntersection(triangle);
System.out.println("三角形顶点:");
for (Point p : result2) {
System.out.printf("(%.4f, %.4f)%n", p.x, p.y);
}
System.out.printf("三角形面积: %.4f%n", polygonArea(result2));
// 示例3:无解情况(矛盾的半平面)
System.out.println("\n=== 示例3:无解情况 ===");
Line[] impossible = {
new Line(new Point(0, 0), new Point(1, 0)), // y >= 0
new Line(new Point(0, 1), new Point(1, 1)) // y <= 0(矛盾,内部方向需调整)
};
// 实际上需要正确设置方向...
Line[] impossibleFixed = {
new Line(new Point(1, 0), new Point(0, 0)), // y >= 0
new Line(new Point(1, -1), new Point(0, -1)) // y <= -1(与y>=0矛盾)
};
List<Point> result3 = halfplaneIntersection(impossibleFixed);
System.out.println("交点数量(应为0): " + result3.size());
}
}
算法正确性分析
双端队列的维护基于以下关键观察:
- 极角单调性:按极角排序后,相邻直线在凸多边形上也是相邻的
- 冗余判定:如果队列中某条直线的两个相邻交点都在新半平面外,则该直线不再构成边界
- 凸性保证:半平面的交集是凸集,其边界由按极角有序的直线段组成
复杂度分析
| 步骤 | 时间复杂度 | 说明 |
|---|---|---|
| 极角排序 | O(n log n) | 主要开销 |
| 双端队列处理 | O(n) | 每条直线入队出队各一次 |
| 交点计算 | O(n) | 相邻直线求交 |
| 总复杂度 | O(n log n) | 由排序决定 |
| 空间复杂度 | O(n) | 存储直线和交点 |
应用场景
半平面交在多个领域有重要应用:
- 线性规划的可行域计算:求约束条件的交集
- 机器人可见性:计算光源或摄像头的可见区域
- 运动规划:计算机器人可移动的安全区域
- 计算机图形学:裁剪多边形的Weiler-Atherton算法的辅助工具
- 核(Kernel)计算:计算点集的核,即可以看到所有点的区域
注意事项
- 精度问题:计算几何中浮点精度至关重要,使用EPS进行相等判断
- 边界情况:空集、无界区域、退化(直线/点)需要特殊处理
- 方向一致性:确保所有半平面的”保留侧”定义一致(本文统一使用左侧)
- 无界区域处理:本文实现返回有界多边形,处理无界区域需要添加一个大包围盒
掌握半平面交后,可以解决大量计算几何中的区域求交问题,是算法竞赛中的高级但实用的工具。