[计算几何] 闵可夫斯基和
(凸包的)闵可夫斯基和:P4557 [JSOI2018]战争
定义
官方定义:两个点集 \(A,B\) 的闵可夫斯基和为: \(A+B={a+b\ |\ x\in A,y\in B}\)
通俗理解:从原点向图形 \(A\) 内部的每一个点做向量,将图形 \(B\) 向每个向量移动,所有的最终位置便是闵可夫斯基和。
(图片是 \(copy\) 的)
一些结论
-
闵可夫斯基和一定是凸集
-
从形态上看,一条出现在任意凸包上的边,一定会出现在闵可夫斯基和上,或者说,每一条凸包的边都对应一条闵可夫斯基和的边
求法
无非就是合并两个凸包的过程 我却调了一个下午
现在有两个凸包,显然他们的点是经过极值排序的。
只需要重新把两份边集按照极角排序,是一个归并排序的过程(两边已经单调了)
所以我把所有边代表的向量用两个数组 \(s1,s2\) 存起来,用于方便比较极角,更重要的是把边连接起来(向量加法)。
我准备把代码拆分来纪念 \(Debug\) 的一个下午,因为题解这份代码巨方便 就是容易写错
代码
首先是判断一个点是否在凸包内:
这样写起来比较方便,前提是 \(A[1]\) 是坐标原点
可以理解为:
-
先判与边界的关系
-
找到极角 \(<\) 它的极角的最大极角的点
-
利用叉积进行比较
inline LL in(node a){
if(a*A[1]>0 || A[tot]*a>0)return 0;
LL pos=lower_bound(A+1,A+tot+1,a,cmp2)-A-1;
return (a-A[pos])*(A[pos%tot+1]-A[pos])<=0;//=0时在边上
}
结构体的定义:
struct node{
LL x,y;
node operator - (node A) {return (node){x-A.x,y-A.y};}
node operator + (node A) {return (node){x+A.x,y+A.y};}
LL operator * ( node A) const {return x*A.y-y*A.x;}
LL len() const {return x*x+y*y;}
};
一种扫描法求凸包的写法:
void Tubao(node *p,LL &n){//返回一个凸包,写法很巧妙
sort(p+1,p+1+n,cmp1);
p0=p[1];
stk[top=1]=1;
for(int i=1;i<=n;i++)p[i]=p[i]-p0;
sort(p+2,p+1+n,cmp2);
for(int i=2;i<=n;i++){
while(top>1 && (p[i]-p[stk[top-1]])*(p[stk[top]]-p[stk[top-1]])>=0)top--;
stk[++top]=i;
}
for(int i=1;i<=top;i++)p[i]=p[stk[i]]+p0;//这个算法的麻烦的地方就是需要还原点的坐标
// for(int i=1;i<=top;i++)printf("A[%d]: %d %d\n",i,p[i].x,p[i].y);
// puts("");
n=top;p[n+1]=p[1];
}
合并凸包(归并排序+向量):
void Minkowski(){
for(LL i=1;i<n;i++)s1[i]=p1[i+1]-p1[i];
s1[n]=p1[1]-p1[n];
for(LL i=1;i<m;i++)s2[i]=p2[i+1]-p2[i];
s2[m]=p2[1]-p2[m];
A[tot=1]=p1[1]+p2[1];
LL cnt1=1,cnt2=1;
//把两个凸包进行合并(p1,p2)
while(cnt1<=n && cnt2<=m)++tot,A[tot]=A[tot-1]+(s1[cnt1]*s2[cnt2]>=0?s1[cnt1++]:s2[cnt2++]);
while(cnt1<=n)++tot,A[tot]=A[tot-1]+s1[cnt1++];
while(cnt2<=m)++tot,A[tot]=A[tot-1]+s2[cnt2++];
}
因此总体代码长这样:
#include <iostream>
#include <cstdio>
#include <cstring>
#include <algorithm>
#include <cmath>
using namespace std;
const int maxn = 1e5 + 100;
#define LL long long
struct node{
LL x,y;
node operator - (node A) {return (node){x-A.x,y-A.y};}
node operator + (node A) {return (node){x+A.x,y+A.y};}
LL operator * ( node A) const {return x*A.y-y*A.x;}
LL len() const {return x*x+y*y;}
}A[maxn],p1[maxn],p2[maxn],s1[maxn],s2[maxn],p0;
LL stk[maxn],n,m,q,top,tot;
LL cmp1(node A,node B) {return A.y==B.y?A.x<B.x:A.y<B.y;}
LL cmp2(node A,node B) {return A*B==0?A.len()<B.len():A*B>0;}
void Tubao(node *p,LL &n){//返回一个凸包,写法很巧妙
sort(p+1,p+1+n,cmp1);
p0=p[1];
stk[top=1]=1;
for(int i=1;i<=n;i++)p[i]=p[i]-p0;
sort(p+2,p+1+n,cmp2);
for(int i=2;i<=n;i++){
while(top>1 && (p[i]-p[stk[top-1]])*(p[stk[top]]-p[stk[top-1]])>=0)top--;
stk[++top]=i;
}
for(int i=1;i<=top;i++)p[i]=p[stk[i]]+p0;//这个算法的麻烦的地方就是需要还原点的坐标
// for(int i=1;i<=top;i++)printf("A[%d]: %d %d\n",i,p[i].x,p[i].y);
// puts("");
n=top;p[n+1]=p[1];
}
void Minkowski(){
for(LL i=1;i<n;i++)s1[i]=p1[i+1]-p1[i];
s1[n]=p1[1]-p1[n];
for(LL i=1;i<m;i++)s2[i]=p2[i+1]-p2[i];
s2[m]=p2[1]-p2[m];
A[tot=1]=p1[1]+p2[1];
LL cnt1=1,cnt2=1;
//把两个凸包进行合并(p1,p2)
while(cnt1<=n && cnt2<=m)++tot,A[tot]=A[tot-1]+(s1[cnt1]*s2[cnt2]>=0?s1[cnt1++]:s2[cnt2++]);
while(cnt1<=n)++tot,A[tot]=A[tot-1]+s1[cnt1++];
while(cnt2<=m)++tot,A[tot]=A[tot-1]+s2[cnt2++];
}
inline LL in(node a){
if(a*A[1]>0 || A[tot]*a>0)return 0;
LL pos=lower_bound(A+1,A+tot+1,a,cmp2)-A-1;
return (a-A[pos])*(A[pos%tot+1]-A[pos])<=0;//=0时在边上
}
int main(){
//freopen("1.in","r",stdin);
scanf("%lld%lld%lld",&n,&m,&q);
for(int i=1;i<=n;i++){
scanf("%lld%lld",&p1[i].x,&p1[i].y);
}
Tubao(p1,n);
for(int i=1;i<=m;i++){
scanf("%lld%lld",&p2[i].x,&p2[i].y);
p2[i].x=-p2[i].x;p2[i].y=-p2[i].y;
}
Tubao(p2,m);
Minkowski();
Tubao(A,tot);
p0=A[1];
for(int i=tot;i>=1;i--)A[i]=A[i]-p0;
while(q--){
scanf("%lld%lld",&A[0].x,&A[0].y);
printf("%lld\n",in(A[0]-p0));
}
return 0;
}