计算几何
非常幼稚和不严谨的学习笔记。
前置
板板
const double eps=1e-8;
struct node{
double x,y;node(double xx=0,double yy=0){x=xx,y=yy;}
inline void in(){scanf("%lf%lf",&x,&y);}
inline void print(char wA,char wB){printf("%.4f%c%.4f%c",x,wA,y,wB);}
inline node operator +(node s1){return (node){x+s1.x,y+s1.y};}
inline node operator -(node s1){return (node){x-s1.x,y-s1.y};}
inline node operator *(double s1){return (node){x*s1,y*s1};}
inline node operator /(double s1){return (node){x/s1,y/s1};}
inline bool operator ==(node s1){return x==s1.x&&y==s1.y;}
inline double len(){return sqrt(x*x+y*y);}
};
inline int sgn(double x){return x<-eps?-1:(x>eps?1:0);}
inline double Abs(double x){return x*sgn(x);}
点积
点积比较简单，即 \(a\cdot b=|a||b|\cos\theta\)，计算方式是 \(\text{dot}=x_1x_2+y_1y_2\)。然后常用的是用点积除以其中一个的模长，得到的就是另一个向量到这个向量的投影长度。根据这个投影方向和原向量的方向这玩意有正有负。这玩意显然是符合交换律的，毕竟 \(\cos\) 的符号圆也是这样的嘛。然后根据 \(\cos\) 的符号性质可以通过点积来判断 \(\theta\) 和 \(90^\circ\) 的关系，如果正的话就小于，否则大于，刚好是 \(90^\circ\) 就刚好是 \(0\)。
inline double dot(node x,node y){
return x.x*y.x+x.y*y.y;
}
叉积
然后是叉积。\(a\times b=|a||b|\sin\theta\)，计算方式是 \(\text{dog}=x_1y_2-y_1x_2\)。数值上等于从 \(a\) 到 \(b\) 形成的平行四边形的面积，由于有方向（也可以参照 \(\sin\) 的符号圆），得出来的数值的符号取决于两个向量的相对位置，如果 \(a\) 到 \(b\) 是逆时针更近那么就是正的，否则就是负的。当然如果叉积是 \(0\) 的话说明两个向量是平行的。
inline double dog(node x,node y){
return x.x*y.y-x.y*y.x;
}
inline bool onRight(node p,node a,node b){
return sgn(dog(a-p,b-p))==-1;
}
点
旋转
点 \((x,y)\) 绕着原点顺时针旋转 \(\theta\) 后的点坐标应该是 \((x\cdot\cos\theta+y\cdot\sin\theta,-x\cdot\sin\theta+y\cdot\cos\theta)\)。如果不是以原点为中心的话直接重新建系还原回去即可。一定要记得是顺时针。
inline node turn(node x,double y){
double xx=x.x*cos(y)+x.y*sin(y);
double yy=-x.x*sin(y)+x.y*cos(y);
return (node){xx,yy};
}
inline node turn(node x,node y,double q){return turn(x-y,q)+y;}
信用卡凸包。首先发现带圆弧的凸包可以分成两个部分，一个部分可以直接投影到圆心形成的凸包上，另一种则是圆弧，而且所有圆弧加起来刚好是一整个圆。所以要做的就是整理出圆心的点集，对这个点集求凸包之后求出凸包的长度加上圆的长度即可。
相似的还有 UVA10652。做这种题的时候一定要注意的是旋转方向，毕竟公式和代码里给的都是顺时针旋转某个弧度的结果。而且要注意给的是弧度还是角度。理论上来说 sin 和 cos 都是定义在实数上的函数，所以不用特别注意旋转角的范围和正负问题。
点是否在线段上
可以直接用 \(dis_{x,p}+dis_{p,y}\) 与 \(dis_{x,y}\) 做比较，但是精度可能会损失比较多。另一个方法是直接看 \(\frac{a\cdot b}{|a||b|}\) 和 \(-1\) 的关系，但同样会损失精度。总之呢就是尽量不要去算线段长度（毕竟 \(\text{sqrt}\) 函数会大大地损失精度）。所以通过判断共线且外开的方式来判断，即 \(a\times b=0\) 并且 \(a\cdot b<0\)。之前的写法貌似在这个点是线段端点的时候会锅，所以稍微特判一下吧。
inline bool onLineR(node p,node a,node b){
if(p==a||p==b)return true;
return sgn(dog(a-p,b-p))==0&&sgn(dot(a-p,b-p))==-1;
}
点到线段的距离
首先算出来 \(p\) 分别到 \(x\) 和 \(y\) 的距离，然后看两个底角有没有钝角，就是看 \((p,x)\cdot (x,y)\) 和 \(0\) 的关系。如果有小于 \(0\)，也就是有钝角的话那么就取前面的那个，否则直接返回三角形的高，即 \(|\frac{(p,x)\times(p,y)}{\text{dis}(x,y)}|\)。然后似乎要特判 \(x,y\) 重合的情况，但是什么样的傻逼才会卡这种细节啊喂。
inline double workDis(node p,node a,node b){
if(a==b)return (a-p).len();
node x=p-a,y=p-b,xx=b-a,yy=a-b;
if(sgn(dot(x,xx))==-1)return (p-a).len();
if(sgn(dot(y,yy))==-1)return (p-b).len();
return (p-findFoot(p,a,b)).len();
}
点是否在直线上
和线段上同理，不需要判断内张还是外张，充要条件就是 \(a\times b=0\)。
inline bool onLine(node p,node a,node b){
return sgn(dog(a-p,b-p))==0;
}
点到直线的垂足
先求出 \(-a\) 在 \(l\) 上的投影 \(a'\)，同理求出 \(b'\)，发现垂足到 \(x\) 的距离是 \(l\) 长度的 \(\frac{a'}{a'-b'}\) 倍，计算即可。别忘了计算投影长度的时候要把地面除掉。
inline node findFoot(node p,node a,node b){
node x=p-a,y=p-b,xx=b-a,yy=a-b;
double lenA=dot(x,xx)/xx.len(),lenB=dot(y,yy)/yy.len();
return a+xx*(lenA/(lenA+lenB));
}
点关于直线的对称点
关于垂足对称。使用初中学的什么中点坐标公式反求一下即可。
inline node findRef(node p,node a,node b){
return findFoot(p,a,b)*2.0-p;
}
例题有 折纸，由于数据规模比较小只需要暴力模拟即可。会发现一个点打下去是有纸的当且仅当在整个过程中它都在折痕的左边，按照这个标准衍生这个点即可。主要只有在最后一步才检验这个点是否在正方形范围内，因为在折纸的过程中有可能会出现正方形纸外的合法点。
点是否在多边形以内
有非常经典的想法是说从这个点引一条向右的水平线看和多边形有多少交点。为了防止一些不必要事情的发生，我们认为一个线段在 \(y\) 轴上的映射是上闭下开的（也就是说完全水平的线会被忽略），判断的时候注意一下符号。然后就是求水平线和线段的交点，传说中的相似三角形即可。最后就是看交点的横坐标和要求的点的横坐标的大小关系即可，累加之后看其奇偶性就可以知道这玩意在不在多边形之中。
inline void solve(node p){
int nowNum=0;
for(int i=1;i<=m;i++){
node x=a[i],y=a[i%m+1];
if(onLine(x,y,p))return puts("YES"),void();
if(x.y>y.y)swap(x,y);
if(sgn(p.y-y.y)==1||sgn(x.y-p.y)>=0)continue;
double xx=(x+(y-x)*((p.y-x.y)/(y.y-x.y))).x;
nowNum+=sgn(xx-p.x)==1;
}
if(nowNum&1)return puts("YES"),void();
else return puts("NO"),void();
}
平面最近点对
好像混进来了一个奇怪的东西。平面最近点对利用了分治的思想，先计算左右半平面最近的点对，然后计算跨越中间线的点对，而跨越中间线的点对不会太多，因为对于一个点而言，河对岸和它靠得比较近的点不会太多，因为既然那么多点和它都靠得那么近，那么它们之间也会靠得很近；既然它们之间都靠得那么近，那么先期的答案就应该更小。然后就矛盾了。复杂度是有保障的线性对数，过程中要同时归并排序。
int solve(int l,int r){
if(l==r)return inf;int mid=l+r>>1;int xx=a[mid].x;
int d=min(solve(l,mid),solve(mid+1,r));int tl=l,tr=mid+1,tt=l;
while(tl<=mid&&tr<=r)
if(a[tl].y<a[tr].y)b[tt++]=a[tl++];
else b[tt++]=a[tr++];
while(tl<=mid)b[tt++]=a[tl++];while(tr<=r)b[tt++]=a[tr++];
for(int i=l;i<=r;i++)a[i]=b[i];tt=0;
for(int i=l;i<=r;i++)if(qq(a[i].x-xx)<=d)b[++tt]=a[i];
for(int i=1;i<=tt;i++)
for(int j=i-1;j;j--){if(qq(b[i].y-b[j].y)>d)break;d=min(d,dis(b[i],b[j]));}
return d;
}
板板：P7883，UVA10245。话说 UVA 的题怎么都那么骨骼清奇，输入格式又傻逼，输出格式又傻逼。
线&多边形&圆
两条直线的交点
仍然使用面积之比来做。画个图就可以大概理解到它求交点的原理。
inline node findPt(node a,node b,node c,node d){
node x=b-a,y=d-c,z=a-c;
return a+x*(dog(y,z)/dog(x,y));
}
例题 切割多边形 由于数据范围很小，只需要暴力枚举或者用一点状压 DP 来做就可以了。有一些细节要注意一下。
线段是否有交点
求直线的交点，然后看交点是否在两个线段上。
inline node havePt(node a,node b,node c,node d){
node an=findPt(a,b,c,d);
return onLineR(an,a,b)&&onLineR(an,c,d);
}
一个点是否在圆上
看这个点和圆心的距离是否小于等于圆的半径即可。
inline bool inCir(node p,Cir c){
return sgn(c.r-(p-c.p).len())>=0;
}
三点确定一个圆
用小学的方法来搞，做三角形两条边的中垂线，然后求两个中垂线的交点，这个交点就是圆的圆心。半径直接求距离即可。
struct Cir{node p;double r;};
inline node turnPl(node a){
return (node){-a.y,a.x};
}
inline Cir findCir(node a,node b,node c){
node p1=(a+b)/2,p2=(a+c)/2;
node p=findPt(p1,p1+turnPl(a-b),p2,p2+turnPl(a-c));
Cir an;an.p=p,an.r=(a-p).len();
return an;
}
板板：UVA190。UVA 傻逼输出格式的活化石和集大成者。
然后是经典老番 CF1C，有谁想到 CF 第一次比赛就开始搞计算几何了。
最小圆覆盖
比较强的想法。我用的随机增量法，基本思路是说假如 \(p\) 不在 \(S\) 的最小覆盖圆上，那么它一定会在 \(S\cup p\) 的最小覆盖圆上。基本流程就是对于当前的点，如果在圆上就直接忽略，否则找前面的不在圆上的点，并以这两个点为基础去尝试包括更多的点。经过一些妙的分析可以得出这个算法是线性的。然后就是在计算之前要把点集随机打乱一下，防止毒瘤出题人卡你。
Cir ans=(Cir){(node){-inf,-inf},0};
for(int i=1;i<=m;i++){
if(inCir(a[i],ans))continue;
ans=(Cir){a[i],0};
for(int j=1;j<i;j++){
if(inCir(a[j],ans))continue;
ans=(Cir){(a[i]+a[j])/2,(a[i]-a[j]).len()/2};
for(int k=1;k<j;k++){
if(inCir(a[k],ans))continue;
ans=findCir(a[i],a[j],a[k]);
}
}
}
板板：最小圆覆盖。然后有个 双倍经验（信号塔）。都比较好写。
凸包
求凸包
有个非常自然的想法，对凸包上半部分和下半部分分别求。显然排序之后第一个点和最后一个点都应该在凸包之中，所以可以先把这两个点定下来，然后从前往后扫所有点，同时维护一个栈，通过弹栈的方式来动态地维护一个凸包，求即可。这种方法可以控制求出来的凸包的边上有没有点。下半部分的代码：
sort(a+1,a+m+1);st[top=1]=1;
for(int i=2;i<=m;i++){
node x=a[i];
while(top>1){
node y=a[st[top]],z=a[st[top-1]];
if(sgn(dog(x-y,z-y))==1)break;top--;
}
st[++top]=i;
}
P2742 板子，不多说。
保护出题人：可以画出一个阶梯状的图，向右爬升的那种。每个台阶的宽都是 \(d\)。发现一个僵尸吃掉植物当且仅当由这个点引一条斜率为攻击速度的直线，与 \(y\) 轴的交点在原点下方。而每个僵尸对应的点是不变的，可以对于每轮游戏重新建系。具象化地，在第 \(i\) 轮游戏之中，第 \(k\) 个点的距离是 \(x_i+(k-1)d\)。而每个点需要消耗的血量实际上是 \(s_i-s_{k-1}\)，所以每一轮的答案会是 \(\max\limits_{i=1}^x\frac{s_x-s_{i-1}}{x_x+(x-i)d}\)，直接暴力是 \(O(N^2)\) 的，考虑优化。显然可以把柿子看成是点 \((x_x+xd,s_x)\) 和 \((id,s_{i-1})\) 连线的斜率。后者是不变的，所以可以对后面那堆东西维护一个凸包，而且显然最后出来的是一个下凸包，最后二分即可。我二分写挂了，所以写的三分，三分还是很显然的。而且跑得也不算太慢。
P2924 发现一些边能形成向量当且仅当它们首尾顺次连接并且向量对应的弧度值是单调的。所以可以把所有边求出来，按 atan2 值排序，然后枚举凸包上的第一个点，DP 即可。大概形式是顺次枚举边，用入点去更新出点。复杂度 \(O(N^3)\)。
动态凸包
按照思路，对上凸包和下凸包分别用 set 维护点集。判断一个点是否在凸包里等价于这个点在下凸包的上面、上凸包的下面。所以只需要在 set 中找到第一个横坐标大于等于当前点的点，特判这个迭代器是 begin 和 end 的情况。然后就只需要看当前点是否在这个点和上一个点形成的向量的左边即可。插入部分，首先判断当前的点是否在凸包以内，如果在的话直接丢掉；如果不在的话说明它会在新凸包上，也就是说旧凸包它两边的点就有可能被删除。相当于是遍历当前点左右两边的点，如果当前遍历到的点和左右两边形成了一个凹陷的话就删掉当前的点并遍历下一个点，否则就退出遍历过程。细节比较多，涉及很多迭代器相关的内容，看代码。
struct tank{
set<node>s;int op;
//up:1 down:-1
inline bool inn(node p){
auto it=s.lower_bound((node){p.x,-inf});
if(it==s.end())return false;
if(sgn((*it).x-p.x)==0)return (p.y-it->y)*op<=0;
if(it==s.begin())return false;
auto j=*prev(it),k=*it;return dog(p-j,k-j)*op>=0;
}
inline int judge(IT it){
IT j=it,k=it;
if(j==s.begin())return 0;--j;
if(++k==s.end())return 0;
node it_=*it,j_=*j,k_=*k;
return dog(it_-j_,k_-j_)*op>=0;
}
inline void insert(node P){
if(inn(P))return;
IT tmp=s.lower_bound((node){P.x,-inf});if(tmp!=s.end()&&tmp->x==P.x)s.erase(tmp);
s.insert(P);IT it=s.find(P),p=it;
if(p!=s.begin()){--p;while(judge(p)){IT tmp=p--;s.erase(tmp);}}
if((p=++it)!=s.end())while(judge(p)){IT tmp=p++;s.erase(tmp);}
}
}up,down;
半平面交
最后交出来的一定是一个凸包，这很显然。把所有向量按照极角从大到小排序，极角相同的向量则贪心选择覆盖平面更大的哪一个。维护一个双端队列，表示当前被选择的向量集合。最后一个向量应该只贡献了一个点，如果这个点在当前向量的右边的话说明这个向量就是没用的，就应该丢出去。同理，队首也应该通过相同的方式丢掉。完成以上两步之后再把当前的向量压入队尾。在处理完所有向量之后再检查一下剩下来的向量，如果按照上面说的相同的方式弹掉队首和队尾，剩下的向量挨个求交点就是最后半平面交的凸包了。有一些需要注意的：如果题目没有保证交出来面积有限的话可以搞一个大框框框起来。然后就是检查的时候要先检查队尾再检查队首。然后就是 atan2 函数先传纵坐标再传横坐标（老生常谈啦）。
inline double getAng(Vec a){
return atan2((a.b-a.a).y,(a.b-a.a).x);
}
inline bool cmp(Vec a,Vec b){
double pa=getAng(a),pb=getAng(b);
if(sgn(pa-pb)!=0)return pa<pb;
return onRight(a.a,a.b,b.a);
}
inline bool check(Vec a,Vec b,Vec p){
node pp=getPt(a,b);
return onRight(p.a,p.b,pp);
}
signed main(){
sort(c+1,c+n+1,cmp);
for(int i=1;i<=n;i++){
if(i!=1&&sgn(getAng(c[i])-getAng(c[i-1]))==0)continue;
while(r-l>=1&&check(q[r-1],q[r],c[i]))r--;
while(r-l>=1&&check(q[l+1],q[l],c[i]))l++;
q[++r]=c[i];
}
while(r-l>=1&&check(q[r],q[r-1],q[l]))r--;
while(r-l>=1&&check(q[l],q[l+1],q[r]))l++;
s[num=1]=getPt(q[l],q[r]);double ans=0;
for(int i=l;i<r;i++)s[++num]=getPt(q[i],q[i+1]);
return 0;
}
板板：P4196。
最小矩形覆盖
显然出来的矩形和凸包有至少一条边重合。枚举这条边，然后用旋转卡壳相似的方式求出最上面的那个点，用相似的方法求出最左边和最右边的点。我求最左边和最右边点求复杂了，我用的什么投影到底线离点最远的方法，当投影点回去的时候（比如类似三角形的情况）会出问题，所以还需要判断一下当前的点是否在射线上（如果是的话就丢掉，射线是开的），就很智障。题解有方法，只需要对底线和当前直线的夹角和九十度比较即可，显然好写很多很多。
