[HNOI2007][BZOJ1185] 最小矩形覆蓋 [凸包+旋轉卡殼]

題面

BZOJ題面php

前置芝士

建議先學習向量相關的計算幾何基礎html

計算幾何基礎戳這裏node

思路

用這道題學習一下凸包和旋轉卡殼ios

首先是凸包部分git

凸包

求凸包用的算法是graham算法算法

算法流程以下:學習

找到$y$座標最小的一點做爲原點spa

對原點以外的全部點按照到原點的極角排序(這裏由於選取了最靠下的,因此極角範圍在$[0,\pi]$).net

依次遍歷全部排序後的點,加入一個單調棧中:每次判斷(棧頂元素和棧頂第二元素之間的斜率)是否大於(當前點和棧頂第二元素之間的斜率)code

注意一旦這個大於成立了,棧頂元素就會在當前元素和棧頂第二元素的連線的「下面」,也就是在凸包裏面了

由於咱們事先按照極角排序了,因此這一單調棧能夠不重複不遺漏地記錄凸包上全部點

注意這樣求出來的凸包上的點是逆時針排序的(根本緣由是由於極角排序就是逆時針繞圈)

graham算法的複雜度是$O(n\log n)$,瓶頸是排序

旋轉卡殼

首先,我在這道題裏面用的不是標準的旋轉卡殼算法......可是也是「旋轉+卡殼」的思路

標準版的旋轉卡殼戳這裏,這裏標準版指的是用4條邊去卡殼,我寫的是一條邊和3個極值點

對於最小面積矩形,咱們有結論:這一外接矩形必定有一條邊和凸包的一條邊重合

注意這個結論對於最小周長矩形依然成立

證實嘛......我愣是沒找到。感性理解一下就是對於一個四邊都接在凸包端點上的矩形,把它旋轉一下必定更優

核心算法流程以下:

咱們遍歷凸包上的每一條邊,並對於每一條邊求出以這條邊爲x軸時,最靠左的點、最靠右的點、最靠上的點

設求出來的凸包上的點是$q[0...m-1]$

假設咱們當前的邊的兩個端點是凸包上的$q[i],q[i+1]$,並且是有向的(i指向i+1),那麼上述三個端點有以下性質:

對於最靠左的點$q[l]$,$vec(q[i],q[l])\ast vec(q[i],q[i+1])$是全部$l$中最小的

對於最靠右的點$q[r]$,$vec(q[i],q[r])\ast vec(q[i],q[i+1])$是全部$r$中最大的

對於最靠上的點$q[u]$,$vec(q[i],q[i+1])$叉乘$vec(q[i],q[u])$是全部$u$中最大的

其中$vec(u,v)$表示從點$u$指向$v$的向量

前兩個的證實,利用點乘的性質:由於點乘的被投影向量長度相等,因此決定點乘結果大小的就是投影的大小

那麼顯然投影最小最靠左,投影最大最靠右

第三個的證實,利用叉乘的性質:叉乘等於兩個向量逆時針旋轉構成的有向平行四邊形面積

由於平行四邊形底邊長度相同,並且$vec(q[i],q[i+1])$必定在全部從$q[i]$出發的向量的順時針方向,因此反過來旋轉必定是正的,叉乘最大就是最高

圖示以下:

知道了這三個點之後,咱們就能夠知道這個矩形的長寬,進而求出面積了

又有性質:咱們每次從$vec(q[i],q[i+1])$旋轉到$vec(q[i+1],q[i+2])$的時候,$l,r,u$也會跟着逆時針旋轉,因此只要枚舉便可

這樣,整個旋轉卡殼就是「遍歷全部邊」+「順序求出端點」的過程,總複雜度是$O(n)$的

Code

代碼中node是向量結構體,對於點咱們用位矢表示,也是向量

標$\ast$的是叉乘,標\的是點乘

