【AHOI2013复仇】动态凸包

Posted on 2013-02-28 18:29 Mato_No1 阅读(2701) 评论(0)  编辑 收藏 引用 所属分类: 经典问题的模型动态规划几何
原题地址

写了几天终于写出来了……(显然,我太弱了,请各位神犇不要鄙视)

在有加点的情况下,动态地维护凸包,有以下两种方法:
<1>维护上、下凸壳(本沙茶采用的方法):
凸包可以拆成上、下凸壳,对它们分别维护。两个凸壳均按照下面定义的<关系(即先x增、再y增)排序,注意,两个凸壳的两端是相同的,均为整个凸包的最小点与最大点,除两端外,它们没有公共定点。
以上凸壳为例,设目前加进去的点是P,则有以下三种情况:
1)P小于上凸壳的最小点(这里,对于点"(x1, y1)<(x2, y2)"定义为:x1<x2或x1=x2且y1<y2),此时将P插入上凸壳,并从P开始从小到大遍历新的上凸壳,将那些旋转方向不对的点删掉;
2)P大于上凸壳的最大点(这里,对于点"(x1, y1)>(x2, y2)"定义为:x1>x2或x1=x2且y1>y2),此时将P插入上凸壳,并从P开始从大到小遍历新的上凸壳,将那些旋转方向不对的点删掉;
3)P位于上凸壳的最小点与最大点之间:此时找到上凸壳中,小于等于P的最大点L和大于P的最小点R,判断PL转向PR是否为逆时针,若为逆时针,则P在上凸壳外,插入上凸壳,并分别向左右两个方向遍历,将旋转方向不对的点删掉;
对于下凸壳,1)和2)一样,3)反过来搞即可。注意在找前趋L和后继R的时候,可以顺便判断出来P是否在上、下凸壳上,若在凸壳上则不插入。
显然,这当中有插入、删除、找值的前趋和后继等操作,因此需要用平衡树。由于每个点最多只会被插入一次、删除一次,且遍历就意味着删除,因此总时间复杂度为O(NlogN)。为了加快速度,可以在平衡树(Splay Tree)中维护子树的最小、最大结点。
问题是,这种方法写起来是很复杂的,因为要维护两棵平衡树。

<2>极角法(Orz!!):
如果凸包面积大于0,可以在其内部选一个点作为定点,然后将凸包上所有点按照到这个定点的极角递增排序,用一棵平衡树维护。
加入新点P时,首先得出P到定点的极角,在平衡树中找到其前趋和后继,然后作叉积判断P是否在外部,若不在外部则不插入,否则插入,并且从这个前趋与这个后继开始,分别向两端删掉旋转方向不对的多余点。注意,虽然凸包是一个环,但由于按照极角排序了,在平衡树中仍只能作为链来维护,这样就会出现,在找前趋和后继的时候,如果在两端找不到,需要循环,到另一端去找,而且在遍历的时候,到头了也要回到另一端。
这种方法的时间复杂度也是O(NlogN),且由于只要维护一棵平衡树,代码量减小了很多(木有必要像下面的烂代码一样每个操作都加上bool _和N多的[_]了,虽然多了点循环的特判)。

当然,本题是可以不用平衡树的,而用线段树——因为可以离线。
还有一个问题就是本题的维护面积问题。面积的2倍可以用各边两端点关于任意定点的叉积之和得出,前提是端点按照逆时针排列(最后结果的正负只与端点排列顺序有关,与这个定点在哪无关),因此,对于有边加入或删除的时候,就可以顺便维护出来,当然在有三角形变化的时候可以直接维护三角形。

动态凸包的最主要应用就是那些用凸壳或曲线凸壳优化的DP,虽然这种DP大多数时候只需要维护上或下凸壳就行了囧……比较典型的应用有:
NOI2007 cash(动态维护凸壳,但本题可以用分治转成离线,从而直接Graham构造解决)
NOI2009 poet(维护曲线凸壳,容易证明这种绝对值P次函数的曲线符合交点只有一个的性质,且对称轴是越来越右,所以可以直接用队列维护……MS它是一种有用模型的代表)
BZOJ2149 拆迁队(由于合法决策要满足三个条件:决策阶段较小、A值较小、F值严格比目前阶段F值小1,因此需要变序,按照A的递增序插入决策射线,维护可以离散化后用线段树,也可以根据A值即为斜率且递增,直接用栈维护凸壳。当然,本题也需要分治转成离线)
CEOI2011 ballons(本题的效果值曲线为顶点在X轴上的开口向上的抛物线,由于对称轴x是递增的,所以可以只维护右半部分,一个决策插入后影响的是一个区间,可以离散化后用线段树解决。当然,本题如果对称轴x不递增才好玩呢囧……)

