「OI之路」07其他-6计算几何

计算几何
more1
more2

辅助函数

1
2
double atan2(double y,double x)
long double atan2l(long double y,long double x)

向量相关

主要见线性代数一章

向量旋转:这个很好推,考虑两个基向量即可,$(a,b)逆时针旋转 \theta ->(acos\theta -bsin\theta ,asin\theta +bcos\theta )$

三角形的五心

外心:垂直平分线交点,外接圆圆心

重心:中线交点,重心到顶点=2*重心到对边中点,平分面积,到三个顶点距离的平方和最小,坐标为平均数,到三个顶点的向量和=零向量

垂心:垂线交点,三角形垂心到任一顶点=2*其外心到对边,单位圆上为三点坐标之和

内心:角平分线交点,直角三角形内心的$r=\frac{a+b-c}{2}$,圆内三点的三角形内心=那三段圆弧的中点形成三角形的垂心

欧拉线:外心O,重心G,垂心H,三点共线,2OG=GH

平面图转对偶图

如果是表格(如经典的BJOI狼抓兔子)很好转,但假如是个普通的平面图呢?例题:「HNOI2016」矿区

主要是给面编号,并且能快速得到一条边的两侧面的编号。一个小技巧是将无向边拆成两条相向的有向边,不妨用若干条顺时针连接的边(右侧有效的半平面)表示出一个面。于是枚举每个点,逆时针枚举每条还没用过的边$(x \to y)$,在y的出边中找到这条边逆时针的下一条边走(删除保证复杂度),这样走回原点就能得到一个面。一条边两侧的面就是两条有向边所属的面。另外还会得到一个无限大的面,其叉积和为负。

将例题讲完吧:因为依次将边断开后肯定能把图断开,随便求个生成树,只需要关心断掉的树边就能得到面的点集(具体而言,人云亦云,逆时针看,如果右侧半平面是父亲就加子树,否则减子树),$O(\sum d*log)$,实现时因为是整数所以面积乘二就是整数,最后除二即可。另外要注意对偶图中可能有重边,loj倒二code

其他零散知识

求多边形(可以凹)的面积:坐标原点到相邻两点的叉积和(逆时针求最后是正数),例题HDU2036。

适当的条件下(主要是考虑边界上咋办),还可以人为定义面积,例如下面的PKUSC2018D2T3 PKUSC

随机任意多边形的数据:据说,「WF2017」Airport Construction

随机恰n边的凸多边形的数据:随机n个和为0的向量,按角度排序,依次连接

欧拉公式

对于一个多边形,$顶点数V-边数E+面数F=2$,适用范围

推论:对于一个平面图,一定满足$m \leq 3n-6$(每条边最多给两个面使用,每个面至少有3条边,$r \leq \frac{2m}{3}$)

「HNOI2010」planar用过,这个推论,就背下不多就行了,应该还挺有用的

极角排序:

1
2
3
4
vc<int> L,R;fo(i,1,n) if(i!=st and i!=mid) (pt[st].cross(pt[mid],pt[i])>0?L:R).PB(i);
sort(all(L),[&](int a,int b){ return pt[st].cross(pt[a],pt[b])>0; });
sort(all(R),[&](int a,int b){ return pt[st].cross(pt[a],pt[b])>0; });
for(auto x:R) sot[st].PB(x);sot[st].PB(mid);for(auto x:L) sot[st].PB(x);//由mid分两半

直线交点
建议画图,一条直线与另外两点的叉积面积比=高的比=另一条直线被交点所分的两部分之比

1
2
3
4
5
6
Pt check(Pt a,Pt b,Pt c,Pt d)//尚未验证正确性
{
Pt t1=b-a,t2=d-c,t3=c-a;
ll A=cross(t1,t2),B=cross(t3,t2);
return a+t1*(double)B/A;
}

过圆外一点的切线:求点与圆心的直线然后用asin旋转

Pick定理:$S=内部整点+边界整点/2-1$,证明(还讲到用其证明Farey序列的性质)

例题CF559D Randomizer:给出n个整点的凸多边形,问所有构成的多边形不退化的点集得到的多边形严格内部整点数,$n \le 1e5,精度1e-9$

做法:拆成求面积、边界整点;面积拆分到每条边上,于是都可以用期望线性性摊到每条边$i \to j$上(多边形看做逆时针),$(x_1y_2-x_2y_1)/2$;至于边界整点显然是$\sum gcd(|x_1-x_2|,|y_1-y_2|)$,概率就是$\frac{2^{n-a}-1}{2^n-(1+n+C_n^2)}=(2^{-a}-2^{-n})/(1-2^{-n}(1+n+C_n^2))$,故a在一定范围内即可,code

基础题

雅礼冬令营集训2019 D7t3

PKUSC2018D2T3 PKUSC

板子

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
//计算几何
namespace Geo
{
const db PI=acos(-1);
db sqr(db x){return x*x;}
struct V//向量
{
db x,y;V(db _x=0,db _y=0){x=_x,y=_y;}
V operator + (V b){return V(x+b.x,y+b.y);}
V operator - (V b){return V(x-b.x,y-b.y);}
V operator - (){return V(-x,-y);}
V operator * (db b){return V(x*b,y*b);}
db operator * (V b){return x*b.y-y*b.x;}
db len(){return sqrt(x*x+y*y);}
V Unit(){return V(x/len(),y/len());}
db cross(V a,V b){return V(a.x-x,a.y-y)*V(b.x-x,b.y-y);}
V rotate(db rad){return (V){x*cos(rad)-y*sin(rad),x*sin(rad)+y*cos(rad)};}
bool operator == (V oth){return fabs(x-oth.x)+fabs(y-oth.y)<=eps;}
}; typedef V pt;
typedef pair<pt,db> cir;
typedef pair<pt,pt> line;
line cir_jd(cir a,cir b)//圆的交点
{
V oo=b.FR-a.FR;assert(oo.len()>eps);
db x=(sqr(a.SE)-sqr(b.SE)+sqr(oo.len()))/oo.len()/2;
pt mid=a.FR+oo*(x/oo.len());
V go=oo.Unit().rotate(PI/2)*sqrt(sqr(a.SE)-sqr(x));return {mid-go,mid+go};
}
db sp(V a,V b)//向量夹角
{
if(fabs(a.x)<=eps and fabs(a.y)<=eps) debug("err\n");
if(fabs(b.x)<=eps and fabs(b.y)<=eps) debug("err\n");
db c=atan2(a.y,a.x),d=atan2(b.y,b.x);
return min(fabs(c-d),PI*2-fabs(c-d));
}
V chui_zu(V a,V b)//向量垂足
{
db len=(b-a).len(),dis=fabs(a*b)/len;
return (a-b).rotate(PI/2)*(dis/len)*(a*b>0?1:-1);
}
};

本文基于 知识共享署名-相同方式共享 4.0 国际许可协议发布
本文地址:http://zory.ink/posts/7158.html
转载请注明出处,谢谢!