平面几何

向量（Vector）
在几乎所有的几何问题中，向量（有时也称矢量）是一个基本点。向量的定义包含方向和一个数（长度）。在二维空间中，一个向量可以用一对x和y来表示。例如由点（1,3）到（5,1的向量可以用(4,-2）来表示。这里大家要特别注意，我这样说并不代表向量定义了起点和终点。向量仅仅定义方向和长度。

向量加法
向量也支持各种数学运算。最简单的就是加法。我们可以对两个向量相加，得到的仍然是一个向量。我们有：
    V1（x1, y1）+V2（x2, y2）=V3(x1+x2, y1+y2)
下图表示了四个向量相加。注意就像普通的加法一样,相加的次序对结果没有影响（满足交换律），减法也是一样的。


点乘（Dot Product）
如果说加法是凭直觉就可以知道的，另外还有一些运算就不是那么明显的，比如点乘和叉乘。
点乘比较简单，是相应元素的乘积的和：
    V1( x1, y1)   V2(x2, y2) = x1*x2 + y1*y2
注意结果不是一个向量，而是一个标量（Scalar）。点乘有什么用呢，我们有：
    A   B = |A||B|Cos(θ)
θ是向量A和向量B见的夹角。这里|A|我们称为向量A的模(norm)，也就是A的长度， 在二维空间中就是|A| = sqrt(x2+y2)。这样我们就和容易计算两条线的夹角：
    Cos(θ) = A   B /(|A||B|)

当然你知道要用一下反余弦函数acos()啦。（回忆一下cos(90)=0 和cos(0) = 1还是有好处的，希望你没有忘记。）这可以告诉我们如果点乘的结果，简称点积，为0的话就表示这两个向量垂直。当两向量平行时，点积有最大值
另外，点乘运算不仅限于2维空间，他可以推广到任意维空间。（译注：不少人对量子力学中的高维空间无法理解，其实如果你不要试图在视觉上想象高维空间，而仅仅把它看成三维空间在数学上的推广，那么就好理解了）


叉乘（cross product）
相对于点乘，叉乘可能更有用吧。2维空间中的叉乘是：
    V1(x1, y1) X V2(x2, y2) = x1y2 – y1x2
看起来像个标量，事实上叉乘的结果是个向量，方向在z轴上。上述结果是它的模。在二维空间里，让我们暂时忽略它的方向，将结果看成一个向量，那么这个结果类似于上述的点积，我们有：
    A x B = |A||B|Sin(θ)
然而角度 θ和上面点乘的角度有一点点不同，他是有正负的，是指从A到B的角度。下图中 θ为负。
另外还有一个有用的特征那就是叉积的绝对值就是A和B为两边说形成的平行四边形的面积。也就是AB所包围三角形面积的两倍。在计算面积时，我们要经常用到叉积。
（译注：三维及以上的叉乘参看维基：http://en.wikipedia.org/wiki/Cross_product）


点-线距离
找出一个点和一条线间的距离是经常遇见的几何问题之一。假设给出三个点，A，B和C，你想找出点C到点A、B定出的直线间距离。第一步是找出A到B的向量AB和A到C的向量AC，现在我们用该两向量的叉积除以|AB|，这就是我们要找的的距离了（下图中的红线）。
    d = (AB x AC)/|AB|


如果你有基础的高中几何知识，你就知道原因了。上一节我们知道(AB X AC)/2是三角形ABC的面积，这个三角形的底是|AB|，高就是C到AB的距离。有时叉积得到的是一个负值，这种情况下距离就是上述结果的绝对值。
当我们要找点到线段的距离时，情况变得稍稍复杂一些。这时线段与点的最短距离可能是点到线段的某一端点，而不是点到直线的垂线。例如上图中点C到线段AB的最短距离应该是线段BC。我们有集中不同的方法来判断这种特殊情况。第一种情况是计算点积AB BC来判定两线段间夹角。如果点积大于等于零，那么表示AB到BC是在-90到90度间，也就是说C到AB的垂线在AB外，那么AB上到C距离最近的点就是B。同样，如果BAAC大于等于零，那么点A就是距离C最近的点。如果两者均小于零，那么距离最近的点就在线段AB中的莫一点。

 

代码

1 3 . 浮点函数
2   // 浮点几何函数库
3 #include < math.h >
4 #define eps 1e-8
5 #define zero(x) (((x)>0?(x):-(x))<eps)
6 struct point{ double x,y;};
7 struct line{point a,b;};
8
9 // 计算cross product (P1-P0)x(P2-P0)
10 double xmult(point p1,point p2,point p0){
11 return (p1.x - p0.x) * (p2.y - p0.y) - (p2.x - p0.x) * (p1.y - p0.y);
12 }
13 double xmult( double x1, double y1, double x2, double y2, double x0, double y0){
14 return (x1 - x0) * (y2 - y0) - (x2 - x0) * (y1 - y0);
15 }
16
17 // 计算dot product (P1-P0).(P2-P0)
18 double dmult(point p1,point p2,point p0){
19 return (p1.x - p0.x) * (p2.x - p0.x) + (p1.y - p0.y) * (p2.y - p0.y);
20 }
21 double dmult( double x1, double y1, double x2, double y2, double x0, double y0){
22 return (x1 - x0) * (x2 - x0) + (y1 - y0) * (y2 - y0);
23 }
24
25 // 两点距离
26 double distance(point p1,point p2){
27 return sqrt((p1.x - p2.x) * (p1.x - p2.x) + (p1.y - p2.y) * (p1.y - p2.y));
28 }
29 double distance( double x1, double y1, double x2, double y2){
30 return sqrt((x1 - x2) * (x1 - x2) + (y1 - y2) * (y1 - y2));
31 }
32
33 // 判三点共线
34 int dots_inline(point p1,point p2,point p3){
35 return zero(xmult(p1,p2,p3));
36 }
37 int dots_inline( double x1, double y1, double x2, double y2, double x3, double y3){
38 return zero(xmult(x1,y1,x2,y2,x3,y3));
39 }
40
41 // 判点是否在线段上,包括端点
42 int dot_online_in(point p,line l){
43 return zero(xmult(p,l.a,l.b)) && (l.a.x - p.x) * (l.b.x - p.x) < eps && (l.a.y - p.y) * (l.b.y - p.y) < eps;
44 }
45 int dot_online_in(point p,point l1,point l2){
46 return zero(xmult(p,l1,l2)) && (l1.x - p.x) * (l2.x - p.x) < eps && (l1.y - p.y) * (l2.y - p.y) < eps;
47 }
48 int dot_online_in( double x, double y, double x1, double y1, double x2, double y2){
49 return zero(xmult(x,y,x1,y1,x2,y2)) && (x1 - x) * (x2 - x) < eps && (y1 - y) * (y2 - y) < eps;
50 }
51
52 // 判点是否在线段上,不包括端点
53 int dot_online_ex(point p,line l){
54 return dot_online_in(p,l) && ( ! zero(p.x - l.a.x) ||! zero(p.y - l.a.y)) && ( ! zero(p.x - l.b.x) ||! zero(p.y - l.b.y));
55 }
56 int dot_online_ex(point p,point l1,point l2){
57 return dot_online_in(p,l1,l2) && ( ! zero(p.x - l1.x) ||! zero(p.y - l1.y)) && ( ! zero(p.x - l2.x) ||! zero(p.y - l2.y));
58 }
59 int dot_online_ex( double x, double y, double x1, double y1, double x2, double y2){
60 return dot_online_in(x,y,x1,y1,x2,y2) && ( ! zero(x - x1) ||! zero(y - y1)) && ( ! zero(x - x2) ||! zero(y - y2));
61 }
62
63 // 判两点在线段同侧,点在线段上返回0
64 int same_side(point p1,point p2,line l){
65 return xmult(l.a,p1,l.b) * xmult(l.a,p2,l.b) > eps;
66 }
67 int same_side(point p1,point p2,point l1,point l2){
68 return xmult(l1,p1,l2) * xmult(l1,p2,l2) > eps;
69 }
70
71 // 判两点在线段异侧,点在线段上返回0
72 int opposite_side(point p1,point p2,line l){
73 return xmult(l.a,p1,l.b) * xmult(l.a,p2,l.b) <- eps;
74 }
75 int opposite_side(point p1,point p2,point l1,point l2){
76 return xmult(l1,p1,l2) * xmult(l1,p2,l2) <- eps;
77 }
78
79 // 点关于直线的对称点 // by lyt
80 // 缺点：用了斜率
81 // 也可以利用"点到直线上的最近点"来做，避免使用斜率。
82 point symmetric_point(point p1, point l1, point l2) {
83 point ret;
84 if (l1.x > l2.x - eps && l1.x < l2.x + eps) {
85 ret.x = ( 2 * l1.x - p1.x);
86 ret.y = p1.y;
87 } else {
88 double k = (l1.y - l2.y ) / (l1.x - l2.x);
89 ret.x = ( 2 * k * k * l1.x + 2 * k * p1.y - 2 * k * l1.y - k * k * p1.x + p1.x) / ( 1 + k * k);
90 ret.y = p1.y - (ret.x - p1.x ) / k;
91 }
92 return ret;
93 }
94
95 // 判两直线平行
96 int parallel(line u,line v){
97 return zero((u.a.x - u.b.x) * (v.a.y - v.b.y) - (v.a.x - v.b.x) * (u.a.y - u.b.y));
98 }
99 int parallel(point u1,point u2,point v1,point v2){
100 return zero((u1.x - u2.x) * (v1.y - v2.y) - (v1.x - v2.x) * (u1.y - u2.y));
101 }
102
103 // 判两直线垂直
104 int perpendicular(line u,line v){
105 return zero((u.a.x - u.b.x) * (v.a.x - v.b.x) + (u.a.y - u.b.y) * (v.a.y - v.b.y));
106 }
107 int perpendicular(point u1,point u2,point v1,point v2){
108 return zero((u1.x - u2.x) * (v1.x - v2.x) + (u1.y - u2.y) * (v1.y - v2.y));
109 }
110
111 // 判两线段相交,包括端点和部分重合
112 int intersect_in(line u,line v){
113 if ( ! dots_inline(u.a,u.b,v.a) ||! dots_inline(u.a,u.b,v.b))
114 return ! same_side(u.a,u.b,v) &&! same_side(v.a,v.b,u);
115 return dot_online_in(u.a,v) || dot_online_in(u.b,v) || dot_online_in(v.a,u) || dot_online_in(v.b,u);
116 }
117 int intersect_in(point u1,point u2,point v1,point v2){
118 if ( ! dots_inline(u1,u2,v1) ||! dots_inline(u1,u2,v2))
119 return ! same_side(u1,u2,v1,v2) &&! same_side(v1,v2,u1,u2);
120 return dot_online_in(u1,v1,v2) || dot_online_in(u2,v1,v2) || dot_online_in(v1,u1,u2) || dot_online_in(v2,u1,u2);
121 }
122
123 // 判两线段相交,不包括端点和部分重合
124 int intersect_ex(line u,line v){
125 return opposite_side(u.a,u.b,v) && opposite_side(v.a,v.b,u);
126 }
127 int intersect_ex(point u1,point u2,point v1,point v2){
128 return opposite_side(u1,u2,v1,v2) && opposite_side(v1,v2,u1,u2);
129 }
130
131 // 计算两直线交点,注意事先判断直线是否平行!
132 // 线段交点请另外判线段相交(同时还是要判断是否平行!)
133 point intersection(line u,line v){
134 point ret = u.a;
135 double t = ((u.a.x - v.a.x) * (v.a.y - v.b.y) - (u.a.y - v.a.y) * (v.a.x - v.b.x))
136 / ((u.a.x - u.b.x) * (v.a.y - v.b.y) - (u.a.y - u.b.y) * (v.a.x - v.b.x));
137 ret.x += (u.b.x - u.a.x) * t;
138 ret.y += (u.b.y - u.a.y) * t;
139 return ret;
140 }
141 point intersection(point u1,point u2,point v1,point v2){
142 point ret = u1;
143 double t = ((u1.x - v1.x) * (v1.y - v2.y) - (u1.y - v1.y) * (v1.x - v2.x))
144 / ((u1.x - u2.x) * (v1.y - v2.y) - (u1.y - u2.y) * (v1.x - v2.x));
145 ret.x += (u2.x - u1.x) * t;
146 ret.y += (u2.y - u1.y) * t;
147 return ret;
148 }
149
150 // 点到直线上的最近点
151 point ptoline(point p,line l){
152 point t = p;
153 t.x += l.a.y - l.b.y,t.y += l.b.x - l.a.x;
154 return intersection(p,t,l.a,l.b);
155 }
156 point ptoline(point p,point l1,point l2){
157 point t = p;
158 t.x += l1.y - l2.y,t.y += l2.x - l1.x;
159 return intersection(p,t,l1,l2);
160 }
161
162 // 点到直线距离
163 double disptoline(point p,line l){
164 return fabs(xmult(p,l.a,l.b)) / distance(l.a,l.b);
165 }
166 double disptoline(point p,point l1,point l2){
167 return fabs(xmult(p,l1,l2)) / distance(l1,l2);
168 }
169 double disptoline( double x, double y, double x1, double y1, double x2, double y2){
170 return fabs(xmult(x,y,x1,y1,x2,y2)) / distance(x1,y1,x2,y2);
171 }
172
173 // 点到线段上的最近点
174 point ptoseg(point p,line l){
175 point t = p;
176 t.x += l.a.y - l.b.y,t.y += l.b.x - l.a.x;
177 if (xmult(l.a,t,p) * xmult(l.b,t,p) > eps)
178 return distance(p,l.a) < distance(p,l.b) ? l.a:l.b;
179 return intersection(p,t,l.a,l.b);
180 }
181 point ptoseg(point p,point l1,point l2){
182 point t = p;
183 t.x += l1.y - l2.y,t.y += l2.x - l1.x;
184 if (xmult(l1,t,p) * xmult(l2,t,p) > eps)
185 return distance(p,l1) < distance(p,l2) ? l1:l2;
186 return intersection(p,t,l1,l2);
187 }
188
189 // 点到线段距离
190 double disptoseg(point p,line l){
191 point t = p;
192 t.x += l.a.y - l.b.y,t.y += l.b.x - l.a.x;
193 if (xmult(l.a,t,p) * xmult(l.b,t,p) > eps)
194 return distance(p,l.a) < distance(p,l.b) ? distance(p,l.a):distance(p,l.b);
195 return fabs(xmult(p,l.a,l.b)) / distance(l.a,l.b);
196 }
197 double disptoseg(point p,point l1,point l2){
198 point t = p;
199 t.x += l1.y - l2.y,t.y += l2.x - l1.x;
200 if (xmult(l1,t,p) * xmult(l2,t,p) > eps)
201 return distance(p,l1) < distance(p,l2) ? distance(p,l1):distance(p,l2);
202 return fabs(xmult(p,l1,l2)) / distance(l1,l2);
203 }
204
205 // 矢量V以P为顶点逆时针旋转angle并放大scale倍
206 point rotate(point v,point p, double angle, double scale){
207 point ret = p;
208 v.x -= p.x,v.y -= p.y;
209 p.x = scale * cos(angle);
210 p.y = scale * sin(angle);
211 ret.x += v.x * p.x - v.y * p.y;
212 ret.y += v.x * p.y + v.y * p.x;
213 return ret;
214 }
215
216