动态凸壳

2023-12-12 19:58:45

我并不会求静态二维凸包，但是在培训中碰到几个需要使用平衡树来维护凸包的计算几何题

二维凸包

凸包：求一个周长最小的，并且能够包含所有给定点的多边形。当多边形表面存在凹陷时，根据三角不等式\(\begin{cases}a+b>c\\b+c>a\\a+c>b\end{cases}\)，一定没有直接把最短边连起来优

维护动态二维凸包

向量叉积

\(\vec{a}=(x_1,y_1),\vec{b}=(x_2,y_2)\)

\(\vec{a}\times\vec{b}=\vec{c}\) 其等价于 \(|a||b|\sin\theta = x_1y_2-x_2y_1\)

其几何本质就是将 \(\vec{a}\) 按照逆时针方向旋转至 \(\vec{b}\)

所以可以根据这个东西判断两点之前的角度关系

atan2(y,x)

把 y 和 x 坐标放入该函数，可以得到这个点和原点形成的直线与 x 轴的夹角，对于向量来说也是同理

平衡树维护凸壳

首先自定义一个点和向量的结构体以及相关运算

struct Node{// 其中还有点乘运算没有定义
	int x,y;
	db Ang;
	friend Node operator+(Node a,Node b){return (Node){a.x + b.x , a.y + b.y,0};}//向量加减法运算
	friend Node operator-(Node a,Node b){return (Node){a.x - b.x , a.y - b.y,0};}
	inline db operator~()const {return sqrt(x*x + y*y);}//计算向量模长
	inline db ang() const {return atan2(y,x);}//算向量的极角
	friend bool operator<(Node a,Node b)//极角排序，中心点按重心即可
	{
		return a.Ang < b.Ang;
	}
}po[np],z;
inline int cross(Node a,Node b)//向量叉积
{
	return a.x * b.y - a.y * b.x;
}
inline db dis(Node a,Node b)
{
	int x_ = (a-b).x;
	int y_ = (a-b).y;
	x_ /= 3;//为了避免算重心/3带来的精度问题，所以给所有坐标都乘上3
	y_ /= 3;
	return sqrt(x_*x_+y_*y_);
}

我们的目的是动态维护凸壳，用平衡树来维护。首先考虑已经有了一个凸壳，往上面加一个点后如何维护。

先计算该点与重心之间的极角，然后在 set 中二分查找第一个极角比自己大的点。此时取该点、第一个极角比自己大的点、最后一个极角比自己小的点这三个点，记为x,suf(x),pre(x)，计算叉积 \(cross(pre(x)-x,suf(x)-x)\)，若叉积大于等于 0，则说明 x 位于凸壳内侧或者凸壳上，set 保持不变。否则说明该点已经超出凸壳，将该点插入凸壳，然后不断更新前驱后继是否继续满足凸包的性质即可。

给出 Insert 代码


inline set<Node>::iterator suf_(set<Node>::iterator x)
{
	if(x == s.end()) return s.begin();
	return (++x) == s.end()?s.begin():x;
}
inline set<Node>::iterator pre_(set<Node>::iterator x){return x == s.begin()?--s.end():(--x);}

inline void Insert(Node x)
{
	sto it = s.lower_bound(x);
	set<Node>::iterator a1 = pre_(it);
	set<Node>::iterator a2 = (it == s.end())?s.begin():it;
	if(cross(*a1 -  x,*a2 - x) >= 0) return;
	s.insert(x);
	it = s.find(x);
	while(s.size() >= 4&&cross(*it - (*suf_(it)) , (*suf_(suf_(it))) - *suf_(it)) >= 0)s.erase(suf_(it));
	while(s.size() >= 4 && cross((*pre_(pre_(it))) - *pre_(it) , *it - (*pre_(it))) >= 0) s.erase(pre_(it));
}

然后就完成了……

防线修建

题目大概就是要求维护一个凸壳，并且随时回答该凸壳的长度。题目要求的凸壳是删点的，但我们的凸壳只能维护加点，所以需要时光倒流将删点变为加点。

detail：注意增量维护凸壳长度时别把点加反了

#include<bits/stdc++.h>

using namespace std;

#define pii pair<int,int>
#define INF 1ll<<30
#define db double

template<typename _T>
inline void read(_T &x)
{
	x=0;char s=getchar();int f=1;
	while(s<'0'||s>'9') {f=1;if(s=='-')f=-1;s=getchar();}
	while('0'<=s&&s<='9'){x=(x<<3)+(x<<1)+s-'0';s=getchar();}
	x*=f;
}

const db Eps = 1e-12;
const int np = 2e5 + 5;

