1.基础的基础
1.点与向量.
点与向量都用x,y两个坐标表示,所以我们在程序上暂时不做区分.
但是注意:
(1).点+向量=点;
(2).向量+向量=向量(向量-向量=向量+(相反向量));
(3).点-点=向量;
(4).点+点没有意义;
其中pt(point)表示点,vt(vector)表示向量.
1 struct pt{ 2 double x,y; 3 point(double x=0,double y=0):x(x),y(y){} 4 }; 5 typedef pt vt; 6 vt operator + (vt a,vt b){ return vt(a.x+b.x,a.y+b.y); } 7 vt operator - (pt a,pt b){ return vt(a.x-b.x,a.y-b.y); } 8 vt operator * (vt a,double x){ return vt(a.x*x,a.y*y); } 9 vt operator / (vt a,double x){ return vt(a.x/x,a.y/x); }
2.点的比较函数
在有些需要排序的时候会用到.
1 bool operator < (const pt &a,const pt &b){ return a.x<b.x||(a.x==b.x&&a.y<b.y); }
3.判断点相等的函数和避免误差的三态函数
浮点误差总是存在的,所以我们用三态函数dcmp来减少精度带来的问题.
1 const double eps=1e-8; 2 int dcmp(double x){ if(fabs(x)<eps) return 0; else return x<0?-1:1; } 3 bool operator == (const pt &a,const pt &b){ return dcmp(a.x-b.x)==0&&dcmp(a.y-b.y)==0; }
4.向量的极角
向量的极角表示从x轴正半轴旋转到该向量所需要的角度.
使用C标准库里面的atan2函数即可,单位是弧度.
1 atan2(y,x);
2.基本运算
1.向量的点积与叉积
$$cross(a,b)=a.x*b.y-a.y*b.x$$
这个是怎么来的呢?我们比较两个向量所在直线的斜率,如果b向量所在直线的斜率大(两个向量的斜率都存在的情况),那么就有$$frac{b.y}{b.x}>frac{a.y}{a.x}$$
然后去分母$$a.x*b.y>a.y*b.x$$
我们注意到,此时b向量所在的直线位于a向量所在直线的右手螺旋方向(a在b的左手螺旋方向).如果是小于号即b向量所在直线的斜率小,那么b向量所在的直线位于a向量所在直线的左手螺旋方向(a在b的右手螺旋方向).如果是等号即斜率一样,那么两条直线平行即向量a,b共线.
但这是两直线斜率都存在的情况.我们来讨论一下斜率不存在的情况:
(1).a所在直线的斜率不存在即(a.x=0).
若原不等式仍成立即[0>a.y*b.x],再讨论(a.y)的正负,可以证明b所在直线仍然位于a所在直线的右手螺旋方向.
(2).b所在直线的斜率不存在即(b.x=0).
若原不等式仍成立即$$a.x*b.y>0$$,再讨论(b.y)的正负,可以证明b所在直线仍然位于a所在直线的右手螺旋方向.
(3).a,b所在直线的斜率都不存在即(a.x=b.x=0).
此时不等式取等,平行.
综上,原不等式对应的结论永远成立.所以我们可以用叉积来判断两条直线的左右位置关系.
另外,叉积的几何意义是以这两个向量为邻边的平行四边形的有向面积(可用余弦定理证明),也就是两倍的三角形面积.
为什么说是有向面积呢?因为根据定义我们很容易发现(cross(a,b)=-cross(b,a)),所以算出来会带正负号!
1 double dot(vt a,vt b){ return a.x*b.x+a.y*b.y; } 2 double cross(vt a,vt b){ return a.x*b.y-a.y*b.x; } 3 double area2(pt a,pt b,pt c){ return cross(b-a,c-a); }
2.向量的模
1 double length(vt a){ return sqrt(dot(a,a)); }
3.两向量的夹角
通过点积来计算夹角的余弦值,再用一个反三角函数.
1 double angle(vt a,vt b){ return acos(dot(a,b)/length(a)/length(b)); }
4.向量的旋转.
将一个向量绕起点逆时针旋转(alpha)弧度.
怎么操作呢?
我们先把这个向量的起点平移到坐标原点.写出向量终点的极坐标:$$x=rcos heta$$
$$y=rsin heta$$
旋转之后终点的极坐标为:
$$x'=rcos( heta+alpha)$$
$$y'=rsin( heta+alpha)$$
打开得到:
$$x'=r(cos{ heta}cos{alpha}-sin{ heta}sin{alpha})=xcos{alpha}-ysin{alpha}$$
$$y'=r(sin{ heta}cos{alpha}+cos{ heta}sin{alpha})=xsin{alpha}+ycos{alpha}$$
即:
$$x'=xcos{alpha}-ysin{alpha}$$
$$y'=xsin{alpha}+ycos{alpha}$$
1 vt rotate(vt a,double rad){ return vt(a.x*cos(rad)-a.y*sin(rad),a.x*sin(rad)+a.y*cos(rad)) }
5.复数
白书上提到了,但是我没有搞懂用来干嘛...貌似也比较慢,FFT里都是手写的,不管了.
3.点和直线
1.线段相交
由于这是出现最频繁的问题,所以我们放在第一个说.
首先来定义一下规范相交:两线段相交且交点不是四个端点中的任何一个.
在解决判断线段是否相交的问题时,最容易想到的是解析几何的办法.但是和斜率相关的点斜式和斜截式都用不了,如果用一般式或两点式需要联立解方程,然后判断是否有解,如果有解还要判断是否在两条线段内部.这样很麻烦,而且误差不好控制.
对于不要求求出交点坐标的问题,我们采用一种神奇的方法.
首先一个引理:两线段(AB)与(CD)规范相交(Longleftrightarrow)(overrightarrow{AC})与(overrightarrow{AD})分别在(overrightarrow{AB})左右两侧且(overrightarrow{CA})与(overrightarrow{CB})分别在(overrightarrow{CD})左右两侧
所以我们用叉积做就好了.
1 bool segment_proper_intersection(pt a,pt b,pt c,pt d){ 2 return(dcmp(cross(b-a,c-a))^dcmp(cross(b-a,d-a))==-2&& 3 dcmp(cross(d-c,a-c))^dcmp(cross(d-c,b-c))==-2); 4 }
注意这里的(dcmp(cross(b-a,c-a)^cross(b-a,d-a))==-2)(Longleftrightarrow)(dcmp(cross(b-a,c-a),cross(b-a,d-a))<0).
因为小于0的情况是一个是1,另一个是-1.在计算的运算中,全部用补码(正数不变,负数取反+1).
这样1的补码是0000...001
-1的补码是1111...111
异或之后结果是1111...110
开头是1,是个负数,-1之后变为了1111...101,然后取反(从反码变回原码),1000...010,也就是-2.
那么对于不规范相交呢?我们发现不规范相交都会有四个端点中至少一个端点在另一个线段上.所以在用叉积进行判断的时候,如果两个向量的叉积为0,就说明有一个点在另一个线段所在的直线上,如果这个点在那个线段上,就是不规范相交啦.
所以我们还需要写一个函数来判断直线上的一个点是否在直线的某一条线段上.
我们可以用坐标比较的方法.还有第二种方法是用点积.叉积可以用来比较左右位置,而点积可以用来比较前后位置.
设这个点为(P),线段为(AB).如果(overrightarrow{PA})与(overrightarrow{PB})的点积小于等于0,就说明(P)在线段(AB)上(当(P)与(A)或(B)重合时取等).
1 int between(pt x,pt a,pt b){ 2 if(fabs(a.x-b.x)>fabs(a.x-b.y)) return dcmp(x.x-a.x)*dcmp(x.x-b.x); 3 else return dcmp(x.y-a.y)*dcmp(x.y-b.y); 4 } 5 int between(pt x,pt a,pt b){ 6 return dcmp(rot(x-a,x-b)); 7 } 8 int segment_intersection(pt a,pt b,pt c,pt d){ 9 d1=dcmp(cross(b-a,c-a)); 10 d2=dcmp(cross(b-a,d-a)); 11 d3=dcmp(cross(d-c,a-c)); 12 d4=dcmp(cross(d-c,b-c)); 13 if(d1==0&&between(c,a,b)<=0|| 14 d2==0&&between(d,a,b)<=0|| 15 d3==0&&between(a,c,d)<=0|| 16 d4==0&&between(b,c,d)<=0) 17 return 2; 18 return 0; 19 }
2.求直线交点
这回是真的要求交点了.
首先我们用参数方程表示直线(一个点+一个向量)
1 struct line{ 2 pt p; vt v; 3 line(){} 4 line(pt a,pt b){ p=a; v=b-a; } 5 };
直线上所有点(Q)都满足(Q=P+toverrightarrow{v}).
我们设两条直线分别为(P+toverrightarrow{v})和(Q+toverrightarrow{w}).设交点在这两条直线上的参数分别为(t_1,t_2),作别为((x,y)).
列出方程
$$x=x_P+t_1x_v=x_Q+t_2x_w$$
$$y=y_P+t_1y_v=y_Q+t_2y_w$$
我们设(overrightarrow{u}=P-Q),则(x_u=x_P-x_Q,y_u=y_P-y_Q)
由1式得到:
$$t_2=frac{x_u+t_1x_v}{x_w}$$
带入2式消元,解得:
$$t_1=frac{cross(overrightarrow{w},overrightarrow{u})}{cross(overrightarrow{v},overrightarrow{w})}$$
这样就可以得到交点:
$$P+v_1t_1$$
1 pt get_line_intersection(pt p,vt v,pt q,vt w){ 2 vt u=p-q; 3 double t=cross(w,u)/cross(v,w); 4 return p+v*t; 5 }
注意到这里的(cross(overrightarrow{v},overrightarrow{w}))不能为0,也就是说(overrightarrow{v},overrightarrow{w})不能共线,也就是说在调用函数之前一定要确保两直线不平行.
3.求线段交点
我们可以用叉积的几何意义来求线段的交点.图嘛...看黑书P357吧...我这里大概推一下.
首先通过三角形等高容易证明:
$$frac{|DP|}{|CP|}=frac{S_{{ riangle}ADB}}{S_{{ riangle}ACB}}=frac{overrightarrow{AD} imes{overrightarrow{AB}}}{overrightarrow{AC} imes{overrightarrow{AB}}}$$
然后我们用向量来证明黑书下面那个公式:
$$overrightarrow{OP}=overrightarrow{OD}+frac{|DP|}{|DP|+|CP|}overrightarrow{DC}$$
我们设(S_1=S_{{ riangle}ADB},S_2=S_{{ riangle}ACB})
$$x_P=x_D+frac{S_1}{S_1+S_2}(x_C-x_D)$$
$$x_P=frac{S_1x_C+S_2x_D}{S_1+S_2}$$
同理可以算出纵坐标.
1 int segment_intersection(pt a,pt b,pt c,pt d,pt &p){ 2 double s1,s2; 3 int d1,d2,d3,d4; 4 d1=dcmp(s1=cross(b-a,c-a)); 5 d2=dcmp(s2=cross(b-a,d-a)); 6 d3=dcmp(cross(d-c,a-c)); 7 d4=dcmp(cross(d-c,b-c)); 8 if((d1^d2)==-2&&(d3^d4)==-2){ 9 p.x=(c.x*s2-d.x*s1)/(s2-s1); 10 p.y=(c.y*s2-d.y*s1)/(s2-s1); 11 return 1; 12 } 13 if(d1==0&&between(c,a,b)<=0|| 14 d2==0&&between(d,a,b)<=0|| 15 d3==0&&between(a,c,d)<=0|| 16 d4==0&&between(b,c,d)<=0) 17 return 2; 18 return 0; 19 }
好了我写不下去了,明天再继续吧...今天晚上回去看凸包,半平面交,不看完是小狗...