不知不觉简单的计算几何就快讲完了,这次我再讲讲旋转卡壳算法(看得我脑袋是有些卡壳QwQ)。这个算法和凸包有着很大的关系,所以建议没看凸包的同学去计算几何分类看一下凸包的求解方法。
先抛出一个问题
在二维平面上给定n个点,求距离最远的两个点之间的距离是多少?
看到这道题,大部分同学第一反应肯定是暴力——暴力枚举所有点对,然后求出最远的那一对相距的距离。这种方法当然没错,可是复杂度已经达到了O(n^2),当n是50000,100000,甚至1000000时,暴力就显得苍白无力。
这时我们进一步思考,这个最远点对一定会出现在这些点集的凸包上,所以在这里我们在凸包上,使用复杂度低至O(n)的旋转卡壳算法进行求解。
什么是旋转卡壳算法
旋转卡(qiǎ)壳算法:是解决一些与凸包有关问题的有效算法,计算时就像一对平行卡壳卡住凸包旋转而得名。
【图片来自:http://cgm.cs.mcgill.ca/~orm/rotcal.html 】
在凸包上,被一对平行卡壳正好卡住的对应点对称为对踵点;可以证明对踵点的个数不超过3n/2个 也就是说对踵点的个数是O(n)的,对踵点的个数也是我们下面解决问题时间复杂度的保证。
我们发现,平行卡壳分为两种情况:
(1)卡住两个点:
(2)卡住一条边和一个点:
由于在处理过程中,卡壳卡住两个点的情况不好处理,所以一般我们都是去关注第二种情况。在这种情况下我们发现,一个对踵点到对应边之间的距离比其他点要大。
也就是说,一个对踵点和对应边所形成的三角形的面积是最大的,所以我们可以使用叉积的第二种用法——求三角形面积来求解对踵点。
经过比对我们发现,对踵点对距离中的最大值,就是在这个二维平面内两点距离的最大值,总的时间复杂度为O(n)。
怎样实现旋转卡壳算法
首先,由于旋转卡壳算法是建立在凸包上的,我们要先对所有点求出一个凸包,放在数组(设为sta数组)内,由于之前已经细细讲过凸包了,这里就不做细讲,不懂凸包的同学请先移步二维凸包。
现在我们需要考虑,如何得到距离每条对应边的的最远点呢?稍加分析,我们可以发现,我们可以把点到直线的距离化解为三角形的面积,再运用叉积进行比较。同时,凸包上的点依次与对应边产生的距离成单峰函数,面积上升到最高点后,又会下降。(具体证明可以从凸包定义入手 用反证法解决)
利用单峰函数这一性质,我们可以根据凸包上点的顺序,枚举对踵点,直到下一个点的距离小于当前点就可以停止了,而且随着对应边的旋转,最远点也只会顺着这个方向旋转,我们可以从上一次的对踵点开始继续寻找这一次的。由于内层while循环的执行次数取决于j增加次数 j最多增加O(n)次,所以求出所有对踵点的时间复杂度为O(n)。
代码如下:
sta[top+1]=sta[1];j=2;ans=0;
for(int i=1;i<=top;i++)
{
while(fabs(multi(sta[i],sta[i+1],sta[j]))<fabs(multi(sta[i],sta[i+1],sta[j+1])))
{
j++;
if(j>top)j=1;
}
//js=dis(sta[i],sta[j]); 错误版本
js=max(dis(sta[i],sta[j]), dis(sta[i+1],sta[j]));
if(js>ans)ans=js;
}
简单的模板题
【caioj 1226】最远点对的距离
题目描述
【题意】
给出n个点的坐标,求最远两点间的距离。
【输入格式】
第一行一个整数n(2 ≤ n ≤ 50000)。
下来n行,每行两个实数x和y表示点坐标(x,y<=|10^6|)。
【输出格式】
一行一个实数,表示最远两点间的距离的平方(保留2位小数)。
【样例输入】
5
0 0
0 5
5 0
5 5
2 0
【样例输出】
50.00
这就是一道经典的模板题。只要凸包再跑一遍旋转卡壳算法就可以完美求解了。
附上代码(同时也是模板),并特别感谢fatzard和fiora_gwx对代码的指正:
#include<cmath>
#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<algorithm>
using namespace std;
struct node
{
double x,y;
}p[510000],sta[510000];
int i,j,k,m,n,o,top;
double ans,jl,js;
double multi(node p1,node p2, node p0)
{
double x1=p1.x-p0.x;
double y1=p1.y-p0.y;
double x2=p2.x-p0.x;
double y2=p2.y-p0.y;
return(x1*y2-x2*y1);
}//计算叉积
double dis(node p1,node p2)
{
return (double((p1.x-p2.x)*(p1.x-p2.x))+double((p1.y-p2.y)*(p1.y-p2.y)));
}//计算任意两点间距离
bool cmp(node p1,node p2)
{
double tt=multi(p1,p2,p[1]);
if(tt<0)return(0);
if(tt==0 && dis(p1,p[1])>dis(p2,p[1]))return(0);
return(1);
}
void graham()
{
sort(p+2,p+n+1,cmp);
for(int i=1;i<=2;i++)sta[i]=p[i];
top=2;
n++;p[n]=p[1];
for(int i=3;i<=n;i++)
{
while(top>1 && multi(sta[top],p[i],sta[top-1])<=0)top--;
top++;
sta[top]=p[i];
}
}
int main()
{
scanf("%d",&n);
for(int i=1;i<=n;i++)
{
scanf("%lf%lf",&p[i].x,&p[i].y);
if(p[i].y<p[1].y ||(p[i].y==p[1].y && p[i].x<p[1].x))
{
node t=p[1];p[1]=p[i];p[i]=t;
}//找最下方的点
}
top=0;//记录凸包上点的个数
graham();
top--;//sta中存储了top+1个点,第top+1个和第一个相同
j=2;ans=0;
for(int i=1;i<top;i++)
{
while(fabs(multi(sta[i],sta[i+1],sta[j]))<fabs(multi(sta[i],sta[i+1],sta[j+1])))
{
j++;
if(j>top)j=1;
}
js=max(dis(sta[i],sta[j]), dis(sta[i+1],sta[j]));
if(js>ans)ans=js;
}//旋转卡壳
printf("%0.2lf",ans);
return 0;
}
结语
通过这篇博客,我相信神犇的你一定已经熟练的掌握旋转卡壳算法了,到这里计算几何也差不多接近尾声了,莫名有些悲伤:cry:。最后感谢你阅读这篇博客!
Orz
是不是哪里有问题
这个数据很明显答案是36,但是您代码是29
4
0.0 0.0
3.0 3.0
-2.0 5.0
-3.0 3.0
您好,非常感谢您对代码的指正。数据结果不对是因为存储凸包上点的数量的变量top多了一导致在轮回一圈时造成错误,现在已经修正,非常感谢!
没明白为什么都去算三角形面积,我得再理解理解,我看了一个大佬的代码,没有算面积,直接算的是距离,这样是对的
double maxDist(vec &t){
vec u=t;//ConvexHull(t);
int N=u.size(),K=0;
double res=0;
for(int i=0;i<N;i++){
while(abs(u[i]-u[K])<abs(u[i]-u[(K+1)%N]))K=(K+1)%N;
res=max(res,abs(u[i]-u[K]));
}
return res;
}
您好,用三角形的面积主要是因为当我们的线卡住一条边的时候,另一边卡住的点一定是离这条边最远的点,即和底边形成的三角形高最长的边,即面积最大的点。当然其他大佬的写法也应该是对的,不过可能和我切入的角度有些区别。
楼主,你这个对于数据
8
5 15
15 10
0 0
10 0
-5 10
0 10
5 5
-1 6
跑出来是325,
但是应该为400
感谢您的指正,代码的问题出现在对于每个三角形,底边的两个点都应该去计算贡献,现在已经改成js=max(dis(sta[i],sta[j]), dis(sta[i+1],sta[j]))!感谢!
讲的太好了
谢谢您!
旋转卡壳算法可以求n个点的最小距离没吗
n个点的最小距离可以看这篇文章,里面有平面最近点对哦 https://www.jvruo.com/archives/162/
看到博主的blog里有一篇是Ziyin/Zayin的题,您是SYSU的?
嗯嗯是的,是大一的新生
希望两年后我能做您的学弟吧~(虽然我现在还在垂死挣扎OI但估计是马上退役了),以我的水平肯定考不上上交那种学校啦,但愿能去双鸭山!(话说Zayin刚给我们上了一个月的课^ ^)
博主好棒QwQ ,现在还在打ACM吗
谢谢您,是的QwQ
请问为什么单对单的对踵点不好处理呢(雾)
在原算法中,我们是通过枚举边和这个边所在三角形的一个单调性(面积先上升后下降)来移动对踵点;而没有边的辅助后,我们的算法就是移动一个点找另一个点,这样就没有以上我们提到的单调性,就可能要枚举完剩下的所有点才能找到这个点对应的对踵点。所以我们说直接找单对单的对踵点是不方便的。
orz
ORZ‘
Orz
太强了,我何时才能这么优秀
您校队C位太强了吧
堇了 才高三
为什么大家都喜欢摸鱼了
不说了 该摸了
哈哈一起摸鱼
%%%%
写的非常好w ୧(๑•̀⌄•́๑)૭
然后有个小问题
这份代码跑这样的数据似乎会出锅
0 0
x 0
y 0
因为为最后栈中剩下的只有两个(0,0)而已
要怎么解决呢∠( ᐛ 」∠)_
确实,旋转卡壳貌似一般都是用于解决二维平面(存在凸包)问题,对于一维的问题(貌似没有一维凸包的说QwQ),由于本身就没有凸包,所以会无法卡壳出现问题,此时应进行特判,最后感谢你喜欢这篇BLOG!
okok!知道啦
感谢答疑!٩(ˊᗜˋ*)و
(话说博主回的好快呀OωO
哈哈感谢你支持蒟蒻的BLOG,博主已经高三啦,所以只有在摸鱼的时候才能上上博客啦QwQ