计算几何学习笔记

以下的图几乎是盗来的,因为我懒得画图了。

如果想学计算几何建议看洛谷日报那篇,这个是自己用来复习的,十分粗略。

向量

向量是既有大小又有方向的量,可以用一个带有箭头的线段表示,表示为 \(\vec{a}\),由于具有平移不变性,所以可以将起点平移到原点,用一个坐标 \(( x , y)\) 来表示一个向量。

向量的模就是向量的长度,对于 \(\vec{a}=( x , y)\) ,表示为 \(|\vec{a}|=\sqrt{x^2+y^2}\)

向量的加减满足交换律,对于 \(\vec{a}=( x_0 , y_0)\ ,\ \vec{b}=( x_1 , y_1)\)\(\vec{a}\pm \vec{b}=( x_0\pm x_1 , y_0\pm y_1)\)

盗图

向量的数量积,也称点积,是一个实数,几何意义为一个向量在另一个向量上的投影再乘上第二个向量的模长,计算式为
\(\vec{a}\cdot\vec{b}=|\vec{a}||\vec{b}|\cos\theta(\theta=\left<\vec{a},\vec{b}\right>)\) ,对于 \(\vec{a}=(x_1,y_1)\ ,\ \vec{b}=(x_2,y_2)\ ,\ \vec{a}\cdot\vec{b}=x_1x_2+y_1y_2\) ,满足交换律。

由点积可以判断并计算两个向量的夹角,大于 \(90^\circ\) 时点积为正,等于 \(90^\circ\) 时点积为零,小于 \(90^\circ\) 时点积为负。

向量的外积,也称叉积,是一个向量,几何意义是两向量由平行四边形法则围成的面积,计算式为 \(\vec{a}\times\vec{b}=|\vec{a}||\vec{b}|\sin\theta\) ,对于 \(\vec{a}=(x_1,y_1)\ ,\ \vec{b}=(x_2,y_2)\ ,\ \vec{a}\times\vec{b}=x_1y_2-x_2y_1\) ,不满足交换律。

向量的叉积的正负可以理解为从 \(a\) 逆时针转向 \(b\) 的夹角的范围,当夹角为平角则为正,否则为负。

所以可以通过叉积的计算来判断两个向量的旋转关系,由于其几何意义以及正负的性质,也可以应用于求多边形的面积,做法是算出所有边的两端的向量的叉积再除以二。

struct Vector{
	int x,y;
	Vector friend operator -(Vector a,Vector b){ return {a.x-b.x,a.y-b.y}; }
	Vector friend operator +(Vector a,Vector b){ return {a.x+b.x,a.y+b.y}; }
	Vector friend operator *(Vector a,int b){ return {a.x*b,a.y*b}; }
	int friend operator *(Vector a,Vector b){ return a.x*b.y-a.y*b.x; }//叉积
	int friend operator ^(Vector a,Vector b){ return a.x*b.x+a.y*b.y; }//点积
};

线段交点

直接盗图和别人的解释吧。

盗图

因为 \(\overrightarrow{AO}=\frac{|AO|}{|AB|}\cdot \overrightarrow{AB}=\frac{|AO|}{|AO|+|BO|}\cdot \overrightarrow{AB}\) ,令原点为 \(O^{'}\),则 \(\overrightarrow{O'O}=\overrightarrow{O'A}+\overrightarrow{AO}\)\(\overrightarrow{O'O}\) 的坐标即为 \(O\) 的坐标。

Vector cross(Line a,Line b){
	return b.pre+b.v*(((b.pre-a.pre)*a.v)/(a.v*b.v));
}

凸包

有一个比较普遍的做法,就是将点按 \(x\) 排序,然后从最左边的一个点出发,不断将点 \(A\) 加入栈,设 \(B,C\) 为栈顶的两个点,如果 \(BA\)\(BC\) 的下面,那么 \(C\) 一定不是凸包上的点,就可以将 \(C\) 弹出栈。

盗图

对于判断,可以通过向量的叉积实现,如果 \(BA\)\(BC\) 的下面,那么从 \(BC\) 转向 \(BA\) 的夹角一定不是平角, \(\vec{BC}\times\vec{BA}\le 0\)

int dis(Vector a){ return a.x*a.x+a.y*a.y; }
bool cmp(Vector a,Vector b){
	int x=(a-d[1])*(b-d[1]);
	return x>0||(x==0&&dis(a-d[1])<dis(b-d[1]));
}

