二维计算几何学习笔记

最新推荐文章于 2021-10-04 20:58:21 发布

Tmotfl

最新推荐文章于 2021-10-04 20:58:21 发布

阅读量323

点赞数 2

本文链接：https://blog.csdn.net/NLt_XXy/article/details/80057278

版权

学习笔记同时被 2 个专栏收录

11 篇文章 0 订阅

订阅专栏

—————计算几何————.....

2 篇文章 0 订阅

订阅专栏

1.基础知识

1.1关于浮点误差

一般的话要设一个 $EPS=1e^{-8}$

//精度控制 
bool Dcmp(double x)
{
    if(fabs(x)<Eps) return 0;
    else return 1;
}

1.2向量基础运算

1.2.1点积

double Dot(Vector A,Vector B) { return A.x*B.x+A.y*B.y; }

点积可以用来判断线段垂直什么的.

1.2.2叉积

double Cross(Vector A,Vector B){ return A.x*B.y-A.y*B.x; }

叉积表示两个向量组成的平行四边形的有向面积.
设两个向量 $v1,v2$
当 $v1$ 在 $v2$ 顺时针方向时，叉积大于0.
当 $v1$ 在 $v2$ 逆时针方向时，叉积小于0.
当 $v1,v2$ 共线时，叉积等于0.

1.2.3叉积点积方向.

这里写图片描述

2.进阶知识

2.1点直线线段关系

2.1.1点到直线的距离

一般的话要三个点，点一个点，直线上两个点。
直接用叉积求出面积来除以底边的长度就可以

//点到直线距离
double PointDisToLine1(Point P,Point A,Point B)
{
    Vector v1=B-A,v2=P-A;
    return fabs(Cross(v1,v2))/Length(v1);
}

2.1.2点到线段的距离

点到线段的距离分为三种情况：
红色为点到线段距离
1.当点位于线段上方时.
这里写图片描述
2.当点位于线段左边时.

3.当点位于线段右边时.

//点到线段距离 
double PointDisToLine2(Point P,Point A,Point B)
{
    if(A==B) 
        return Length(P-A);
    Vector v1=B-A,v2=P-A.v3=P-B;
    if(Dcmp(v1,v2)<0) return Length(v2);
    else if(Dcmp(v1,v3)<0) return Length(v3);
    else return fabs(Cross(v1,v2))/Length(v1); 
}

2.1.3线段相交

我们判断线段相交直接判断一个线段两个端点是否位于另一个线段的两侧就可以.
这里写图片描述
假如线段相交的话，那么 $V_{b_1a_1}$ 叉乘 $V_{a_1a_2}$ 和 $V_{b_2a_1}$ 叉乘 $V_{a_1a_2}$ 的符号一定不同,
同理另一符号也一定不同，所以我们可以通过叉积的符号来判断线段相交.

//线段相交判定
bool LineIntersection(Point a1,Point a2,Point b1,Point b2)
{
    double c1=Cross(a2-a1,b1-a1),c2=Cross(a2-a1,b2-a1),
           c3=Cross(b2-b1,a2-b1),c4=Cross(b2-b1,a2-b1);
    return Dcmp(c1)*Dcmp(c2)<0 && Dcmp(c3)*Dcmp(c4)<0; 
}

2.1.4直线交点的计算

对于直线我们设直线的方程 $l_1=P_0+tv$ ,其中 $P_0表示一个点$ ， $v$ 表示方向向量
我们设两个直线分别是 $P+tv$ 和 $Q+tw$ ，设 $u=P-Q$
则 $t_1=\frac{Cross(w,u)}{Cross(v,w)},t_2=\frac{Cross(v,u)}{Cross(v,w)}$
带入求解即可.

//直线交点的计算
Point GetLineIntersection(Point P,Vector v,Point Q,Vector w)
{
    Vector u=P-Q;
    double t=Cross(w,u)/Cross(v,w);
    return P+v*t;
}

2.2多边形

2.2.1多边形的面积

