计算几何专题
谨以此篇,铭记NOI2022退役的lnzwz。
计算几何
向量
向量,即带有方向的量。一个$n$维向量可以表示为一个$n$维的点。比如$\vec{v} = (a,b)$,一般情况下两维向量比较常用.
运算
向量的模长:$|(a_1,a_2,\dots ,a_n)|=\sqrt{\sum a_i^2}$
向量加减:直接把每一维的值加减即可。$(a_1,a_2,\dots ,a_n) \pm (b_1,b_2,\dots ,b_n) = (a_1 \pm b_1,a_2 \pm b_2,\dots ,a_n \pm b_n)$
向量乘常量: 每一维都乘 $\lambda (a_1,a_2,\dots ,a_n)=(\lambda a_1,\lambda a_2,\dots ,\lambda a_n)$
向量点积: $(a_1,a_2,\dots ,a_n) \cdot (b_1,b_2,\dots ,b_n)=\sum a_i * b_i = |\vec{a}||\vec{b}|cos (\vec{a},\vec{b}) $ 相当于$\vec{a}$在$\vec{b}$上的投影长度乘$|\vec{b}|$。当 $\vec{a} \bot \vec{b}$ 时有 $\vec{a} \cdot \vec{b} =0$ 。
向量叉积:$(x_1,y_1) \times (x_2,y_2) = x_1 y_2-x_2 y_1=|\vec{a}||\vec{b}|sin (\vec{a},\vec{b})$。比较重要,当$\vec{a} \times \vec{b} $ > $0$时$ \vec{a}$在$\vec{b}$顺时针方向,当$\vec{a} \times \vec{b} $ < $0$时$ \vec{a}$在$\vec{b}$逆时针方向。
叉积判断线段相交
若两条线段相交,当且仅当对于一条线段,另一条线段的端点在其所在直线上或者分别在其所在直线两侧。
对于两条线段,分别任取一个端点利用叉积判断另一条线段两端点和他的连线 与其所在线段的向量的叉积的乘积是否小于等于0.
如图,则有$(\vec{AC} \times \vec{AB})*(\vec{AD} \times \vec{AB}) \leq 0 \bigwedge (\vec{CA} \times \vec{CD})*(\vec{CB} \times \vec{CD}) \leq 0 \Longleftrightarrow AB与CD相交$
叉积求三角形面积
$S_{\bigtriangleup ABC}=\frac{1}{2} BC *AD=\frac{1}{2}|\vec{CA}||\vec{CB}|sin C=\frac{1}{2}|\vec{CA} \times \vec{CB}| $
叉积求凸多边形面积
在平面内任取一点O,$$S=\frac{1}{2}|\vec{OA_{N-1}} \times \vec{OA_{0}}+\sum_{i=1}^{N-1} \vec{OA_{i-1}} \times \vec{OA_{i}}|$$
使用向量来表示点和直线
对于一个点,可以用原点到这个点的向量来表示它。
对于一条直线,可以用直线上的一个点和以这个点为起点且与直线平行的向量来表示,即点$A$可以表示为$(\vec{v_1} , \vec{v_2})$。
这样表示直线的好处是可以表示斜率无穷大的直线。
利用向量叉积求直线交点
两条直线$(\vec{p_1} , \vec{v_1})$,$(\vec{p_2} , \vec{v_2})$
令$$ t = \frac{(\vec{p_2}-\vec{p_1}) \times \vec{v_2}}{\vec{v_1} \times \vec{v_2}}$$
则有交点$\vec{u} = \vec{p_1} + t \vec{v_1}$
凸包
二维凸包可以类比为:在二维平面上有若干个点,用橡皮筋紧紧围住他们,橡皮筋形成的图形。
凸包是个凸多边形。
凸包求法
先找出y坐标最小的点O(一定在凸包上),将其他点按照O与其的向量的极角排序,利用单调栈(斜率单调)维护即可。
using namespace std;
#include<bits/stdc++.h>
#define db double
namespace zhu{
int N;
const int maxn=100010;
struct vct{
db x,y;
friend vct operator -(vct a,vct b){return (vct){a.x-b.x,a.y-b.y};}
friend db operator *(vct a,vct b){return a.x*b.y-a.y*b.x;}
}k[maxn];
db dis(vct &a,vct &b){
return sqrt((a.x-b.x)*(a.x-b.x)+(a.y-b.y)*(a.y-b.y));
}
int sta[maxn],top=0;
bool cmp(vct a,vct b){
db tmp=(b-k[1])*(a-k[1]);
if(tmp>0) return 1;
else if(tmp<0) return 0;
else if(dis(a,k[1])<dis(b,k[1])) return 1;
else return 0;
}
int main(){
cin>>N;
int h=1;
for(int i=1;i<=N;i++){
scanf("%lf%lf",&k[i].x,&k[i].y);
if(k[i].y<k[h].y||(k[i].y==k[h].y&&k[i].x<=k[h].x)) h=i;
}
swap(k[1],k[h]);
sta[++top]=1;
sort(k+2,k+N+1,cmp);
for(int i=2;i<=N;i++){
while(top>2&&((k[i]-k[sta[top-1]])*(k[sta[top]]-k[sta[top-1]]))<=0){
top--;
}
sta[++top]=i;
}
db ans=k[sta[1]]+k[sta[top]];
for(int i=1;i<top;i++){
ans+=k[sta[i]]+k[sta[i+1]];
}
printf("%.2f",ans);
return 0;
}
};
int main(){
return zhu::main();
}
旋转卡壳
旋转卡壳是用来求凸包直径的,容易证明,凸包直径一定是某两个顶点连线长度。
这可以用平面最远点对来求,即KD-Tree($O(n\sqrt{n})$).
旋转卡壳时间复杂度是线性的。(贺一个动图)
旋转卡壳直观表示。
如果过凸包上的两个点可以画一对平行直线,使凸包上的所有点都夹在两条平行线之间或落在平行线上,那么这一对点叫做一对 对踵点。
所以我们可以枚举凸包上的点,找出其对应的对踵点求距离即可。
然后我们发现随着枚举,对踵点会旋转,只需要判断一下接下来的对踵点是哪个即可。
using namespace std;
#include<bits/stdc++.h>
#define db double
namespace zhu{
int N;
const int maxn=50010;
struct vct{
int x,y;
friend vct operator -(vct a,vct b){return (vct){a.x-b.x,a.y-b.y};}
friend int operator *(vct a,vct b){return a.x*b.y-a.y*b.x;}
}k[maxn];
int dis(vct &a,vct &b){
return (a.x-b.x)*(a.x-b.x)+(a.y-b.y)*(a.y-b.y);
}
int sta[maxn],top=0;
bool cmp(vct &a,vct &b){
db tmp=(b-k[1])*(a-k[1]);
if(tmp>0) return 1;
else if(tmp<0) return 0;
else if(dis(a,k[1])<dis(b,k[1])) return 1;
else return 0;
}
int main(){
cin>>N;
int h=1;
for(int i=1;i<=N;i++){
scanf("%d%d",&k[i].x,&k[i].y);
if(k[i].y<k[h].y||(k[i].y==k[h].y&&k[i].x<=k[h].x)) h=i;
}
swap(k[1],k[h]);
sta[++top]=1;
sort(k+2,k+N+1,cmp);
for(int i=2;i<=N;i++){
while(top>=2&&((k[i]-k[sta[top-1]])*(k[sta[top]]-k[sta[top-1]]))<=0){
top--;
}
sta[++top]=i;
}
if(top==1) cout<<0<<endl,exit(0);
if(top==2) cout<<dis(k[sta[1]],k[sta[2]])<<endl,exit(0);
sta[top+1]=1;
int ans=0;
for(int i=1,j=2;i<=top;i++){
while(abs((k[sta[i+1]]-k[sta[i]])*(k[sta[j]]-k[sta[i]]))<abs((k[sta[i+1]]-k[sta[i]])*(k[sta[j+1]]-k[sta[i]]))){
j++;
if(j==top+1) j=1;
}
vct &now1=k[sta[i]],&now2=k[sta[i+1]];
vct &fac=k[sta[j]];
ans=max(ans,max(dis(now1,fac),dis(now2,fac)));
}
cout<<ans<<endl;
return 0;
}
};
int main(){
return zhu::main();
}
半平面交
文哲一路顺风!
一条直线将一个平面分成了两个 半平面。
一个凸$N$边形可视为$N$个半平面相交。
我们求出所有的半平面交,最后会得到一个凸包,凸包的边在某些直线上。之后可以对其进行求面积等操作。
首先将这些直线按照极角排序,维护一个队列,利用叉积判断新直线是否更优。
核心部分为:求两条直线交点,找凸包,求面积。
using namespace std;
#include<bits/stdc++.h>
typedef double db;
namespace zhu{
int N;
const int maxn=510;const db eps=1e-6;
struct vct{
db x,y;
db jijiao(){
return atan2(y,x);
}
friend vct operator -(vct a,vct b){return (vct){a.x-b.x,a.y-b.y};}
friend vct operator +(vct a,vct b){return (vct){a.x+b.x,a.y+b.y};}
friend db operator *(vct a,vct b){return a.x*b.y-a.y*b.x;}
friend vct operator *(db b,vct a){return (vct){a.x*b,a.y*b};}
}__tmp[1000];
struct line{
vct poi,v;
line(){}
line(vct a,vct b){poi=a,v=b;}
friend bool operator == (const line &a,const line &b){
return abs(a.v*b.v)<eps&&abs(a.v*(a.poi-b.poi))<eps;
}
}l[maxn];
int cnt=0;
bool cmp(line &a,line &b){
return a.v.jijiao()<b.v.jijiao();
}
vct jd(const line &a,const line &b){
db t=(b.poi-a.poi)*b.v / (a.v * b.v);
return a.poi + t * a.v;
}
line tb[maxn];
vct jp[maxn];
int L=1,R=0;
void work(){
tb[++R]=l[1];
for(int i=2;i<=cnt;i++){
while( L<R && (l[i].v * (jp[R-1] - jd(l[i],tb[R-1]))) <eps) R--;
while( L<R && (l[i].v * (jp[L] - jd(l[i],tb[L] ))) <eps) L++;
tb[++R]=l[i];
if( L<R && abs(tb[R].v * tb[R-1].v) <eps){
R--;
if( (tb[R+1].poi-tb[R].poi) * tb[R].v < eps) tb[R]=l[i];
}
if(L<R) jp[R-1]=jd(tb[R],tb[R-1]);
}
while( L<R && (tb[L].v * (jp[R-1] - jd(tb[L],tb[R-1]))) <eps ) R--;
if(L!=R) jp[R]=jd(tb[R],tb[L]);
}
db solve(){
db ans=jp[R] * jp[L];
for(int i=L;i<R;i++) ans+=jp[i] * jp[i+1];
return 0.5*abs(ans);
}
int main(){
cin>>N;
for(int i=1;i<=N;i++){
int k;
scanf("%d",&k);
for(int i=1;i<=k;i++) scanf("%lf%lf",&__tmp[i].x,&__tmp[i].y);
for(int i=1;i<k;i++){
l[++cnt]=line(__tmp[i],__tmp[i+1]-__tmp[i]);
}
l[++cnt]=line(__tmp[k],__tmp[1]-__tmp[k]);
}
sort(l+1,l+cnt+1,cmp);
work();
cout<<fixed<<setprecision(3)<<solve()<<endl;
return 0;
}
};
int main(){
return zhu::main();
}