计算几何

计算几何初步 - 博客 - root的博客

前置知识

平面直角坐标系

我们的点和向量都是通过坐标来保存的。

处理精度:

1
2
3
4
5
6
7
8
const double eps = 1e-8;
int dcmp(double x)//处理精度，判断正负 
{
	if(fabs(x) < eps) return 0;
	if(x < 0) return -1;
	return 1;
}

向量及其运算

向量

向量的模

相反向量

垂直向量

共线向量

向量的加减

零向量

向量的数乘

向量的数量积

由向量的数量积可以判断它们的夹角和投影。

image-20240222150704008

向量的外积

向量的外积，也叫叉积，几何意义是两向量由平行四边形法则围成的面积。外积是一个向量，垂直于原来两个向量所在的平面。

外积没有交换律。

image-20240222150759176

由向量外积可以判断两向量的旋转关系，同时可以方便地求出点到直线的距离。分别会用于下面的凸包和旋转卡壳。

向量旋转

image-20240222151744355
image-20240222151852560

线性代数的表示方法

image-20240222151954577

三角剖分求面积

三角剖分可以用来求任意多边形的面积。

我们只需要知道一个多边形的各个连续顶点分别在哪里即可。

image-20240222152927361
image-20240222152935967
image-20240222153003591

线的记录

直线与射线

我们记录的是：直线上一点和直线的方向向量。

线段

线段很好记录：只需要记录左右端点即可。

多边形

开数组按一定顺序记录多边形的每个顶点即可。

特殊地，如果矩形的各边均与某坐标轴平行的话，我们只记录左下角和右上角的顶点即可。

曲线

一些特殊曲线，如函数图像等一般记录其解析式。对于圆，直接记录其圆心和半径即可。

基本公式

正弦定理

余弦定理

略

凸包定义

在平面上给定若干点，其 凸包 被定义为能包含所有点的 凸多边形 的 交集。

性质

  1. 在所有包含所有点的凸多边形中：
    • 凸包的面积是最小的（显然）。
    • 凸包的周长也是最小的。
      • 因为每个包含所有点的凸多边形，总可以不断缩小周长，最终得到凸包。
  2. 凸包的每个顶点都是原有点集中的点。

直观理解

在二维欧几里得空间中，凸包可以想象为一条刚好包着所有点的橡皮圈。简单来说，给定二维平面上的点集，凸包就是将最外层的点连接起来构成的凸多边形，它能包含点集中所有的点。

凸包问题

问题描述：给定点集，求构成凸包的点。

求解方法

求凸包的方法有很多，这里介绍两种经典算法：

  1. Graham 扫描法
  2. Andrew 算法

Graham扫描法

本质：

Graham扫描算法维护一个凸壳 通过不断在凸壳中加入新的点和去除影响凸性的点 最后形成凸包

开始操作：

  • 我们还是选择一个y值最小（如有相同选x最小）的点，记为P1
  • 剩下的点集中按照极角的大小逆时针排序，然后编号为P2~Pm
  • 我们按照排序结束时的顺序枚举每一个点，依次连线，这里可以使用一个栈来存储，每次入栈，如果即将入栈的元素与栈顶两个元素所构成了一个类似于凹壳的东西，那么显然处于顶点的那个点一定不在这个点集的凸包上，而他正好在栈顶，所以把它弹出栈，新点入栈。

但是，新来的点有可能既踢走了栈顶，再连接新的栈顶元素后却发现仍然可以踢出，此时就不能忘记判断。

img

由于第一篇题解的讲解实在太好，这里就不要脸的贴上（）

ShineEternal （大佬洛谷名字）

其中p0为起始点

然后p1准备入栈，由于栈中元素过少，所以检验合格，可直接进入。

之后因为p2仍为向左转，符合凸包条件，所以暂时先让它进去

p3出现了右转现象，那么我们就把顶点p2舍去，在检查p3的性质，合格

就是大致这么个过程

动图

在这里插入图片描述

由此，我们的Graham算法的全过程就结束了。

扫描的时间复杂度：O(n)

