[模板] 计算几何2: 自适应Simpson/凸包/半平面交/旋转卡壳/闵可夫斯基和

Posted 2022-02-22 ubospica

//to update

一些基本的定义在这里: [模板] 计算几何1(基础): 点/向量/线/圆/多边形/其他运算

自适应Simpson

凸包

Andrew 算法, 即分别求上, 下凸包. 时间复杂度 \\(O(n \\log n)\\).

struct tvecdb x,y;;

il int dcmp(db a)return fabs(a)<=eps?0:(a>0?1:-1);
il db p2(db a)return a*a;
il db gougu1(db a,db b)return sqrt(p2(a)+p2(b));
il tvec operator+(tvec a,tvec b)return (tvec)a.x+b.x,a.y+b.y;
il tvec operator-(tvec a,tvec b)return (tvec)a.x-b.x,a.y-b.y;
il tvec operator*(tvec a,db b)return (tvec)a.x*b,a.y*b;
il tvec operator*(db a,tvec b)return b*a;
il db operator*(tvec a,tvec b)return a.x*b.y-b.x*a.y;
il db operator^(tvec a,tvec b)return a.x*b.x+a.y*b.y;
il db len(tvec a)return gougu1(a.x,a.y);
bool cmp(tvec a,tvec b)int tmp=dcmp(a.x-b.x);return tmp?tmp<0:dcmp(a.y-b.y)<0;

tvec li[nsz],conv[nsz];
int pc=0;
void getconv()
    sort(li+1,li+n+1,cmp);
    rep(i,1,n)
        while(pc>1&&dcmp((conv[pc]-conv[pc-1])*(li[i]-conv[pc]))<=0)--pc;
        conv[++pc]=li[i];
    
    int tmp=pc;
    repdo(i,n-1,1)
        while(pc>tmp&&dcmp((conv[pc]-conv[pc-1])*(li[i]-conv[pc]))<=0)--pc;
        conv[++pc]=li[i];
    
    if(n>1)--pc;

半平面交

增量法, 时间复杂度 \\(O(n \\log n)\\) (排序的复杂度).

需要保证不是开放的半平面. 否则加上四个 \\(\\pm \\infty\\) 的平面即可.

细节较多. 详见代码...

const int psz=550;
const db eps=1e-9;
int n,m;

db dcmp(db v)return fabs(v)<=eps?0:(v>0?1:-1);
db p2(db v)return v*v;

struct tvecdb x,y;;
tvec operator+(tvec a,tvec b)return (tvec)a.x+b.x,a.y+b.y;
tvec operator-(tvec a,tvec b)return (tvec)a.x-b.x,a.y-b.y;
tvec operator*(tvec a,db b)return (tvec)a.x*b,a.y*b;
tvec operator*(db a,tvec b)return b*a;
db operator*(tvec a,tvec b)return a.x*b.y-a.y*b.x;
db operator^(tvec a,tvec b)return a.x*b.x+a.y*b.y;
db len(tvec a)return sqrt(p2(a.x)+p2(a.y));

struct tl
    tvec p,v;
    db d;
    tl()
    tl(tvec a,tvec b):p(a),v(b-a)d=atan2(v.y,v.x);
li[psz];
int pl=0;
bool operator<(tl a,tl b)return a.d<b.d;
bool isleft(tl a,tvec b)return dcmp(a.v*(b-a.p))>0;
il tvec inters(tl a,tl b)db v=(b.v*(a.p-b.p))/(a.v*b.v);return a.p+a.v*v;

tvec poly[psz];
int ppo=0;

tl que[psz]; //queue
tvec qp[psz]; //intersect points
int qh=1,qt=0;
int hplane()//0 fail, 1 success
    sort(li+1,li+pl+1);
    int pl1=1;//suppose that pl>=1
    rep(i,2,pl)
        if(li[i].d>li[pl1].d)li[++pl1]=li[i];
        else if(isleft(li[pl1],li[i].p))li[pl1]=li[i];
    
    pl=pl1;
    qh=1,qt=0;
    rep(i,1,pl)
        while(qh<qt&&!isleft(li[i],qp[qt-1]))--qt;
        while(qh<qt&&!isleft(li[i],qp[qh]))++qh;
        que[++qt]=li[i];
        if(qh<qt)qp[qt-1]=inters(que[qt-1],que[qt]);
    
    while(qh<qt&&!isleft(que[qh],qp[qt-1]))--qt; //**
    ppo=0;
    if(qt-qh<=1)return 0; //no sol
    qp[qt]=inters(que[qh],que[qt]);
    rep(i,qh,qt)poly[++ppo]=qp[i];
    return 1;

旋转卡壳

这是一种拥有 \\(4\\) 个多音字, \\(2^4 = 16\\) 种读音的优秀算法.