代码:
#include <iostream> 
#include 
<stdio.h> 
#include 
<stdlib.h> 
#include 
<string.h> 
using namespace std; 
#define re(i, n) for (int i=0; i<n; i++) 
#define re1(i, n) for (int i=1; i<=n; i++) 
#define re2(i, l, r) for (int i=l; i<r; i++) 
#define re3(i, l, r) for (int i=l; i<=r; i++) 
#define rre(i, n) for (int i=n-1; i>=0; i--) 
#define rre1(i, n) for (int i=n; i>0; i--) 
#define rre2(i, r, l) for (int i=r-1; i>=l; i--) 
#define rre3(i, r, l) for (int i=r; i>=l; i--) 
#define ll long long 
const int MAXN = 100010
struct poi { 
    ll x, y; 
    
bool operator< (poi p0) const {return x < p0.x || x == p0.x && y < p0.y;} 
    
bool operator> (poi p0) const {return x > p0.x || x == p0.x && y > p0.y;} 
    poi 
operator- (poi p0) {return (struct poi) {x - p0.x, y - p0.y};} 
}; 
struct node { 
    poi v; 
    
int c[2], p, sz, minNo, maxNo; 
    
bool d; 
} T[
2][MAXN]; 
int N[2], root[2]; 
ll res; 
void sc(bool _, int _p, int _c, bool _d) 

    T[_][_p].c[_d] 
= _c; T[_][_c].p = _p; T[_][_c].d = _d; 

void upd(bool _, int No) 

    
int lch = T[_][No].c[0], rch = T[_][No].c[1]; 
    T[_][No].sz 
= T[_][lch].sz + T[_][rch].sz + 1
    T[_][No].minNo 
= lch ? T[_][lch].minNo : No; T[_][No].maxNo = rch ? T[_][rch].maxNo : No; 

void rot(bool _, int No) 

    
int p = T[_][No].p; bool d = T[_][No].d; 
    
if (p == root[_]) T[_][root[_] = No].p = 0else sc(_, T[_][p].p, No, T[_][p].d); 
    sc(_, p, T[_][No].c[
!d], d); sc(_, No, p, !d); upd(_, p); 

void splay(bool _, int No, int r) 

    
int p; while ((p = T[_][No].p) != r) if (T[_][p].p == r) rot(_, No); else if (T[_][p].d == T[_][No].d) rot(_, p), rot(_, No); else rot(_, No), rot(_, No); upd(_, No); 