但是显然不可能做到这么优秀. 于是还有排序的时间复杂度： O(nlogn)

合起来总的时间复杂度： O(nlogn)

[USACO5.1] 圈奶牛Fencing the Cows /【模板】二维凸包

题目背景

upd: 新增一组 hack 数据。

题目描述

农夫约翰想要建造一个围栏用来围住他的奶牛，可是他资金匮乏。他建造的围栏必须包括他的奶牛喜欢吃草的所有地点。对于给出的这些地点的坐标，计算最短的能够围住这些点的围栏的长度。

输入格式

输入数据的第一行是一个整数。表示农夫约翰想要围住的放牧点的数目 \(n\)。

第 \(2\) 到第 \((n + 1)\) 行，每行两个实数，第 \((i + 1)\) 行的实数 \(x_i, y_i\) 分别代表第 \(i\) 个放牧点的横纵坐标。

输出格式

输出输出一行一个四舍五入保留两位小数的实数，代表围栏的长度。

样例 #1

样例输入 #1

1
2
3
4
5
4
4 8
4 12
5 9.3
7 8

样例输出 #1

1
12.00

提示

数据规模与约定

对于 \(100\%\) 的数据，保证 \(1 \leq n \leq 10^5\)，\(-10^6 \leq x_i, y_i \leq 10^6\)。小数点后最多有 \(2\) 位数字。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
#include<bits/stdc++.h>
using namespace std;
int n;
const int N=100005;
struct ben
{
    double x,y;
}p[N],s[N];
double check(ben a1,ben a2,ben b1,ben b2)//检查叉积是否大于0，如果是a就逆时针转到b 
{
    return (a2.x-a1.x)*(b2.y-b1.y)-(b2.x-b1.x)*(a2.y-a1.y);
}
double d(ben p1,ben p2)//两点间距离 
{
    return sqrt((p2.y-p1.y)*(p2.y-p1.y)+(p2.x-p1.x)*(p2.x-p1.x));
}
bool cmp(ben p1,ben p2)//排序函数，这个函数别写错了，要不然功亏一篑 
{
    double tmp=check(p[1],p1,p[1],p2);
    //检测一下叉积
    if(tmp>0) 
		return 1;
    if(tmp==0&&d(p[0],p1)<=d(p[0],p2)) 
		return 1;
    return 0;
}
int main()
{
    cin>>n;
    double mid;
    for(int i=1;i<=n;i++)
    {
       	cin>>p[i].x>>p[i].y;
       if(i!=1&&(p[i].y<p[1].y||(p[i].y==p[1].y&&p[i].x<p[1].x)))//这是是去重 
        {
            mid=p[1].y;p[1].y=p[i].y;p[i].y=mid;
            mid=p[1].x;p[1].x=p[i].x;p[i].x=mid;
        }
        //这里保证p1的值一定是最小的，方便排序
       //别看这么长，只是两个交换x,y
    } 
    sort(p+2,p+1+n,cmp);//系统快排 
    s[1]=p[1];//最低点一定在凸包里 
    int cnt=1;
    //生成凸包过程
    for(int i=2;i<=n;i++)
    {
        while(cnt>1&&check(s[cnt-1],s[cnt],s[cnt],p[i])<=0) 
            //判断前面的会不会被踢走，如果被踢走那么出栈，//没有向左转的意思啦
			cnt--;
        cnt++;
        s[cnt]=p[i];
    }
    s[cnt+1]=p[1];
    //最后一个点回到凸包起点 
    //答案处理过程
    double ans=0; 
    for(int i=1;i<=cnt;i++) 
		ans+=d(s[i],s[i+1]);//然后s里存好了凸包序列，只需要把两两距离累加就行
 printf("%.2lf\n",ans);
}

旋转卡壳

凸包的直径——旋转卡壳 - 小蒟蒻yyb - 博客园 (cnblogs.com)

我们很显然可以得出，平面内最远的点对一定在凸包上面 而凸包的直径也就是凸包上最远点对的距离

显然不一定所有点都会在凸包上，显然比O(n^2)的枚举的效率会更加优秀（但是求解凸包还有一个O(nlogn)）

