計算幾何基礎
Tags：高級算法
前言
Noip爆炸后差點退組。
不取模設小上界的人大概不配學計算幾何吧。
copy一些網上的博客好了。
計算幾何詳解：https://blog.csdn.net/clover_hxy/article/details/53966405
凸包詳解：http://www.cnblogs.com/aiguona/p/7232243.html
旋轉卡殼：https://blog.csdn.net/wang_heng199/article/details/74477738
題目標簽：WFJ
聲明一下：一下代碼片段現已統一格式：*表示叉積，&表示點積，^表示數乘。符合博主現在的碼風。當然對於實數並不存在運算符混淆問題。
向量線段相關
點積
\(a.x*b.x+a.y*b.y\)
叉積
\(a.x*b.y-a.y*b.x\)
借圖：左邊是點積，右邊是叉積
轉角公式
逆時針旋轉角度\(B\)
原：\((cosAr,sinAr)\)
現：\((cos(A+B)r,sin(A+B)r)=((cosAcosB-sinAsinB)r,(sinAcosB+cosAsinB)r)\)
即：\((x,y)->(xcosB-ysinB,xsinB+ycosB)\)
極角排序
兩種方法
\(atan2\)
這是一個給定三角形對邊和鄰邊算角的函數，在cmath庫中
對每個向量求\(angle=atan2(a.y,a.x)\)后直接排序即可
舉例：\(atan2(1,1)=0.78,atan2(1,0)=1.57\)
叉積
叉積的正負可以判斷左右關系，所以在最大跨角不超過180度的時候可以非常好地實現，否則有可能隨機一個起點繞好幾圈
兩直線夾角
若兩個向量都是從原點出發，則可以很方便地用$$\theta =atan2(x*y,x&y)$$叉積是平行四邊形，底乘高；點積是底乘投影。
更特殊點，該夾角是有向夾角，從x向量旋轉到y，逆時針是正角、順時針是負角。
兩直線交點
面積比->線段比->定比分點
Node Crosss(Node a1,Node a2,Node b1,Node b2)
{
Node a=a2-a1,b=b2-b1,c=b1-a1;
if(fabs(b*a)<1e-12) return (Node){-1e9,-1e9};//平行
double t=(b*c)/(b*a);
return a1+(a^t);
}
給個圖形象一點（之前畫錯了orzFlashHu）
判斷兩線段是否相交
\(a1,a2,b1,b2\)為兩條直線的兩個端點
如果滿足\((b1-a1)×(b1-a2)\)和\((b2-a1)×(b2-a2)\)不同號，並且\((a1-b1)×(a1-b2)\)和\((a2-b1)×(a2-b2)\)不同號，則兩線段相交，叫跨立實驗。
多邊形相關
凸包
找到最下的點設為原點，將剩下的點用叉積極角排序
用單調棧維護凸包，當\((A[i]-A[sta[top-1]])×(A[sta[top]]-A[sta[top-1]])\)為非負時需要彈棧
int cmp1(const Node&A,const Node&B) {return A.y<B.y||(A.y==B.y&&A.x<B.x);}
int cmp2(const Node&A,const Node&B) {return A*B>0||(A*B==0&&A.dis()<B.dis());}
//這個表示如果兩向量共線，把短的放前面（方便被踢開）
void Tubao()
{
sort(A+1,A+n+1,cmp1);//找到最底下的點
for(int i=n;i>=1;i--) A[i]=A[i]-A[1];
sort(A+2,A+n+1,cmp2);//極角排序
for(int i=1;i<=n;sta[++top]=i,i++)
while(top>=2&&(A[i]-A[sta[top-1]])*(A[sta[top]]-A[sta[top-1]])>=0) top--;
n=top;for(int i=1;i<=top;i++) A[i]=A[sta[i]];
}
判斷點是否在多邊形內
從該點向右引射線，與多邊形交奇數次即在其內。
順時針+1，逆時針-1，注意與頂點重合的情況
\(Practice：Codeforces375C\ Circling\ Round\ Treasures\)
判斷點是否在凸包內
先判定和邊界的關系
然后找到與其極角相鄰的兩點，憑此判斷
須保證\(A[1]=(0,0)\)
ll in(Node a)
{
if(a*A[1]>0||A[tot]*a>0) return 0;
ll ps=lower_bound(A+1,A+tot+1,a,cmp2)-A-1;
return (a-A[ps])*(A[ps%tot+1]-A[ps])<=0;
}
求多邊形面積
逆時針極角排序后直接用叉積求
double CalcS()
{
double res=0;
for(int i=2;i<node;i++) res+=(A[i]-A[1])*(A[i+1]-A[1]);
return res/2;
}
三角形重心
三點各坐標的平均值\((\frac{x1+x2+x3}{3},\frac{y1+y2+y3}{3})\)
算法
半平面交
用來求解一類線性規划問題？
首先存一個極大的凸多邊形，考慮增加一條直線
那么將其視為一個向量，左邊的為合法半平面
復雜度\(O(n^2)\)
void Cut(Node a,Node b)//增加一個向量
{
A[tt+1]=A[1];c=0;
for(int i=1;i<=tt;i++)//枚舉每一條邊
{
double v1=(A[i]-a)*(A[i]-b);
double v2=(A[i+1]-a)*(A[i+1]-b);
if(v1>=0) C[++c]=A[i];//如果A[i]在合法平面內則存入
if(v1*v2<0) Add(A[i],A[i+1],a,b);//如果邊(A[i],A[i+1])穿過新增直線，把交點加入
}
tt=c;for(int i=1;i<=tt;i++) A[i]=C[i];
}
如果不是凸多邊形的時候，可以采用水平可見直線的方法，先按照斜率排序后直接用棧維護下平面。復雜度為排序復雜度，有一定局限性，但在坐標范圍較大采用上面方法會炸精度時非常好用
struct Line{int k,b;int id;}A[N];
int cmp1(const Line&A,const Line&B) {return A.k<B.k||(A.k==B.k&&A.b>B.b);}
double Node(int i,int j) {return 1.0*(A[j].b-A[i].b)/(A[i].k-A[j].k);}
void work()
{
sort(A+1,A+n+1,cmp1);
for(int i=1;i<=n;sta[++top]=i,i++)
while(top>=2&&Node(i,sta[top])-Node(sta[top],sta[top-1])<=eps) top--;
}
旋轉卡殼
利用凸包上一些奇妙的單調性，求解
- 多邊形直徑
- 多邊形寬度
- 最小面積矩形覆蓋
- 等等
復雜度\(O(n)\)
for(int i=1,p=1;i<=n;i++)
{
//這份代碼只卡了一個點，還可以同時卡多個點
while(H(A[i],A[i+1],A[p%n+1])>=H(A[i],A[i+1],A[p])) p=p%n+1;
Ans=max(Ans,max((A[p]-A[i]).dis(),(A[p]-A[i+1]).dis()));
}
最小圓覆蓋
隨機增量構造，復雜度\(O(n)\)
感覺可以被扁平三角形卡掉，但是實際上並不會，由於\(ijk\)的枚舉順序，三角形一定是在嘗試把三條邊都作為直徑都失敗后、才會有三角形的外接圓。（2019.2.23 byGXY）
pair<Node,Node> Line(int k,int i)//中垂線
{
Node a=A[k]-A[i];a.rot();
Node c=(A[k]+A[i])^0.5;
return mp(c,a+c);
}
int main()
{
random_shuffle(A+1,A+n+1);//先隨機化
for(int i=1,j,k;i<=n;i++)
if((A[i]-C).dis()>R)//某點不在圓內，更新為半徑為0的圓
for(C=A[i],R=0,j=1;j<i;j++)
if((A[j]-C).dis()>R)//某線段不在圓內，更新為以該線段為直線的圓
for(C=(A[i]+A[j])^0.5,R=(A[i]-C).dis(),k=1;k<j;k++)
if((A[k]-C).dis()>R)//某三角形不在圓內，更新為該三角形的外切圓
C=Cross(Line(k,i),Line(k,j)),R=(A[i]-C).dis();//求兩條中垂線交點
}
自適應辛普森積分
用來求定積分。
當函數極其鬼畜、退不出原函數的時候，采用二次函數擬合，在精度達到要求之后便視為相等
即$$\int_L^Rf(x)dx\approx\int_L^R(Ax^2+Bx+C)dx=\frac{(R-L)(f(L)+f(R)+4f(\frac{L+R}{2}))}{6}$$
double f(double x) {return (c*x+d)/(a*x+b);}
double calc(double l,double r) {return (r-l)*(f(l)+f(r)+4*f((l+r)/2))/6;}
double simpson(double l,double r,double ans)
{
double mid=(l+r)/2,ans1=calc(l,mid),ans2=calc(mid,r);
if(fabs(ans1+ans2-ans)<=eps) return ans;
return simpson(l,mid,ans1)+simpson(mid,r,ans2);
}
int main()
{
std::cin>>a>>b>>c>>d>>L>>R;
printf("%.6f\n",simpson(L,R,calc(L,R)));
}
注意\(\frac{1}{x}\)的原函數是\(ln|x|\)！！
三維凸包
出門左拐：https://www.cnblogs.com/xzyxzy/p/10225804.html
閔可夫斯基和
出門左拐第二家：https://www.cnblogs.com/xzyxzy/p/10229921.html
Pick定理、歐拉公式和圓的反演
出門左拐第三家：https://www.cnblogs.com/xzyxzy/p/10241872.html
Voronoi圖與Delaunay三角剖分
出門左拐第四家：https://www.cnblogs.com/xzyxzy/p/10349403.html