我一直以来都挺不想学这个玩意儿的......html
奈何最近常常碰到凸包啊半平面交之类的题,还有平面图node
因此仍是作一下吧算法
平面向量及其坐标表示数据结构
平面几何图形的基础定义、公理、定理.net
向量是个很是方便的东西,能够把不少平面几何空间几何里面用笛卡尔坐标暴力算很麻烦的东西变得很简单,因此必定要熟练运用code
如下约定向量随着字母表的顺序a...z,依次对应坐标下标1...26htm
坐标表示:$\mathbf{a}=(x_1,y_1)$blog
向量的模:就是长度,$|\mathbf{a}|=\sqrt{x_1^2+y_1^2}$排序
$\mathbf{a}\ast\mathbf{b}=x_1x_2+y_1y_2$队列
能够用来算两个直线的夹角:求出两个直线$a,b$的向量$\mathbf{a},\mathbf{b}$
夹角$cos<a,b>=\frac{\mathbf{a}\ast\mathbf{b}}{|\mathbf{a}||\mathbf{a}|}$
实际上叉积这个概念是定义于空间向量里面的:
对于空间向量$\mathbf{a}={x_1,y_1,z_1}$,
空间向量叉积:$\mathbf{a}\times\mathbf{b}=(y_1z_2-y_2z_1,x_2z_1-x_1z_2,x_1y_2-x_2y_1)$
注意这个式子获得的是一个空间向量
那么对于平面向量来讲,咱们把$z$当作0,能够获得一个实数结果:
平面向量叉积:$\mathbf{a}\times\mathbf{b}=x_1y_2-x_2y_1$
这玩意儿有什么用呢?
它等于从$\mathbf{a}$旋转到$\mathbf{b}$的过程当中构成的有向平行四边形**的面积(也就是能够是负数)
这里逆时针旋转为正,顺时针为负
那么它能够算什么呢?咱们后面再说
直线有不少种存储方法:
能够存储直线上两个端点,常见于半平面交中
能够对于直线方程$Ax+By+C=0$存储$A,B,C$,这种比较不直观,相对来讲用的少一点
能够对于直线方程$y=kx+b$存储$k,b$,分别是斜率和$y$轴截距,这种比较直观,常见于斜率相关的处理中
固然了,用哪种最好仍是依照各位本身的习惯,适合本身的才是最好的
最直观的的方式通常是求出$k$和$b$,用直线上两点$A,B$列方程
$k=\frac{y_2-y_1}{x_2-x_1}$
$b=y_1-x_1\ast k$
这个你们都学过,用直线方程$Ax+By+C=0$直接求
$dis=|\frac{Ax_0+By_0+C}{\sqrt{A^2+B^2}}|$
也能够用直线上两点坐标,构成三角形,求高便可
这里咱们就会看到平面向量叉积的第一个使用:用从向量(这里也就是有向直线)上任意两点指向目标点的两个向量的叉积的正负,能够判断点和直线的位置关系
更准确地:设有向直线上按照顺序排列的两个点$A,B$,目标点$C$,则$\overrightarrow{AC}\times\overrightarrow{BC}$为正则$C$在$\overrightarrow{AB}$右侧,不然在左侧
能够画个图体会一下
若是有直线方程的话直接解方程便可
这里讲一个用两个直线上四个点求交点的方法:
设两条直线上四点分别为$A,B$和$C,D$,坐标为$(x_{1...4},y_{1...4})$
咱们把点变成位置向量
令$v1=(A-D)\times(B-D),v2=(A-C)\times(B-C)$,那么交点的位置向量(位矢)为$D+(v1/(v1-v2))*(C-D)$
这个东西也是画个图就出来了:求出的两个平行四边形面积之比正好等于交点到$C,D$的距离之比,只不过其中一个是正的一个是负的,因此那里是$(v1-v2)$
代码以下:
struct p{ long double x,y; p(long double xx=0.0,long double yy=0.0){x=xx;y=yy;} }; inline p operator *(const p &a,const long double &b){return p(a.x*b,a.y*b);} inline long double operator *(const p &a,const p &b){return a.x*b.y-a.y*b.x;}//'x-multiple' of planary vector inline p operator -(const p &a,const p &b){return p(a.x-b.x,a.y-b.y);} inline p operator +(const p &a,const p &b){return p(a.x+b.x,a.y+b.y);} struct seg{ p a,b;long double k; seg(p aa=p(),p bb=p(),long double kk=0.0){a=aa;b=bb;k=kk;} }; inline p cross(const seg &x,const seg &y){//calculate the intersection using planary vector long double v1=(x.a-y.b)*(x.b-y.b); long double v2=(x.a-y.a)*(x.b-y.a); long double c=v1/(v1-v2); p re=(y.b+((y.a-y.b)*c)); return re; }
基本上是很好理解的:一堆直线的右边那一半平面的交就是半平面交【听起来贼好理解是否是】
操做也比较简单,放一个dalao的连接在这里(懒得本身写了23333)
1.以逆时针为正方向,建边(输入方向不肯定时,可用叉乘求面积看正负得知输入的顺逆方向)
2.对线段根据极角排序
3.去除极角相同的状况下,位置在右边的边
4.用双端队列储存线段集合,遍历全部线段
5.判断该线段加入后对半平面交的影响(对双端队列的头部和尾部进行判断,由于线段加入是有序的)
(这里判断的方式是看最前面两个或者最后面两个的交点是否是在新加入线的右边,能够画个图理解一下)
6.若是某条线段对于新的半平面交没有影响,则从队列中剔除掉(双端队列头尾删除)
7.最后剩下的线段集合,即便最后要求的半平面交
代码给一下:
这份代码用双端队列求出半平面交的直线集合,以及集合内直线的交点,再用交点位矢求出面积
这个代码有问题,看下面的吧
update 2019/4/6
前两天写了[HNOI2012]射箭这道题,从新认识了半平面交的写法
具体参考上面博客,这里放一个比较全的模板代码
必定要看注释!!!!!!!!!!!!
inline bool sign(long double x){//判断符号 if(x>eps) return 1; if(x<-eps) return -1; return 0; } int n,m; struct node{ long double x,y; node(long double xx=0.0,long double yy=0.0){x=xx;y=yy;} //这里是向量的基本运算,注意这是个通用的数据结构,点和二维向量都集成进去了 //注意到点和二维向量在作大部分运算的时候是同样的,因此这么作可行 //标*的是数乘和平面向量叉乘,标/的是平面向量点乘 //slope是求两个点之间的斜率 inline friend node operator +(const node &a,const node &b){return node(a.x+b.x,a.y+b.y);} inline friend node operator -(const node &a,const node &b){return node(a.x-b.x,a.y-b.y);} inline friend node operator *(const node &a,const long double &b){return node(a.x*b,a.y*b);} inline friend long double operator *(const node &a,const node &b){return a.x*b.y-a.y*b.x;} inline friend long double operator /(const node &a,const node &b){return a.x*b.x+a.y*b.y;} inline friend long double slope(const node &a,const node &b){return atan2l(a.y-b.y,a.x-b.x);} }rt[300010]; struct seg{ //这里线段(直线、半平面)的定义是从a开始到b结束的向量,是有方向的 node a,b;long double k;int id; seg(node aa=node(),node bb=node()){a=aa;b=bb;k=slope(aa,bb);id=0;} seg(node aa,node bb,long double kk){a=aa;b=bb;k=kk;id=0;} inline friend bool operator <(const seg &a,const seg &b){return a.k<b.k;}//按照斜率排序 inline friend node cross(const seg &a,const seg &b){//求两个线段的交点,讲解见上面 long double v1=(a.a-b.b)*(a.b-b.b); long double v2=(a.a-b.a)*(a.b-b.a); return b.b+(b.a-b.b)*(v1/(v1-v2)); } inline friend bool right(const node &a,const seg &b){ //判断一个点是否是在一条线的右边 //注意这里是大于eps,由于半平面交可能出现最后的交是一个点的状况 //有的题目须要排除上述状况,就写大于-eps(这样包括了点在线段上) return ((a-b.b)*(a-b.a))>eps; } }lis[300010],a[300010],q[300010]; inline bool solve(int lim){//重要!!!!!这份代码的半平面交是每一个有向直线(seg)的左侧半平面的交 int i,head=1,tail=0,flag,tot=0; for(i=1;i<=m;i++) if(lis[i].id<=lim) a[++tot]=lis[i]; for(i=1;i<=tot;i++){ flag=0; while((head<=tail)&&(!sign(a[i].k-q[tail].k))){ if(right(q[tail].a,a[i])) tail--; else{flag=1;break;} } if(flag) continue; while(head<tail&&right(rt[tail],a[i])) tail--; while(head<tail&&right(rt[head+1],a[i])) head++; q[++tail]=a[i]; if(head<tail) rt[tail]=cross(q[tail],q[tail-1]); } while(head<tail&&right(rt[tail],q[head])) tail--; while(head<tail&&right(rt[head+1],q[tail])) head++; return (tail-head>1); } const long double pi=acos(-1.0); int main(){ //这里是边界条件,这份模板里面要加入 //这里也能够看出来求的是线段左侧的半平面的交 lis[++m]=seg(node(-1e12,1e12),node(-1e12,-1e12),pi/2.0); lis[++m]=seg(node(1e12,1e12),node(-1e12,1e12),0); lis[++m]=seg(node(1e12,-1e12),node(1e12,1e12),-pi/2.0); lis[++m]=seg(node(-1e12,-1e12),node(1e12,-1e12),pi); }
咕咕咕(主要是我没见过纯凸包的......)周末补
我错了,凸包真是博大精深
凸包能够简单地理解为一条绳子,绕在一个点集(木桩集合)的最外面,造成的凸多边形
凸包有以下性质:
1.全部点都在凸包内部或者凸包上
2.凸包的端点必定是点集中的点
求凸包经常使用graham算法,时间复杂度$O(n\log n)$
流程以下:
找到$y$坐标最小的一点做为原点
对原点以外的全部点按照到原点的极角排序(这里由于选取了最靠下的,因此极角范围在$[0,\pi]$)
依次遍历全部排序后的点,加入一个单调栈中:每次判断(栈顶元素和栈顶第二元素之间的斜率)是否大于(当前点和栈顶第二元素之间的斜率)
注意一旦这个大于成立了,栈顶元素就会在当前元素和栈顶第二元素的连线的“下面”,也就是在凸包里面了
由于咱们事先按照极角排序了,因此这一单调栈能够不重复不遗漏地记录凸包上全部点
注意这样求出来的凸包上的点是逆时针排序的(根本缘由是由于极角排序就是逆时针绕圈)
示例代码:
struct node{ double x,y; node(double xx=0.0,double yy=0.0){x=xx;y=yy;} inline bool operator <(const node &b){return ((fabs(y-b.y)<eps)?(x<b.x):(y<b.y));} inline friend bool operator ==(const node &a,const node &b){return ((fabs(a.x-b.x)<eps)&&(fabs(a.y-b.y)<eps));} inline friend bool operator !=(const node &a,const node &b){return !(a==b);} inline friend node operator +(const node &l,const node &r){return node(l.x+r.x,l.y+r.y);} inline friend node operator -(const node &l,const node &r){return node(l.x-r.x,l.y-r.y);} inline friend node operator *(node l,double r){return node(l.x*r,l.y*r);} inline friend double operator *(const node &l,const node &r){return l.x*r.y-l.y*r.x;} inline friend double operator /(const node &l,const node &r){return l.x*r.x+l.y*r.y;} inline friend double dis(const node &a){return sqrt(a.x*a.x+a.y*a.y);} }a[100010],q[100010],x[10]; inline bool cmp(node l,node r){ double tmp=(a[1]-l)*(a[1]-r); if(fabs(tmp)<eps) return dis(a[1]-l)<dis(a[1]-r); else return tmp>0; } void graham(){//get a counter-clockwise convex int i; for(i=2;i<=n;i++){ if(a[i]<a[1]) swap(a[1],a[i]); } sort(a+2,a+n+1,cmp); q[++top]=a[1]; q[++top]=a[2]; for(i=3;i<=n;i++){ while(top>1&&((q[top]-q[top-1])*(a[i]-q[top])<eps)) top--; q[++top]=a[i]; } q[0]=q[top]; }
旋转卡壳,顾名思义,就是“旋转着”+“卡qia在凸壳qiao上”
简单来说,就是拿平行直线卡在凸包外面,而后让凸包在里面转,这样的感受
旋转卡壳能够解决凸包最大直径、凸包最小直径(凸包宽)、两个凸包之间最大最小距离等问题
它也能够解决凸包外接最小面积、最小周长凸$n$边形问题
详细的讲解请戳这里
旋转卡壳+graham的一道例题:BZOJ1185
凸包有很是很是很是多的应用
斜率DP中用到的就是上下凸壳的单调栈求法
对于决策最大化问题,经常能够经过考虑不一样决策之间的关系,来导出凸包问题或者半平面交问题
也能够在求最优决策的时候,把答案带入式子中,并最大化包含答案那一项
和半平面交一块儿使用的时候,一般是两种:
1.断定形问题,断定是否点集在某个半平面内,此时在直线上方须要下凸壳,直线下方须要上凸壳
2.极值形问题,求出点集关于半平面的某个极值,此时在直线上方须要上凸壳,直线下方须要下凸壳(最大化),最小化的方法同类型一
线段树能够维护区间询问的凸包,可是要考虑到凸包上传更新的时间复杂度问题
cdq分治能够解决原本须要动态凸包的问题,比较好写好调
纯动态凸包可使用set或者平衡树维护,例题在这里
所谓平面图,就是一个图,全部的边能够画在平面上,互相之间不在定点以外的地方相交
通常会给出平面图每一个点的坐标、每一个边是直线
能够看到,一个平面图会把一个平面划分红一个无限大的区域和若干面积有限的区域
咱们可使用最左转线算法,把每个区域变成点、每个原图边变成新图中相邻两个区域之间的边,这样就构建出了平面图的对偶图
先把全部边改为双向的(不过通常题目都会给双向的)
对于原图每个点的出边,按照其出边的斜率进行排序。
以以下方式遍历图中全部的边:
找一条没有访问过的边(u,v),标记之
对于点v,找到v的出边中极点序在(v,u)后面的那一个,并重复上一行的过程,直到某一次找到的下一条边是标记过的
每作一次这个过程,咱们就会找到一个 对偶图定点(亦即原图中的一个区域)
一次过程当中访问到的全部边就是这个区域的边界(若是是无限大的那个区域,就是分割无限大和其余人的全部边)
由于咱们全部的边都是双向的,而每一条双向边分割了两个区域,因此最终咱们对于原图边和区域的一一对应,能够不重复不遗漏
由于“找到极点序中的上一个”这个操做就是找到这条边的下一条中“最向左转”的一条边,因此算法称做最左转线
若想要图片解释,请戳这里
贴个代码,来自这道题:
for(i=1;i<=n;i++){ tot=e[i].size();ee.clear(); for(j=0;j<tot;j++){ ee.push_back(mp(atan2(y[e[i][j].first]-y[i],x[e[i][j].first]-x[i]),e[i][j].first)); } sort(ee.begin(),ee.end()); for(j=0;j<tot-1;j++){//最左转线预处理:标记每个点的后继 suf[i][ee[j].second]=ee[j+1].second; } if(tot) suf[i][ee[tot-1].second]=ee[0].second; } for(i=1;i<=n;i++){ for(j=0;j<e[i].size();j++){ pos=e[i][j].first;from=i; if(col[i][pos]) continue; cnt++; col[i][pos]=cnt; suml[cnt]+=light[pos]; sumd[cnt]+=dark[pos]; while(1){//求出一个区域 next=suf[pos][from]; if(col[pos][next]) break; from=pos;pos=next; col[from][pos]=cnt; suml[cnt]+=light[pos]; sumd[cnt]+=dark[pos]; } } }
点定位,顾名思义,就是定位平面上的点位于平面图分割出来的那个区域里面
能够用排序+平衡树实现$O(n\log n)$,也是周末补上