怎样计算多边形的面积呢？
我们可以把多边形来划分成数个三角形，然后求叉积就可以了.
这里写图片描述
对于上图，我们要求 $S_{n_1n_2n_3n_4n_5}$ 的话，我们只需要求出 $(n_1n_3 × n_1n_2 )+(n_1n_4 × n_1n_3)+(n_1n_5 ×n_1n_4)$ 即可
这个可以通过 $for(····)$ 来实现，实现起来比较简单.
上面的例子是凸多边形，假如我们要求的不是凸多边形呢？
如下图：
这里写图片描述
我们设 $S_{n_1n_2n_3}=S_1,S_{n_1n_4n_3}=S_2,S_{n_1n_5n_4}=S_3$
我们发现由于我们计算的时候是按照顶点的顺序来计算的，所以我们计算叉积的时候会存在正负的关系，所以对于不是凸多边形的情况，我们依然可以按照上面的方法计算比如说上图 $S=S_1-S_2+S_3$ ，我们计算 $S_2$ 的时候恰好叉积为负。
最后不要忘记除2哦,因为计算的是平行四边形的有向面积

//多边形面积
double S_ploy(Point p[],int n)
{
    double S=0;
    for(int i=3;i<=n;i++)
        S+=Cross(p[i]-p[1],p[i-1]-p[1]);
    return S/2;
}

2.2.2关于点在多边形内的判断.

一般的话采用射线法,我们只需要判断穿过多边形的数目就可以

//判断点是否在多边形内
int Point_In_Poly(Point P,Point poly[],int num)
{
    int ans=0,k,d1,d2;
    for(int i=1;i<=num;i++)
    {
        if(!Dcmp(Point_Dis_To_1(P,poly[i],poly[i%num+1]))) return -1;
        k=Dcmp(Cross(ploy[i%num+1]-poly[i],p-poly[i]));
        d1=Dcmp(ploy[i].y-p.y);
        d2=Dcmp(poly[i%num+1].y-p.y);
        if(k>0 && d1<=0 && d2>0 ) ++ans;
        if(k<0 && d1>=0 && d2<0) --ans;   
    }
    if(ans) 
        return 1;
    return 0;
}

3.高级知识

3.1凸包

3.1.1定义

凸包是啥呢？
给定二维平面上的点集，凸包就是将最外层的点连接起来构成的凸多边形，它能包含点集中所有的点.
这里写图片描述
上面红色的就是一个凸包。

3.1.2Andrew算法

$Andrew$ 是 $Graham$ 算法的变种，但是 $Andrew$ 更快.
算法流程：
1.就是先对所有点进行排序，以 $x$ 为第一关键字，以 $y$ 为第二关键字
2.按照顺序依次遍历所有点，看是否符合要求：
这里写图片描述
怎么看是否符合要求呢，我们假设 $A,B,C,D$ 已经是凸包中的点了，现在我们将要插入 $E$ ，这是我们要判断凸包中倒数第一个点与倒数第二个点组成的向量，与即将插入的点与倒数第二个点组成的向量的叉积的正负，若为负，也就是上图中的情况 $((D,C)×(E,C))<0$ ,我们就退出凸包中的最后一个点，因为这个点已经不可能在凸包内了，然后我们不停的判断直达叉积大于0，就可以把当前要插入的点插入了。
搞完这一个点之后就是这个样子了：
这里写图片描述
3.这样按照排序遍历一遍之后我们只是求得了下凸包，也就是凸包的一半，所以我们还要再倒序枚举一遍，求出上凸包来
4.凸包就求完了。

//凸包
int ConvexHull(Point P[],int n)
{
    sort(p,p+n);
    int m=0;
    for(int i=;i<n;i++)
    {
        while(m>1 && Cross(ch[m-1]-ch[m-2],p[i]-ch[m-2])<=0) m--;
        ch[m++]=p[i];
    }
    int k=m;
    for(int i=n-2;i>=0;i--)
    {
        while(m>k) && Cross(ch[m-1]-ch[m-2],p[i]-ch[m-2]<=0) m--;
        ch[m++]=p[i];
    }
    if(n>1) m--;
    return m;
}

3.2旋转卡壳

不要问我怎么读，我也不会读..

3.2.1旋转卡壳

旋转卡壳是用来求平面内，最远点的距离的(就是距离最远的两个点的距离)。
有一个性质非常显然：最远点对一定在凸包上，这种需要我们先求出凸包来.
考虑暴力枚举: $O(n^2)$
然后我们来看旋转卡壳:

//旋转卡壳
double RotatingCaliper(Point ch[],int n)
{
    int q=1; double ans=0;
    ch[n]=ch[0];
    for(int p=0;p<n;p++)
    {
        while(fabs(Cross(ch[q+1]-ch[p+1],ch[p]-ch[p+1])) > fabs(Cross(ch[q]-ch[p+1],ch[p]-ch[p+1]))) q=(q+1)%n;
        ans=max(ans,max(Length2(ch[p]-ch[q]),Length2(ch[p+1]-ch[q+1])));
        ans=max(ans,max(Length2(ch[p]-ch[q+1]),Length2(ch[p+1]-ch[q])));
    }
    return ans;
}

3.3半平面交

半平面交
极角：从x正半轴旋转到当前直线的角度函数cmath atan2()
极角排序：
用双端队列维护
然后每一个考虑加入一条直线，
如果新加入的直线的与倒数第二条直线的交点在倒数第一条直线与倒数第二条直线的交点在新加入直线的右边那么把倒数第一条直线退队.还要和队首判断一下，假如队首第一条直线的交点和队首第二条直线的交点在新加入直线的右边退队..
然后拿队首第一条直线和倒数第二条直线判断一下就可以.

//半平面交
struct Line()
{
    Point P;
    Vector v;
    double ang;
    Line() {}
    Line(Point P,Vector v) : P(P),v(v){ ang=atan2(v.y,v.x); }
    bool operator < (const Line &L) const
    {
        return ang<L.ang;
    }
} 
bool OnLeft(Line L,Point p)
{
    return Cross(L.v,p-L.P)>0;
}
int HalfPlane(Line L[],int n,Point poly[])
{
    sort(L,L+n);
    int first,last;
    Point p[MAXN];
    Line q[MAXN];
    q[first=last=0]=L[0];
    for(int i=1;i<n;i++)
    {
        while(first<last && !OnLeft(L[i],p[last-1]))) last--;
        while(first<last && !OnLeft(L[i],p[first])) first++;
        q[++last]=L[i];
        if(!Dcmp(Cross(q[last].v,q[last-1].v)))
        {
            last--;
            if(OnLeft(q[last],L[i].P)) 
                q[last]=L[i];
        }
        if(first<last) p[last-1]=GetLineIntersection(q[last-1],q[last]);
    }
    if(last-first<=1) return 0;
    p[last]=GetLineIntersection(q[last],q[first]);
    int m=0;
    for(int i=first;i<=last;i++) poly[m++]=p[i];
    return m;
    ···; 
}

3.4最小圆覆盖

这里写图片描述

3.5最近点对

就是求平面中最近的两个点，怎么求呢可以 $n^2$ 暴力，还可以用分治的方法来求，怎么分治，每一次分成左边和右边，最近点对只有三种情况一种是在左边，一种是在右边，一种是在中间，所以我们先求出左边和右边的最近点对，跨过中间的最近点对只能在左右的距离，肯定小于左右最近点对距离然后把在范围的点按y排序，求出距离来就可以了. $O(nlogn)$

double Closest_Pair(int left,int right)
{
    double dis=INF;
    if(left==right)
        return dis;
    if(left+1==right)
        return dist(p[left],p[right]);
    int mid=(left+right)>>1;
    double dis1=Closest_Pair(left,mid);
    double dis2=Closest_Pair(mid,right);
    dis=min(dis1,dis2);
    int num=0;
    for(int i=left;i<=right;i++)
    {
        if(fabs(p[mid].x-p[i].x<=dis))
            now[++num]=p[i];    
    }
    sort(now+1,now+num+1,Y);
    for(int i=1;i<=num;i++)
    {
        for(int j=i+1;j<=num && now[j].y-now[i].y<dis;j++)
        {
            double dis3=Dis(now[i],now[j]);
            dis=min(dis,dis3);
        }   
    }   
    return dis;
}

Tmotfl

关注

2
点赞
踩
1

收藏

觉得还不错? 一键收藏
0
评论
二维计算几何学习笔记

1.基础知识1.1关于浮点误差一般的话要设一个EPS=1e−8EPS=1e−8EPS=1e^{-8}//精度控制 bool Dcmp(double x){ if(fabs(x)&lt;Eps) return 0; else return 1;}1.2向量基础运算1.2.1点积double Dot(Vector A,Vector B)...
复制链接

扫一扫