#include<iostream>
#include<cstdio>
#include<cstring>
#include<algorithm>
#include<cassert>
#include<cmath>
#define eps 1e-9
#define ll long long
using namespace std;
inline int read(){
    int re=0,flag=1;char ch=getchar();
    while(!isdigit(ch)){
        if(ch=='-') flag=-1;
        ch=getchar();
    }
    while(isdigit(ch)) re=(re<<1)+(re<<3)+ch-'0',ch=getchar();
    return re*flag;
}
struct node{
    double x,y;
    node(double xx=0.0,double yy=0.0){x=xx;y=yy;}
    inline bool operator <(const node &b){return ((fabs(y-b.y)<eps)?(x<b.x):(y<b.y));}
    inline friend bool operator ==(const node &a,const node &b){return ((fabs(a.x-b.x)<eps)&&(fabs(a.y-b.y)<eps));}
    inline friend bool operator !=(const node &a,const node &b){return !(a==b);}
    inline friend node operator +(const node &l,const node &r){return node(l.x+r.x,l.y+r.y);}
    inline friend node operator -(const node &l,const node &r){return node(l.x-r.x,l.y-r.y);}
    inline friend node operator *(node l,double r){return node(l.x*r,l.y*r);}
    inline friend double operator *(const node &l,const node &r){return l.x*r.y-l.y*r.x;}
    inline friend double operator /(const node &l,const node &r){return l.x*r.x+l.y*r.y;}
    inline friend double dis(const node &a){return sqrt(a.x*a.x+a.y*a.y);}
}a[100010],q[100010],x[10];
int n,top;double ans=1e60;
inline bool cmp(node l,node r){
    double tmp=(a[1]-l)*(a[1]-r);
    if(fabs(tmp)<eps) return dis(a[1]-l)<dis(a[1]-r);
    else return tmp>0;
}
void graham(){//get a counter-clockwise convex
    int i;
    for(i=2;i<=n;i++){
        if(a[i]<a[1]) swap(a[1],a[i]);
    }
    sort(a+2,a+n+1,cmp);
    q[++top]=a[1];
    q[++top]=a[2];
    for(i=3;i<=n;i++){
        while(top>1&&((q[top]-q[top-1])*(a[i]-q[top])<eps)) top--;
        q[++top]=a[i];
    }
    q[0]=q[top];
}
void RC(){//RotatingCalipers
    int l=1,r=1,p=1,i;
    double L,R,D,H,tmp;
    for(i=0;i<top;i++){
        D=dis(q[i]-q[i+1]);
        while((q[i+1]-q[i])*(q[p+1]-q[i])-(q[i+1]-q[i])*(q[p]-q[i])>-eps) p=(p+1)%top;
        while((q[i+1]-q[i])/(q[r+1]-q[i])-(q[i+1]-q[i])/(q[r]-q[i])>-eps) r=(r+1)%top;
        if(i==0) l=r;
        while((q[i+1]-q[i])/(q[l+1]-q[i])-(q[i+1]-q[i])/(q[l]-q[i])<eps) l=(l+1)%top;
        L=(q[i+1]-q[i])/(q[l]-q[i])/D;
        R=(q[i+1]-q[i])/(q[r]-q[i])/D;
        H=(q[i+1]-q[i])*(q[p]-q[i])/D;
        tmp=(R-L)*H;
        if(tmp<ans){
            ans=tmp;
            x[0]=q[i]+(q[i+1]-q[i])*(R/D);
            x[1]=x[0]+(q[r]-x[0])*(H/dis(x[0]-q[r]));
            x[2]=x[1]-(x[0]-q[i])*((R-L)/dis(q[i]-x[0]));
            x[3]=x[2]-(x[1]-x[0]);
        }
    }
}
int main(){
    n=read();int i,j;
    for(i=1;i<=n;i++) scanf("%lf%lf",&a[i].x,&a[i].y);
    graham();
    RC();
    printf("%.5lf\n",ans);
    j=0;
    for(i=1;i<4;i++) if(x[i]<x[j]) j=i;
    for(i=0;i<4;i++) printf("%.5lf %.5lf\n",x[j].x,x[j].y),j=(j+1)%4;
}
相關文章
相關標籤/搜索