void ins(bool _, poi v0) 

    
if (!root[_]) { 
        T[_][root[_] 
= ++N[_]].p = 0; T[_][N[_]].c[0= T[_][N[_]].c[1= 0; T[_][N[_]].sz = 1; T[_][N[_]].minNo = T[_][N[_]].maxNo = N[_]; T[_][N[_]].v = v0; return
    } 
    poi v; 
int i = root[_], j; 
    
while (1) { 
        v 
= T[_][i].v; T[_][i].sz++; j = T[_][i].c[v0 > v]; if (j) i = j; else break
    } 
    T[_][
++N[_]].c[0= T[_][N[_]].c[1= 0; T[_][N[_]].sz = 1; T[_][N[_]].v = v0; T[_][N[_]].minNo = T[_][N[_]].maxNo = N[_]; sc(_, i, N[_], v0 > v); 
    splay(_, N[_], 
0); 

int uni(bool _, int A, int B) 

    
if (!A) return B; else if (!B) return A; else if ((A + B) & 1) { 
        
int S = uni(_, A, T[_][B].c[0]); 
        sc(_, B, S, 
0); upd(_, B); return B; 
    } 
else { 
        
int S = uni(_, T[_][A].c[1], B); 
        sc(_, A, S, 
1); upd(_, A); return A; 
    } 

int L(bool _, poi v0) 

    
int i = root[_], j, res0 = -1; poi v; 
    
while (1) { 
        v 
= T[_][i].v; 
        
if (v0 < v) j = T[_][i].c[0]; else if (v0 > v) {j = T[_][i].c[1]; res0 = i;} else return i; 
        
if (j) i = j; else break
    } 
    
return res0; 

int R(bool _, poi v0) 

    
int i = root[_], j, res0 = -1; poi v; 
    
while (1) { 
        v 
= T[_][i].v; 
        
if (v0 < v) {j = T[_][i].c[0]; res0 = i;} else if (v0 > v) j = T[_][i].c[1]; else return i; 
        
if (j) i = j; else break
    } 
    
return res0; 

ll cr(poi p0, poi p1) 

    
return p0.x * p1.y - p0.y * p1.x; 

void solve(poi p) 

    
int _L = L(0, p); if (_L >= 0 && !(T[0][_L].v < p || p < T[0][_L].v)) return
    
int _R = R(0, p), _L0, _R0, _; bool FF; ll cr0; 
    
if (_L == -1) { 
        res 
+= cr(T[0][_R].v, p); FF = 1
    } 
else if (_R == -1) { 
        res 
+= cr(p, T[0][_L].v); FF = 1
    } 
else { 
        cr0 
= cr(T[0][_L].v - p, T[0][_R].v - p); 
        
if (cr0 > 0) {res += cr0; FF = 1;} else FF = 0
    } 
    
if (FF) { 
        ins(
0, p); 
        
while (_L = T[0][root[0]].c[0]) { 
            _L0 
= T[0][_L].maxNo; splay(0, _L0, root[0]); _L = T[0][root[0]].c[0]; 
            
if (T[0][_L].c[0]) _L0 = T[0][T[0][_L].c[0]].maxNo; else break
            cr0 
= cr(T[0][_L].v - p, T[0][_L0].v - p); 
            
if (cr0 <= 0) {res -= cr0; _ = uni(0, T[0][_L].c[0], T[0][_L].c[1]); sc(0, root[0], _, 0); upd(0, root[0]);} else break
        } 
        
while (_R = T[0][root[0]].c[1]) { 
            _R0 
= T[0][_R].minNo; splay(0, _R0, root[0]); _R = T[0][root[0]].c[1]; 
            
if (T[0][_R].c[1]) _R0 = T[0][T[0][_R].c[1]].minNo; else break
            cr0 
= cr(T[0][_R].v - p, T[0][_R0].v - p); 
            
if (cr0 >= 0) {res += cr0; _ = uni(0, T[0][_R].c[0], T[0][_R].c[1]); sc(0, root[0], _, 1); upd(0, root[0]);} else break
        } 
    } 
    _L 
= L(1, p); if (_L >= 0 && !(T[1][_L].v < p || p < T[1][_L].v)) returnelse _R = R(1, p); 
    
if (_L == -1) { 
        res 
+= cr(p, T[1][_R].v); FF = 1
    } 
else if (_R == -1) { 
        res 
+= cr(T[1][_L].v, p); FF = 1
    } 
else { 
        cr0 
= cr(T[1][_L].v - p, T[1][_R].v - p); 
        
if (cr0 < 0) {res -= cr0; FF = 1;} else FF = 0
    } 
    
if (FF) { 
        ins(
1, p); 
        
while (_L = T[1][root[1]].c[0]) { 
            _L0 
= T[1][_L].maxNo; splay(1, _L0, root[1]); _L = T[1][root[1]].c[0]; 
            
if (T[1][_L].c[0]) _L0 = T[1][T[1][_L].c[0]].maxNo; else break
            cr0 
= cr(T[1][_L].v - p, T[1][_L0].v - p); 
            
if (cr0 >= 0) {res += cr0; _ = uni(1, T[1][_L].c[0], T[1][_L].c[1]); sc(1, root[1], _, 0); upd(1, root[1]);} else break
        } 
        
while (_R = T[1][root[1]].c[1]) { 
            _R0 
= T[1][_R].minNo; splay(1, _R0, root[1]); _R = T[1][root[1]].c[1]; 
            
if (T[1][_R].c[1]) _R0 = T[1][T[1][_R].c[1]].minNo; else break
            cr0 
= cr(T[1][_R].v - p, T[1][_R0].v - p); 
            
if (cr0 <= 0) {res -= cr0; _ = uni(1, T[1][_R].c[0], T[1][_R].c[1]); sc(1, root[1], _, 1); upd(1, root[1]);} else break
        } 
    } 

int main() 

    
int m; poi p0, p1, p2, _; 
    scanf(
"%lld%lld%lld%lld%lld%lld%d"&p0.x, &p0.y, &p1.x, &p1.y, &p2.x, &p2.y, &m); 
    
if (p1 < p0) {_ = p0; p0 = p1; p1 = _;} if (p2 < p0) {_ = p0; p0 = p2; p2 = _;} if (p2 < p1) {_ = p1; p1 = p2; p2 = _;} 
    ins(
0, p0); ins(0, p2); ins(1, p0); ins(1, p2); ll __ = cr(p0 - p1, p2 - p1); 
    
if (__ > 0) ins(0, p1); else if (__ < 0) ins(1, p1); 
    res 
= __ >= 0 ? __ : -__; 
    re(i, m) { 
        scanf(
"%lld%lld"&_.x, &_.y); 
        solve(_); 
        printf(
"%lld\n", res >= 0 ? res : -res); 
    } 
    
return 0
}

只有注册用户登录后才能发表评论。
【推荐】超50万行VC++源码: 大型组态工控、电力仿真CAD与GIS源码库
网站导航: 博客园   IT新闻   BlogJava   知识库   博问   管理