假设我现在考虑的是AB边和他距离最远的点

image-20240222162900794

很显然，D点距离他最远

换条边，换成BC边

image-20240222212516193

这个时候变成了E点

再继续？ 看看CD边

image-20240222212532022

我们的边依次是： AB->BC->CD 而点依次是：D->E->A

发生了一件很神奇的事情 我们的边是逆时针的枚举 而此时距离他们最远的点也是逆时针的出现！！

那么，这个时候，利用这个性质， 我们只需要逆时针枚举边，然后从上次枚举到的最远点继续逆时针向后枚举就能求解。

现在我们需要考虑，如何得到距离每条对应边的的最远点呢?稍加分析,我们可以发现,我们可以把点到直线的距离化解为三角形的面积，再运用叉积进行比较。同时，凸包上的点依次与对应边产生的距离成单峰函数,面积上升到最高点后，又会下降。(具体证明可以从凸包定义入手 用反证法解决)

利用单峰函数这一性质，我们可以根据凸包上点的顺序，枚举对踵点，直到下一个点的距离小于当前点就可以停止了，而且随着对应边的旋转，最远点也只会顺着这个方向旋转，我们可以从上一次的对踵点开始继续寻找这一次的。由于内层while循环的执行次数取决于j增加次数 j最多增加O(n)次，所以求出所有对踵点的时间复杂度为O(n)。

xz6.png

【模板】旋转卡壳 | [USACO03FALL] Beauty Contest G

题目描述

给定平面上 \(n\) 个点，求凸包直径。

输入格式

第一行一个正整数 \(n\)。
接下来 \(n\) 行，每行两个整数 \(x,y\)，表示一个点的坐标。保证所有点的坐标两两不同。

输出格式

输出一行一个整数，表示答案的平方。

样例 #1

样例输入 #1

1
2
3
4
5
4
0 0
0 1
1 1
1 0

样例输出 #1

1
2

提示

【数据范围】
对于 \(100\%\) 的数据，\(2\le n \le 50000\)，\(|x|,|y| \le 10^4\)。

\(\text{upd 2022.7.22}\)：新增加四组 Hack 数据。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
#include<bits/stdc++.h>
using namespace std;
#define MAXN 50005

struct Point
{
    int x, y;
    bool operator < (const Point& _P) const
    {
        return y<_P.y||(y==_P.y&&x<_P.x);
    };
}pset[MAXN],ch[MAXN];


int cross(Point a,Point b,Point o)      
{
    return (a.x - o.x) * (b.y - o.y) - (b.x - o.x) * (a.y - o.y);
}

void convex_hull(Point *p,Point *ch,int n,int &len)
{
    sort(p, p+n);
    ch[0]=p[0];
    ch[1]=p[1];
    int top=1;
    for(int i=2;i<n;i++)
    {
        while(top>0&&cross(ch[top],p[i],ch[top-1])<=0)
            top--;
        ch[++top]=p[i];
    }
    int tmp=top;
    for(int i=n-2;i>=0;i--)
    {
        while(top>tmp&&cross(ch[top],p[i],ch[top-1])<=0)
            top--;
        ch[++top]=p[i];
    }
    len=top;    
}


int dist2(Point a,Point b)
{
    return (a.x - b.x) * (a.x - b.x) + (a.y - b.y) * (a.y - b.y);
}

int rotating_calipers(Point *ch,int n)
{
    int q=1,ans=0;
    ch[n]=ch[0];
    for(int p=0;p<n;p++)
    {
        while(cross(ch[p+1],ch[q+1],ch[p])>cross(ch[p+1],ch[q],ch[p]))
            q=(q+1)%n;
        ans=max(ans,max(dist2(ch[p],ch[q]),dist2(ch[p+1],ch[q+1])));            
    }
    return ans; 
}

int main()
{
    int n, len;
		cin>>n;
        for(int i = 0;i < n;i++)
        {
           cin>>pset[i].x>>pset[i].y;
        }
        convex_hull(pset,ch,n,len);
      cout<<rotating_calipers(ch,len);
 
    return 0;
}