int k=1,top=2;
void Get(){
	for(int i=2;i<=n;i++) if(d[i].y<d[k].y||(d[i].y==d[k].y&&d[i].x<d[k].x)) k=i;
	swap(d[1],d[k]);sort(d+2,d+n+1,cmp);s[1]=d[1],s[2]=d[2];
	for(int i=3;i<=n;i++){
		while(top>=2&&(s[top]-s[top-1])*(d[i]-s[top-1])<=0) top--;
		s[++top]=d[i];
	}
}

旋转卡壳

先了解对踵点,如果过凸多边形上两点作一对平行线,使得整个多边形都在这两条线之间,那么这两个点被称为一对对踵点。

再了解凸多边形的直径,就是凸包上距离最远的两个点的距离,通过反证法可以发现这样的点对一定是对踵点。

至于旋转卡壳。。。

盗图

可以发现平行线同时接触的两个点就是对踵点,所以旋转卡壳就是在枚举对踵点。

由于是凸多边形,一个点要求另外一个点组成对踵点可以通过三分法,单个查询 \(\mathcal{o(\log n)}\) ,总复杂度 \(\mathcal{o(n\log n)}\)

不过可以发现,对踵点不会回退,画个图大概就可以理解,所以可以使用双指针来遍历凸包,判断可以通过计算三角形的面积,由于对踵点的定义是关于平行线的,所以对于凸包上的点的一条边,另外一个对踵点一定是所有点中到那条边的垂直距离最大的,至于三角形面积可以通过向量的叉积来计算。

int GetMax(){//求直径
	int res=0;
	if(top==2) return dis(s[1]-s[2]);
	s[++top]=s[1];
	for(int i=0,j=2;i<top;i++){
		while((s[i]-s[j])*(s[i+1]-s[j])<(s[i]-s[j+1])*(s[i+1]-s[j+1])) j=(j+1==top)?1:j+1;
		res=max(res,max(dis(s[i]-s[j]),dis(s[i+1]-s[j])));
	}
	return res;
}

半平面交

半平面的定义是一条直线和直线的一侧,是一个点集,包含直线时称闭半平面,不包含时称开半平面,在计算几何中用向量表示,整个题统一以向量的左侧或右侧为半平面。

半平面交就是多个半平面交的交集,如果面积有限,显然这是个凸包。

半平面交的求解需要先极角排序,通过C++自带的函数 atan2(double y,double x) 可以返回 \(θ=arctan(\frac{y}{x})\) ,如果用向量来表示边的话大概就是向量的斜率。

如果要求更高的精度,可以自己写一个比较函数,比较两个向量时可以将交点看作原点,利用点积判断夹角就可以了。

由于一个凸包连续的边的斜率是单调的,所以极角排序后的半平面组成的半平面交从第一个开始连续的边在原来的排序中的下标也是单调的,那么就可以保证后面操作的正确性。

先考虑像做凸包一样来做半平面交,利用一个单调栈,然后往里头加半平面,大概像是这样:

红色的为加入的半平面,显然第三个半平面会被踢出单调栈,直观的判断是与围成的半平面交找到一个交点,在这个大于这个交点横坐标的半平面交上的点对应的半平面都要踢出,显然除第一个外一个半平面只会贡献一个点。

考虑这个东西怎么写比较好,仔细观察,发现被踢出的半平面和加入的半平面的交点处于加入的半平面的异侧,这个可以自行证明。

//[CQOI2006]凸多边形 部分代码
int dcmp(db x){ return fabs(x)<=eps?0:(x>0?1:-1); }//精度专属判断
struct Vector{
	db x,y;
	Vector friend operator -(Vector a,Vector b){ return {a.x-b.x,a.y-b.y}; }
	Vector friend operator +(Vector a,Vector b){ return {a.x+b.x,a.y+b.y}; }
	Vector friend operator *(Vector a,db b){ return {a.x*b,a.y*b}; }
	db friend operator *(Vector a,Vector b){ return a.x*b.y-a.y*b.x; }
	db friend operator ^(Vector a,Vector b){ return a.x*b.x+a.y*b.y; }
}b[M];
struct Line{
	Vector pre,suf,v;db pos;
	friend bool operator <(Line x,Line y){ return dcmp(x.pos-y.pos)==0?dcmp(x.v*(y.suf-x.pre))>0:dcmp(x.pos-y.pos)<0;}
}a[M],que[M];
Vector cross(Line a,Line b){ return b.pre+b.v*(((b.pre-a.pre)*a.v)/(a.v*b.v)); }
bool pd(Line i,Line j,Line k){
	Vector p=cross(i,j);
	return dcmp(k.v*(p-k.pre))<0;
}

