[计算几何]POJ2079 求点集中面积最大的三角形

http://acm.pku.edu.cn/JudgeOnline/problem?id=2079

如果枚举3个点再算面积,必然超时。

很显然的最大面积的三角形的三个顶点必然是这个点集的凸包上的点,因此先求出凸包。

如果求出凸包仍然枚举,一样会超时。

这个可以借助求凸包直径类似的方法来求最大面积的三角形,使用旋转卡壳方法。

枚举三角形的第一个顶点i,

然后初始第二个顶点j=i+1,第三个顶点k=j+1,

循环k+1直到Area(i,j,k)>Area(i,j,k+1)

更新面积的最大值,下面就开始旋转卡壳了(旋转j,k两个点)

(1)如果Area(i,j,k)<Area(i,j,k+1)且k!=i,则k=k+1,否则转(2)

(2)更新面积,j=j+1,如果j=i,跳出循环

这样旋转一圈后,求得的面积就是以i为顶点的最大三角形的面积了。

时间复杂度是O(n^2)

//author:chenkun
//旋转卡壳!
#include <cstdio>
#include 
<cmath>

const int MAXN=50000;

struct point {
    
int x, y;
} p[MAXN
+1],h[MAXN+1];

double distance(const point& p1,const point& p2) {
   
return sqrt( (p1.x - p2.x) * (p1.x - p2.x) + (p1.y - p2.y) * (p1.y - p2.y));
}

double multiply(const point& sp,const point& ep,const point& op) {
      
return((sp.x-op.x)*(ep.y-op.y)-(ep.x-op.x)*(sp.y-op.y));
}

int partition(point a[],int p,int r) {
    
int i=p,j=r+1,k;
    
double ang,dis;
    point R,S;
    k
=(p+r)/2;
    R
=a[p];
    a[p]
=a[k];
    a[k]
=R;
    R
=a[p];
    dis
=distance(R,a[0]);
    
while(1) {
        
while(1) {
            
++i;
            
if(i>r) {
                i
=r;
                
break;
            }
            ang
=multiply(R,a[i],a[0]);
            
if(ang>0)
                
break;
            
else if(ang==0) {
                
if(distance(a[i],a[0])>dis)
                
break;
            }
        }
        
while(1) {
             
--j;
            
if(j<p) {
                j
=p;
                
break;
            }
            ang
=multiply(R,a[j],a[0]);
            
if(ang<0)
                
break;
            
else if(ang==0) {
                
if(distance(a[j],a[0])<dis)
                    
break;
            }
        }
        
if(i>=j)break;
        S
=a[i];
        a[i]
=a[j];
        a[j]
=S;
    }
    a[p]
=a[j];
    a[j]
=R;
    
return j;
}

void anglesort(point a[],int p,int r) {
   
if(p<r) {
      
int q=partition(a,p,r);
      anglesort(a,p,q
-1);
      anglesort(a,q
+1,r);
   }
}

//对PointSet求凸包,点数为n,凸包上的点保存在ch中,点数位len
void Graham_scan(point PointSet[],point ch[],int n,int &len) {
    
int i,k=0,top=2;
    point tmp;
    
for(i=1;i<n;i++)
        
if ( PointSet[i].x<PointSet[k].x ||
            (PointSet[i].x
==PointSet[k].x) && (PointSet[i].y<PointSet[k].y) )
               k
=i;
    tmp
=PointSet[0];
    PointSet[
0]=PointSet[k];
    PointSet[k]
=tmp;
    anglesort(PointSet,
1,n-1);
    
if(n<3) {
        len
=n;
        
for(int i=0;i<n;i++) ch[i]=PointSet[i];
        
return ;
    }
    ch[
0]=PointSet[0];
    ch[
1]=PointSet[1];
    ch[
2]=PointSet[2];
    
for (i=3;i<n;i++) {
        
while (multiply(PointSet[i],ch[top],ch[top-1])>=0) top--;
            ch[
++top]=PointSet[i];
    }
    len
=top+1;
}

double Area(int i,int j,int k) {
    
return multiply(h[j],h[k],h[i])/2;
}

inline 
double max(double a,double b){return a>b?a:b;}
int main() {
    
int n;
    
while(scanf("%d",&n)!=EOF&&n!=-1) {
        
for(int i=0;i<n;i++)
            scanf(
"%d %d",&p[i].x,&p[i].y);
        
int len;
        Graham_scan(p,h,n,len);
        
//for(int i=0;i<len;i++)
        
//    printf("%d %d\n",h[i].x,h[i].y);

        
double ans=0;
        
for(int i=0;i<len;i++) {
            
int j=(i+1)%len;
            
int k=(j+1)%len;
            
while(k!=i&&Area(i,j,k)<Area(i,j,(k+1)%len)) {
                k
=(k+1)%len;
            }
            
if(k==i) continue;
            
int kk=(k+1)%len;
            
while(j!=kk&&k!=i) {
                ans
=max(ans,Area(i,j,k));
                
while(k!=i&&Area(i,j,k)<Area(i,j,(k+1)%len)) {
                    k
=(k+1)%len;
                }
                j
=(j+1)%len;
            }
        }

        printf(
"%.2f\n",ans);
    }
    
return 0;
}

posted on 2008-09-01 21:48  woodfish  阅读(982)  评论(3编辑  收藏  举报

导航