【筆記篇】最良心的計算幾何學習筆記（七）
動態凸包
本文的github傳送門在這裏~
======================================================================
不會凸包的趕緊去學一下哦~
======================================================================
好的我們已經會求凸包了. 那我們來看這樣一道題.
題目大意(英文題必須要有的翻譯過程OvO):
寫一個程序支持一下兩種操作:
- 1 x y 添加一個點\((x,y)\)
- 2 x y 詢問點\((x,y)\)是否在當前點集的凸包中
這就不是很好搞了吧?
所以我們就要動態維護.
我們來考慮一下加入一個點的時候會對之前求好的凸包造成什麽影響.
很容易從圖上看出:
- Case 1. 插入的點在凸包內, 不對凸包造成影響. (藍)
- Case 2. 凸包中其他點不動, 將新點插入. (紅)
- Case 3. 插入新點使得原凸包上的一些點不再在凸包上. (綠)
而考慮到我們的Graham掃描法, 我們可以知道這個就與新加入的點與兩側的邊的凹凸有關.
所以我們要找到新點的位置, 然後往兩邊掃. 因為只有插入, 所以點被刪了就不會再回來, 復雜度顯然就是\(O(1)\)的.
那麽如何找新點的位置呢? 我們在Graham掃描法中用到了極角排序. 那麽這裏我們也可以用極角來搞.
但是這裏的極角序還是和Graham掃描不太一樣... 因為Graham中可以找到\(y\)坐標最小的點, 但這裏顯然不能這麽做.(因為不一定在某些時間的凸包上).
我們要用做半平面交時候用的分\(x\)軸上下的極角排序(我就是不用atan2啦啦啦...)
那麽基準點就選最開始的凸包上的點就好了. 因為這個點一直都會在凸包裏.
由於題目中說開始一定會給一個不退化的三角形, 我們用重心就好啦.
這樣我們加入一個點的時候就:
- 插入這個點
- 找到這個點的前驅和後繼
- 分別往前和往後判斷是否滿足凸包性質(就是叉積嘛), 不是就刪點.
時間復雜度是查找\(O(logn)\)
而查詢的時候就找到前驅和後繼然後求一下叉積看看是不是在裏面就完了, 時間復雜度也是查找的\(O(logn)\)的, 所以最後的復雜度就是\(O(qlogn)\)的..
但是非常不想碼splay....
然後驚奇的發現這個題的操作只有插入、刪除、前驅、後繼, 那不就可以用set水過去了麽= =
於是就愉快地碼set.... (其實並不怎麽愉快, 因為STL還是有些地方挺反人類的)但是碼力太弱set版都寫不對還調了好久.... 顯然是藥丸.
要註意的地方就是我們要手動把set搞成環, 就是把end放在begin的前面, 因為極角排序本身是繞圈圈的.
大約就這樣吧.
代碼:
#include <set>
#include <cmath>
#include <cstdio>
using namespace std;
const double eps=1e-9;
inline int dcmp(const double &a){
if(fabs(a)<eps) return 0;
return a<0?-1:1;
}
struct vec{
double x,y;
vec(double X=0,double Y=0):x(X),y(Y){}
}p[4],sq;
inline vec operator-(const vec&a,const vec&b){
return vec(a.x-b.x,a.y-b.y);
}
inline double operator^(const vec&a,const vec&b){
return a.x*b.x+a.y*b.y;
}
inline double operator*(const vec&a,const vec&b){
return a.x*b.y-a.y*b.x;
}
inline double len(const vec&a){
return sqrt(a^a);
}
inline bool aboveX(const vec&a){
int d=dcmp(a*vec(1,0));
if(!d) return dcmp(a*vec(0,1))>0;
return d<0;
}
inline bool operator<(const vec&A,const vec&B){
vec a=A-sq,b=B-sq;
if(!dcmp(len(a))) return 1;
bool x1=aboveX(a),x2=aboveX(b);
if(x1!=x2) return x1;
int d=dcmp(a*b);
if(!d) return dcmp(len(b)-len(a))>0;
return d>0;
}
inline int gn(int a=0,char c=0,int f=1){
for(;(c<48||c>57)&&c!='-';c=getchar());if(c=='-')f=-1,c=getchar();
for(;c>='0'&&c<='9';c=getchar()) a=a*10+c-'0'; return a*f;
}
typedef set<vec> pset;
typedef pset::iterator pt;
pset s;
pt pre(pt x){
if(x==s.begin()) x=s.end();
return --x;
}
pt suc(pt x){
++x;
if(x==s.end()) return s.begin();
return x;
}
inline bool query(const vec&a){
pt it=s.lower_bound(a);
if(it==s.end()) it=s.begin();
int d=dcmp((a-*it)*(*pre(it)-*it));
if(!d) return dcmp(len(a-*it)-len(*pre(it)-*it))<0;
return d>0;
}
inline void inser(const vec&a){
if(query(a)) return;
pt it=s.lower_bound(a);
if(it==s.end()) it=s.begin();
s.insert(a);
while(s.size()>3&&dcmp((*it-a)*(*suc(it)-*it))<1){
s.erase(it); it=suc(s.find(a));
}
it=pre(s.find(a));
while(s.size()>3&&dcmp((*it-a)*(*pre(it)-*it))>-1){
s.erase(it); it=pre(s.find((a)));
}
}
int main(){
int n=gn(),m;
for(int i=1;i<=3;++i){
m=gn(),p[i].x=gn(),p[i].y=gn();
sq=sq-p[i];
}
sq=vec(-sq.x/3,-sq.y/3);
for(int i=1;i<=3;++i)
s.insert(p[i]);
for(int i=4;i<=n;++i){
int t=gn(),x=gn(),y=gn();
if(t==1) inser(vec(x,y));
else puts(query(vec(x,y))?"YES":"NO");
}
}你看一整道題碼完都比一個splay板子短= =
非常感謝寫STL的大牛們...
然後我們就可以去水題了.
當然一般的題目是沒有這麽裸的, 一般要維護一些東西, 比如周長啊 面積啊什麽的..
你看這道題就是這樣.
咦? 這題不是刪除麽= =
刪除是不好做的, 刪掉的點如果在凸殼上就要重新求一遍, 復雜度就上天了, 但是我們可以把操作倒著做, 視為不斷的加點和查詢凸殼長...
這道題目的出題人非常的良心, 只需要維護一個上凸殼, 而且還做了各種各樣的保證, 這樣就不用判邊界了.. 但是我還是WA了好多次, 原因也很扯淡我們一會再說.
這樣的話我們就考慮周長的維護.
很明顯地, 加入一個點勢必會導致打叉的邊被去除, 因為至少有一個左邊的點和一個右邊的點分別要連向新點. 要刪掉的點一定與要刪掉的邊是一一對應的關系. 這個從圖上也能很顯然地看出來. 這樣我們就可以復雜度不變的維護周長, 最後記得加上兩條藍色邊的長度即可.
由於這個題不需要判x軸上下而且特殊情況比較少, 所以寫起來會比較舒服.
起碼在我交之前是有這麽個感覺的.
#include <set>
#include <cmath>
#include <cstdio>
using namespace std;
const int N=101010;
const double eps=1e-9;
int dcmp(const double &a){
if(fabs(a)<eps) return 0;
return a<0?-1:1;
}
struct vec{
double x,y;
vec(double X=0,double Y=0):x(X),y(Y){}
}sq,p[N];
inline vec operator -(const vec &a,const vec &b){return vec(a.x-b.x,a.y-b.y);}
inline double operator ^(const vec &a,const vec &b){return a.x*b.x+a.y*b.y;}
inline double operator *(const vec &a,const vec &b){return a.x*b.y-a.y*b.x;}
inline double len(const vec &a){return sqrt(a^a);}
inline bool operator <(const vec &A,const vec &B){
if(A.x==B.x&&A.y==B.y) return 0;
vec a=A-sq,b=B-sq;
int d=dcmp(a*b);
if(!d){
int dd=dcmp(len(a)-len(b));
if(!dd) return dcmp(a*vec(0,1))>0;
return dd>0;
} return d>0;
}
inline int gn(int a=0,char c=0){
for(;c<'0'||c>'9';c=getchar());
for(;c>47&&c<58;c=getchar())a=a*10+c-48;return a;
}
set<vec> s;
typedef set<vec>::iterator pt;
int qu[N<<1]; double ans[N<<1],cur;
bool del[N];
pt pre(pt x){return --x;}
pt suc(pt x){return ++x;}
void inser(vec x){
pt it=s.lower_bound(x);
if(dcmp((x-*it)*(*pre(it)-*it))>-1) return;
cur-=len(*it-*pre(it));
while(suc(it)!=s.end()&&dcmp((x-*it)*(*it-*suc(it)))<1){
cur-=len(*it-*suc(it));
it=suc(it); s.erase(pre(it));
} it=pre(s.lower_bound(x));
while(it!=s.begin()&&dcmp((x-*it)*(*it-*pre(it)))>-1){
cur-=len(*it-*pre(it));
it=pre(it); s.erase(suc(it));
} s.insert(x); it=s.find(x);
cur+=len(*it-*pre(it))+len(*it-*suc(it));
}
int main(){
int n=gn(),x=gn(),y=gn(),m=gn(); sq=vec(n*0.5,0);
vec p1=vec(0,0),p2=vec(n,0),p3=vec(x,y); cur=len(p3-p2)+len(p3-p1);
s.insert(p1); s.insert(p2); s.insert(p3);
for(int i=1;i<=m;++i) p[i].x=gn(),p[i].y=gn(); int q=gn(),tot=0;
for(int i=1;i<=q;++i){
int t=gn();
if(t==1){
int w=gn();
del[w]=1; qu[i]=w;
}
else qu[i]=0;
}
for(int i=1;i<=m;++i) if(!del[i]) inser(p[i]);
for(int i=q;i;--i)
if(!qu[i]) ans[++tot]=cur;
else inser(p[qu[i]]);
for(int i=tot;i;--i) printf("%.2lf\n",ans[i]);
}然後有一個要註意的地方就是set的弱排序原理, 就是重載<號之後, 出現了\(x<y=false且y<x=false\)的情況, 那麽會視為\(x==y\)
所以把重載<的地方改得那麽鬼畜...
寫完過了個樣例交上去WA了... 然後就調唄.. 但是試了各種邊界數據都沒事...
甚至還把上面cf70D的代碼貼過來但還是WA..
最後發現: 竟然是
for(int i=1;i<=m;++i) if(!del[i]) inser(p[i]);裏面的\(p[i]\)打成了\(i\) !!!! 為什麽
void inser(vec p);的函數傳一個int變量不會報錯啊 !!!!!
C++的強制轉換太高端了吧...
不過這個題只維護一個上凸殼, 所以還有一種看上去好像更簡單一點的做法是用平衡樹來維護坐標, 這樣坐標的排序就是以\(x\)為第一關鍵字, 以\(y\)為第二關鍵字排序了..
但是最近覺得極角序非常地好用, 所以就寫極角序咯...
用極角序這種判個<就要慢好多的 竟然還能跑進第一頁也是挺神奇的...
有些時候凸包上的信息滿足一些特殊條件的時候還可以用cdq分治來維護什麽的...
據說代碼會比平衡樹好寫, 但是平衡樹如果寫set就不一樣了吧2333.....
不過這個應該很快就會去學...
那麽學完了再說...
【筆記篇】最良心的計算幾何學習筆記（七）