void half_plane(){
	sort(a+1,a+cnt+1);
	for(int i=1;i<=cnt;i++){
		if(dcmp(a[i].pos-a[i-1].pos)!=0) tot++;
		a[tot]=a[i];
	}
	h=1,t=0;que[++t]=a[1],que[++t]=a[2];
	for(int i=3;i<=tot;i++){
		while(h<t&&pd(que[t-1],que[t],a[i])) t--;
		while(h<t&&pd(que[h+1],que[h],a[i])) h++;
		que[++t]=a[i];
	}
	while(h<t&&pd(que[t-1],que[t],que[h])) t--;
	while(h<t&&pd(que[h+1],que[h],que[t])) h++;
	que[t+1]=que[h];cnt=0;
	for(int i=h;i<=t;i++) b[++cnt]=cross(que[i],que[i+1]);
}

void solve(){
	scanf("%d",&n);
	for(int i=1,mi;i<=n;i++){
		scanf("%d",&mi);
		for(int j=1;j<=mi;j++) scanf("%lf %lf",&b[j].x,&b[j].y);
		b[mi+1]=b[1];
		for(int j=1;j<=mi;j++) a[++cnt].pre=b[j],a[cnt].suf=b[j+1],a[cnt].v=b[j+1]-b[j];
	}
	for(int i=1;i<=cnt;i++) a[i].pos=atan2(a[i].v.y,a[i].v.x);
	half_plane();b[cnt+1]=b[1];
	if(cnt>2) for(int i=1;i<=cnt;i++) Ans+=b[i]*b[i+1];
	printf("%.3f\n",fabs(Ans)/2.0);
}

动态半平面交

给定一个平面,要求支持一下几种操作:

  1. 插入一个半平面

  2. 删除一个还在当前平面的半平面

  3. 询问一个点是否在当前半平面的交集内

因为最近有搞过这类的东西,就放上来了。

正常的半平面交是不能删除的,所以要依靠一个黑科技——线段树分治。

假如没有删除操作,可以使用 cdq 分治,按时间排序后,每次维护左区间的插入的半平面交的上凸包和下凸包,再依次判断右区间的点是否在左边的半平面交中,一旦在每层判断时有一次不在,答案就为否。

如果存在撤销操作,那一个半平面只会覆盖一段时间上的询问点,可以将询问点按时间做成线段树,插入和删除的半平面看做区间修改,像一般的线段树上的区间操作一样在节点上打上标记,每个节点开个 vector 存储半平面,然后遍历整棵线段树,对于一个节点,维护标记的半平面交的上凸包和下凸包,然后枚举对应区间的所有点进行判断。

由于线段树的性质,标记的个数的 \(\mathcal{O(n\log n)}\) 级别的,所有节点枚举对应区间的点的次数也是 \(\mathcal{O(n\log n)}\) 级别的,精细实现复杂度可以做到 \(\mathcal{O(n\log n)}\)

最小圆覆盖

APIO 讲课就听懂了这一部分,wtcl (用线性规划讲解计算几何好阴间啊)。

随机增量法是线性规划中的一个重要方法,主要是通过随机排序后,对加入的规约条件和已有的最优解进行讨论,修改最优解。

用在此处的主要思想是如果加入的点不在已有的最小圆内,那么该点一定在加入该点后的最小圆的边上。

分为三个部分,第一部分是开始无确定点,随机增量到一个不在最小圆的点时,考虑确定一个点如何求解最小圆,转到第二部分。

第二部分与第一部分类似,找到后变为确定两个点如何求最小圆,转到第三部分。

由于第三部分已经确定两个点了,如果找到一个不在之前最小圆上的点,三个点可以得出唯一的圆,可以直接求解。

然后可以发现期望复杂度是 \(\mathcal{O(n)}\) 的。

线性代数和计算几何的一些联系

最近搞了一下线性代数,感觉没啥用,后来去 THUSC 做学习题时发现在工程方面太有用了(当然对 OI 还是没啥用)!

有一个影响比较深的,是将一个在 \(n\) 维空间内的表达式转成 \(n\) 维线性无关的表达式,搞到线性代数里就是将一个矩阵转成对角矩阵,众所周知可以用拉普拉斯变换得出特征值后搞出基础解系再正交化,直接就机械化了!

THUSC 的学习题是光线追踪,感觉挺有意思的,就是我被 Linux 给杀了(

posted @ 2020-08-11 16:31  Dabuliuzp  阅读(635)  评论(0编辑  收藏  举报
/* */ 返回顶端