struct Node{
	int x,y;
	db Ang;
	friend Node operator+(Node a,Node b){return (Node){a.x + b.x , a.y + b.y,0};}
	friend Node operator-(Node a,Node b){return (Node){a.x - b.x , a.y - b.y,0};}
	inline db operator~()const {return sqrt(x*x + y*y);}
	inline db ang() const {return atan2(y,x);}
	friend bool operator<(Node a,Node b)
	{
		return a.Ang < b.Ang;
	}
}po[np],z;
Node sta[np];
int top;
set<Node> s;
inline int cross(Node a,Node b)
{
	return a.x * b.y - a.y * b.x;
}

inline db dis(Node a,Node b)
{
	int x_ = (a-b).x;
	int y_ = (a-b).y;
	x_ /= 3;
	y_ /= 3;
	return sqrt(x_*x_+y_*y_);
}

struct query{
	int typ,pos;
}q[np];
int vis[np];
db ans[np];
pii city[np];

inline set<Node>::iterator suf_(set<Node>::iterator x)
{
	if(x == s.end()) return s.begin();
	return (++x) == s.end()?s.begin():x;
}
inline set<Node>::iterator pre_(set<Node>::iterator x){return x == s.begin()?--s.end():(--x);}
db Edge;

inline void Insert(Node x)
{
	set<Node>::iterator it = s.lower_bound(x);
	set<Node>::iterator a1 = pre_(it);
	set<Node>::iterator a2 = it == s.end()?s.begin():it;
	if(cross(*a1 - x,*a2 - x) >= 0) return ;
	s.insert(x);
	it = s.find(x);
	top = 0;
	
	while(s.size() >= 4 && cross(*it - (*suf_(it)), (*suf_(suf_(it))) - (*suf_(it))) >= 0)sta[++top] = *suf_(it),s.erase(suf_(it));
	reverse(sta + 1,sta + 1 + top);
	while(s.size() >= 4 && cross(*pre_(pre_(it)) - (*pre_(it)),*it - (*pre_(it))) >= 0)sta[++top] = *pre_(it),s.erase(pre_(it));
	if(top)
	{
		Edge -= dis(*suf_(it),sta[1]);
		for(int i=2;i<=top;i++) Edge -= dis(sta[i - 1],sta[i]);
		Edge -= dis(*pre_(it),sta[top]);
	}
	else Edge -= dis(*pre_(it),*suf_(it));
	Edge += dis(*suf_(it),*it);
	Edge += dis(*pre_(it),*it);
}

signed  main()
{
	int n,xx,yy,m,q_;
	read(n); 
	read(xx);
	read(yy);
	int zx,zy;
	po[1] = (Node){0,0,0};
	po[2] = (Node){3 * n,0,0};
	po[3] = (Node){3 * xx,3 * yy,0};
	zx = (n + xx);
	zy = yy;
	z = (Node){zx,zy,0};
	Node k = po[1] - z;
	po[1].Ang = k.ang();
	k = po[2] - z;
	po[2].Ang = k.ang();
	k = po[3] - z;
	po[3].Ang = k.ang();
	s.insert(po[1]);
	s.insert(po[2]);
	s.insert(po[3]);
	db c = dis(po[1],po[2]);
	Edge = 0;
	for(set<Node>::iterator it = s.begin();it != s.end();it++) 
	{
		set<Node>::iterator ier = ++it;
		if(ier == s.end()) ier = s.begin();
		it--;
		Edge += dis(*it,*ier);
	}
	read(m);
	for(int i=1,x,y;i<=m;i++)
	{
		read(x);
		read(y);
		city[i] = (pii){3 * x,3 * y};
	}
	read(q_);
	for(int i=1,opt,pos;i<=q_;i++)
	{
		read(opt);
		if(opt == 1) read(pos),vis[pos] = 1;
		else pos = 0;
		q[i] = (query){opt,pos};
	}
	for(int i=m;i>=1;i--)
	{
		if(!vis[i])
		{
			Node p = (Node){city[i].first,city[i].second,0};
			k = p - z;
			p.Ang = k.ang();
			Insert(p);
		}
	}
	for(int i=q_;i>=1;i--)
	{
		if(q[i].typ == 1)
		{
			int pos = q[i].pos;
			Node p = (Node){city[pos].first,city[pos].second,0};
			k = p - z;
			p.Ang = k.ang();
			Insert(p);
		}
		else ans[i] = Edge;
	}
	
	for(int i=1;i<=q_;i++)
	{
		if(q[i].typ == 2)
		{
			printf("%.2lf\n",ans[i] - c);
		}
	}
 }

\(End\)

计算几何的第一个知识

码农公寓