最近去hust oj 上挂了一套题目计算几何,里面有几道半平面交的题,于是学了一下半平面交。。
所谓半平面,通俗一点就是就是一些形如 a*x+b*y+c <=0 (或 >=0)的区域,实际上就是平面的一半,成为半平面,
而半平面交,就是一些这样的半平面解的交集,反映上直角坐标系上就是半平面的共同区域
而求半平面的的方法一般有3种,
1)联机算法,这是一种朴素的做法,即对新加入的边判断原来的交点是否在合法的区域内,不合法就该点,并加入新点。。
最终的点就是合法的点。复杂度o(n^2)
2)分治法。即利用二分的思想 将 N条边划分成2个n/2,利用递归求解。算法复杂度o(nlogn)
3)增量排序法。此方法是zzy大神2006国家集训队论文中提出的。。下面具体操作
Step1:对于所有的直线进行极角排序,如果极角相同,保留需要的一个(这几个极角相同的全部为公共解的那一个)
step2:使用双端队列,加入初始的两个半平面
step3:依次加入新的半平面,并判断原来队尾的2个半平面交点符不符合当前的半平面,不符合删除队尾的最后一个元素
判断原来队头的2个半平面交点符不符合当前的半平面,不符合删除队头的第一一个元素
将半平面加到队尾,重复step3操作。。
step4:判断当前队头跟队尾有木有矛盾(有矛盾就是队头的前两个元素交点不再队尾班平面区域内),删除队头的第一个元素
队尾同样操作。。
step5,判断有木有解。。(一般判断剩下的直线,少于3条就无解)
解决的问题
求多边形的核:就是求看多边形里面存不存在可以看到边上所有点的的点
poj3335 poj2451等
自己写了一下,感觉很蛋疼。。就去看了别人模板。。
Orz Non_Cease写得十分优美易懂的模板
模板(poj1279):
1 #include<iostream> 2 #include<cstring> 3 #include<cstdlib> 4 #include<cstdio> 5 #include<cmath> 6 #include<algorithm> 7 8 #define maxn 2005 9 #define eps 1e-9 10 using namespace std; 11 struct point{ double x,y; } p[maxn]; //点 12 struct line{ 13 point a, b; 14 double angle; 15 } l[maxn], l1[maxn]; //直线 16 17 int ln, ord[maxn], T, q[maxn], tot, n, bot, top;// ord数组为直线标号 18 19 void add_line(double x1, double y1, double x2, double y2){//加入边,有向的,有向的话保证他的某一方向为半平面 20 l[ln].a.x = x1; 21 l[ln].a.y = y1; 22 l[ln].b.x = x2; 23 l[ln].b.y = y2; 24 l[ln].angle = atan2(y2-y1, x2-x1); 25 ord[ln] = ln; 26 ++ln; 27 } 28 29 void init(){ 30 scanf("%d",&n); 31 for (int i = 0; i < n; ++i) 32 scanf("%lf%lf", &p[i].x, &p[i].y); 33 p[n] = p[0]; 34 ln = 0; 35 for (int i = 0; i < n; ++i) 36 add_line(p[i+1].x, p[i+1].y, p[i].x, p[i].y); 37 } 38 39 int zero(double k){ 40 if (fabs(k) < eps) return 0; 41 return k < 0? -1 : 1; 42 } 43 44 double multi(point p0, point p1, point p2){//叉积,判断或求面积 45 return (p1.x-p0.x)*(p2.y-p0.y) - (p2.x-p0.x)*(p1.y-p0.y); 46 } 47 48 bool cmp(const int u, const int v){//排序 49 double d = l[u].angle - l[v].angle; 50 if (d == 0) return zero(multi(l[u].a, l[u].b, l[v].b)) < 0; 51 return d < 0; 52 53 } 54 55 void get_point(line l1, line l2, point &p){//求交点,采用比例求,而比例直接比较面积 56 double dot1 = multi(l1.a, l1.b, l2.a); 57 double dot2 = multi(l1.a, l2.b, l1.b); 58 p.x = (l2.a.x*dot2 + l2.b.x*dot1) / (dot1 +dot2); 59 p.y = (l2.a.y*dot2 + l2.b.y*dot1) / (dot1 +dot2); 60 } 61 62 bool judge(line l0, line l1, line l2){//判断l1,l2交点是否在l0的解范围内 63 point p; 64 get_point(l1, l2, p); 65 return zero(multi(l0.a, l0.b , p)) < 0; 66 67 } 68 void SAI(){ 69 sort(ord, ord + ln, cmp); 70 /*for (int i = 0; i < ln; ++i){ 71 printf("l[%d].angle = %0.2lf\n ", ord[i], l[ord[i]].angle); 72 73 } */ 74 int j = 0; 75 for (int i = 1; i < ln; ++i) //去除相同极角的,留下符合的 76 if (zero(l[ord[i]].angle - l[ord[j]].angle) > 0) 77 ord[++j] = ord[i]; 78 ln = ++j; 79 q[0] = ord[0]; 80 q[1] = ord[1]; 81 bot = 0; 82 top = 1; 83 for (int i = 2; i < ln; ++i){ //加入半平面 84 while (top > bot && judge(l[ord[i]], l[q[top - 1]], l[q[top]])) --top; 85 while (top > bot && judge(l[ord[i]], l[q[bot + 1]], l[q[bot]])) ++bot; 86 q[++top] = ord[i]; 87 } 88 while (top > bot && judge(l[q[bot]], l[q[top - 1]], l[q[top]])) --top; 89 while (top > bot && judge(l[q[top]], l[q[bot + 1]], l[q[bot]])) ++bot; 90 91 q[++top] = q[bot]; 92 tot = 0; 93 for (int i = bot; i < top; ++i) 94 get_point(l[q[i]], l[q[i+1]], p[tot++]); 95 96 } 97 98 double squre(){ 99 if (tot < 2) return 0; 100 double ans = 0; 101 for (int i = 1; i < tot -1; ++i) //求面积 102 ans += multi(p[0], p[i], p[i+1]); 103 return fabs(ans / 2); 104 105 } 106 107 void solve(){ 108 SAI(); 109 printf("%.2lf\n", squre()); 110 } 111 112 int main(){ 113 freopen("poj1279.in","r", stdin); 114 freopen("poj1279.out","w",stdout); 115 scanf("%d", &T); 116 for (int i = 1; i <= T; ++i){ 117 init(); 118 solve(); 119 } 120 fclose(stdin); fclose(stdout); 121 }