题解 #9919 Concave Hull

· · 题解

Concave Hull

题意是给定一个大小为 n 的点集，求这个点集的凹包的总面积。

凹包定义为 有且仅有一个内角大于 \pi 且能包含点集中所有点的简单多边形。

先求凸包。枚举内角大于 \pi 对应的那个点，再枚举凸包凹进去的那条边，可以扣掉一个三角形。按照凹进去的点做极角排序后发现在三角形范围内，任意两个极角序上相邻的点都可以唯一地确定凹包的形态。

再根据大于 \pi 的点的唯一性，我们可以发现 凹包 相对于 凸包扣掉三角形 的面积增量是两个凸包。

这两个凸包的面积是对称的，考虑如何求一边，发现极角序同样有优秀的性质。按照极角序往凸包加点，每个点只会朝一个方向弹栈，所以可以同样用单调栈维护，用增量法求前缀面积。

另一个方向同理，预处理完之后就可以在枚举极角序相邻点时 \mathcal{O}(1) 查询对应的面积了。

时间复杂度 \mathcal{O}\big(n^2\log n\big)。

#pragma GCC optimize(2)
#pragma GCC optimize(3)
#pragma GCC optimize("Ofast")
#pragma GCC optimize("unroll-loops")

#include<bits/stdc++.h>
#define fi first
#define se second
#define mp make_pair
#define pb push_back
typedef long long ll;
typedef long double ld;
using namespace std;
const int N=2010,M=20010,lpw=1e9+7;

struct node{
    ll x,y,id;
    node(ll xx=0,ll yy=0,ll idd=0){x=xx,y=yy,id=idd;}
    void in(){cin>>x>>y;}
    void out(){cout<<'('<<x<<','<<y<<')'<<'\n';}
    bool operator <(const node &a)const{
        return x!=a.x?x<a.x:y<a.y;
    }
    bool operator ==(const node &a)const{
        return x==a.x&&y==a.y;
    }
};

node operator +(const node &a,const node &b){return node(a.x+b.x,a.y+b.y);}
node operator -(const node &a,const node &b){return node(a.x-b.x,a.y-b.y);}
ll operator *(const node &a,const node &b){return a.x*b.x+a.y*b.y;}
ll operator ^(const node &a,const node &b){return a.x*b.y-a.y*b.x;}

node O(0,0),P(-1,0);
bool cmp(const node &a,const node &b){
    if(!((a-O)^(b-O))){
        if(a*b>0)return a*a<b*b;
        else if(!(a^P))return a.x<b.x;
        else return (a^P)<(b^P);
    }
    if(!((P-O)^(a-O))&&(P.x-O.x)*(a.x-O.x)>0)return 1;
    if(!((P-O)^(b-O))&&(P.x-O.x)*(b.x-O.x)>0)return 0;
    if((((P-O)^(a-O))>0)!=(((P-O)^(b-O))>0))return ((P-O)^(a-O))>((P-O)^(b-O));
    return ((a-O)^(b-O))>0;
}

node p[N],b[N];

int n,cnt,v[N],s[N],t;

ll ans,S,pre[N],suf[N];

vector<node> h;

int prv(int i){return i==0?n-2:i-1;}
int nxt(int i){return i==n-2?0:i+1;}

int main(){
    ios::sync_with_stdio(false);cin.tie(0);cout.tie(0);
    cin>>n;
    for(int i=1;i<=n;i++)
        p[i].in();
    sort(p+1,p+n+1);
    s[t=1]=1;
    for(int i=2;i<=n;i++){
        while(t>1&&((p[i]-p[s[t]])^(p[i]-p[s[t-1]]))>0)t--;
        s[++t]=i;
    }
    for(int i=1;i<=t;i++){
        v[s[i]]=1;
        b[++cnt]=p[s[i]];
    }
    s[t=1]=1;
    for(int i=2;i<=n;i++){
        while(t>1&&((p[i]-p[s[t]])^(p[i]-p[s[t-1]]))<0)t--;
        s[++t]=i;
    }
    for(int i=t-1;i>1;i--){
        v[s[i]]=1;
        b[++cnt]=p[s[i]];
    }
    for(int i=1;i<=cnt;i++)
        S+=(b[i]^b[i==cnt?1:i+1]);
    S=abs(S)%lpw;
    for(int i=1;i<=n;i++){
        if(v[i])continue;
        O=node(p[i].x,p[i].y);
        P=node(p[i].x-1,p[i].y);
        h.clear();
        for(int j=1;j<=n;j++)
            if(j!=i)h.pb(node(p[j].x,p[j].y,j));
        sort(h.begin(),h.end(),cmp);
        for(int l=0;l<=n-2;l++){
            if(!v[h[l].id])continue;
            int r=nxt(l);
            for(;;r=nxt(r))
                if(v[h[r].id])break;
            ll T=abs((p[h[l].id]-p[i])^(p[h[r].id]-p[i]))%lpw;
            if(r==nxt(l)){
                ans=(ans+S-T+lpw)%lpw;
                continue;
            }
            pre[l]=suf[r]=0;
            s[t=1]=l;
            for(int j=nxt(l);j!=r;j=nxt(j)){
                while(t>1&&((p[h[j].id]-p[h[s[t]].id])^(p[h[j].id]-p[h[s[t-1]].id]))>0)t--;
                pre[j]=pre[s[t]]+abs((p[h[j].id]-p[i])^(p[h[s[t]].id]-p[i]));
                pre[j]%=lpw;
                s[++t]=j;
            }
            s[t=1]=r;
            for(int j=prv(r);j!=l;j=prv(j)){
                while(t>1&&((p[h[j].id]-p[h[s[t]].id])^(p[h[j].id]-p[h[s[t-1]].id]))<0)t--;
                suf[j]=suf[s[t]]+abs((p[h[j].id]-p[i])^(p[h[s[t]].id]-p[i]));
                suf[j]%=lpw;
                s[++t]=j;
            }
            for(int j=l;j!=r;j=nxt(j))
                ans=(ans+S-T+pre[j]+suf[nxt(j)]+lpw)%lpw;
        }
    }
    cout<<ans<<'\n';
    return 0;
}