计算几何--二维几何基础练习
内容参考书籍——《算法竞赛入门经典训练指南》
例题1 题目链接：https://uva.onlinejudge.org/index.php?option=com_onlinejudge&Itemid=8&page=show_problem&problem=2119
莫利定理（Morley's theorem），也称为莫雷角三分线定理。将三角形的三个内角三等分，靠近某边的两条三分角线相交得到一个交点，则这样的三个交点可以构成一个正三角形。这个三角形常被称作莫利正三角形。
本题无算法可言，暴力求解即可(根据夹角根据旋转三角形的边，求出交点即为答案)。
代码如下：
1 #include <bits/stdc++.h> 2 using namespace std; 3 struct Point 4 { 5 double x,y; 6 Point(double x=0, double y=0):x(x),y(y) {}//构造函数，方便编写代码 7 }; 8 typedef Point Vector;//别名 9 //向量+向量=向量，点+向量=点 10 Vector operator + (Vector A, Vector B){return Vector(A.x+B.x,A.y+B.y);} 11 //点-点=向量 12 Vector operator - (Vector A, Vector B){return Vector(A.x-B.x,A.y-B.y);} 13 //向量*数=向量 14 Vector operator * (Vector A, double p){return Vector(A.x*p,A.y*p);} 15 //向量/数=向量 16 Vector operator / (Vector A, double p){return Vector(A.x/p,A.y/p);} 17 18 bool operator < (const Point& a, const Point& b){return a.x<b.x||(a.x==b.x&&a.y<b.y);} 19 20 const double eps = 1e-10; 21 int dcmp(double x) 22 { 23 if (fabs(x) < eps) return 0; 24 else return x < 0 ? -1 : 1; 25 } 26 27 bool operator == (const Point& a, const Point &b){return dcmp(a.x-b.x) == 0 && dcmp(a.y-b.y) == 0;} 28 29 //计算点积、向量长度、夹角函数 30 double Dot(Vector A, Vector B){return A.x*B.x+A.y*B.y;} 31 double Length(Vector A){return sqrt(Dot(A,A));} 32 double Angle(Vector A, Vector B){return acos(Dot(A,B)/Length(A)/Length(B));} 33 34 //计算叉积 35 double Cross(Vector A,Vector B){return A.x*B.y-A.y*B.x;} 36 double Area2(Point A, Point B, Point C){return Cross(B-A,C-A);} 37 38 //计算旋转和单位法向量 39 //rad为弧度 40 Vector Rotate(Vector A, double rad){return Vector(A.x*cos(rad)-A.y*sin(rad),A.x*sin(rad)+A.y*cos(rad));} 41 Vector Normal(Vector A) 42 { 43 double L =Length(A); 44 return Vector(-A.y/L,A.x/L); 45 } 46 47 //交点，注意精度问题 48 //调用该函数前请确保两条直线P+tv和Q+tw有唯一交点。当且仅当Cross(v,w)非0 49 Point GetLineIntersection(Point P, Vector v, Point Q, Vector w) 50 { 51 Vector u = P-Q; 52 double t = Cross(w,u) / Cross(v,w); 53 return P+v*t; 54 } 55 //点到直线的距离。用叉积计算。平行四边形的面积除以底 56 double DistanceToLine(Point P, Point A, Point B) 57 { 58 Vector v1 = B-A,v2= P-A; 59 return fabs(Cross(v1,v2)) / Length(v1); 60 } 61 //点到线段的距离 62 double DistanceToSegment(Point P, Point A, Point B) 63 { 64 if (A==B) return Length(P-A); 65 Vector v1 = B-A,v2 = P-A,v3=P-B; 66 if (dcmp(Dot(v1,v2))<0) return Length(v2); 67 else if (dcmp(Dot(v1,v3))>0) return Length(v3); 68 else return fabs(Cross(v1,v2)) / Length(v1); 69 } 70 //计算投影点 71 Point GetLineProjection(Point P, Point A, Point B) 72 { 73 Vector v = B-A; 74 return A+v*(Dot(v,P-A) / Dot(v,v)); 75 } 76 //按“相交规范”判断相交 77 bool SegmentProperIntersection(Point a1, Point a2, Point b1, Point b2) 78 { 79 double c1 = Cross(a2-a1,b1-a1),c2 = Cross(a2-a1,b2-a1), 80 c3 = Cross(b2-b1,a1-b1),c4 = Cross(b2-b1,a2-b1); 81 return dcmp(c1)*dcmp(c2)<0 && dcmp(c3)*dcmp(c4)<0; 82 } 83 //允许在端点处相交，判断是否发生加入下面一段判断一个点是否在一条线段上的代码 84 bool OnSegment(Point P,Point a1, Point a2) 85 { 86 return dcmp(Cross(a1-P,a2-P)) == 0 && dcmp(Dot(a1-P,a2-P))<0; 87 } 88 //多边形的有向面积 89 double PolygonArea(Point* p, int n) 90 { 91 double area = 0; 92 for (int i = 1; i < n-1; ++i) 93 { 94 area +=Cross(p[i]-p[0],p[i+1]-p[0]); 95 } 96 return area/2; 97 } 98 Point GetDaan(Point A,Point B,Point C) 99 { 100 //bc逆时针转 101 Vector v1 = C-B;//直线方向向量 102 double a1=Angle(A-B,v1); 103 v1 = Rotate(v1,a1/3); 104 //cb顺时针转 105 Vector v2 = B-C; 106 double a2=Angle(A-C,v2); 107 v2 = Rotate(v2,-a2/3);//负数表示顺时针旋转 108 109 return GetLineIntersection(B,v1,C,v2); 110 } 111 int main(int argc, char const *argv[]) 112 { 113 int t; 114 cin>>t; 115 while(t--) 116 { 117 Vector a,b,c; 118 Vector d,e,f;//答案 119 cin>>a.x>>a.y>>b.x>>b.y>>c.x>>c.y; 120 d=GetDaan(a,b,c); 121 e=GetDaan(b,c,a); 122 f=GetDaan(c,a,b); 123 printf("%.6f %.6f %.6f %.6f %.6f %.6f\n",d.x,d.y,e.x,e.y,f.x,f.y); 124 } 125 return 0; 126 }
例题2题目链接：https://uva.onlinejudge.org/index.php?option=com_onlinejudge&Itemid=8&page=submit_problem&problemid=2896&category=0
已知两狗以恒定速度奔跑但速度未知，同时出发同时达到，已知各自起点折线跑过程中的拐点和终点，求狗距。狗距=（奔跑过程中的最大距离减去最小距离）。
首先，我们考虑两条狗直线奔跑，那么我们可以认为一条狗不动，另外一条狗相对于这条狗奔跑。即问题转化为，求一个点到一条线段的最短距离和最大距离。
现在我们来模拟整个过程，用pa和pb来记录当前两条狗的位置，sa和sb为各自拐点的编号，每一次我们计算这两条狗谁先到达下一个拐点，那么在他到达下一个拐点的这段时间里便是上面说的直线奔跑的情况，不断跟新pa，pb，sa，sb，以及狗距最后的的结果便是答案。
代码如下：
1 #include <bits/stdc++.h> 2 using namespace std; 3 struct Point 4 { 5 double x,y; 6 Point(double x=0, double y=0):x(x),y(y) {}//构造函数，方便编写代码 7 }; 8 typedef Point Vector;//别名 9 //向量+向量=向量，点+向量=点 10 Vector operator + (Vector A, Vector B){return Vector(A.x+B.x,A.y+B.y);} 11 //点-点=向量 12 Vector operator - (Vector A, Vector B){return Vector(A.x-B.x,A.y-B.y);} 13 //向量*数=向量 14 Vector operator * (Vector A, double p){return Vector(A.x*p,A.y*p);} 15 //向量/数=向量 16 Vector operator / (Vector A, double p){return Vector(A.x/p,A.y/p);} 17 18 bool operator < (const Point& a, const Point& b){return a.x<b.x||(a.x==b.x&&a.y<b.y);} 19 20 const double eps = 1e-10; 21 int dcmp(double x) 22 { 23 if (fabs(x) < eps) return 0; 24 else return x < 0 ? -1 : 1; 25 } 26 27 bool operator == (const Point& a, const Point &b){return dcmp(a.x-b.x) == 0 && dcmp(a.y-b.y) == 0;} 28 29 //计算点积、向量长度、夹角函数 30 double Dot(Vector A, Vector B){return A.x*B.x+A.y*B.y;} 31 double Length(Vector A){return sqrt(Dot(A,A));} 32 double Angle(Vector A, Vector B){return acos(Dot(A,B)/Length(A)/Length(B));} 33 34 //计算叉积，叉积等于两向量组成的有向面积的两倍 35 double Cross(Vector A,Vector B){return A.x*B.y-A.y*B.x;} 36 double Area2(Point A, Point B, Point C){return Cross(B-A,C-A);} 37 38 //计算旋转和单位法向量 39 //rad为弧度 40 Vector Rotate(Vector A, double rad){return Vector(A.x*cos(rad)-A.y*sin(rad),A.x*sin(rad)+A.y*cos(rad));} 41 Vector Normal(Vector A) 42 { 43 double L =Length(A); 44 return Vector(-A.y/L,A.x/L); 45 } 46 47 //交点，注意精度问题 48 //调用该函数前请确保两条直线P+tv和Q+tw有唯一交点。当且仅当Cross(v,w)非0 49 Point GetLineIntersection(Point P, Vector v, Point Q, Vector w) 50 { 51 Vector u = P-Q; 52 double t = Cross(w,u) / Cross(v,w); 53 return P+v*t; 54 } 55 //点到直线的距离。用叉积计算。平行四边形的面积除以底 56 double DistanceToLine(Point P, Point A, Point B) 57 { 58 Vector v1 = B-A,v2= P-A; 59 return fabs(Cross(v1,v2)) / Length(v1); 60 } 61 //点到线段的距离 62 double DistanceToSegment(Point P, Point A, Point B) 63 { 64 if (A==B) return Length(P-A);//线段退化为点 65 Vector v1 = B-A,v2 = P-A,v3=P-B; 66 if (dcmp(Dot(v1,v2))<0) return Length(v2);//垂直 67 else if (dcmp(Dot(v1,v3))>0) return Length(v3);//P在AB的右边靠近B故为PB 68 else return fabs(Cross(v1,v2)) / Length(v1);//平行四边形面积处以底 69 } 70 //计算投影点 71 Point GetLineProjection(Point P, Point A, Point B) 72 { 73 Vector v = B-A; 74 return A+v*(Dot(v,P-A) / Dot(v,v)); 75 } 76 //按“相交规范”判断相交 77 bool SegmentProperIntersection(Point a1, Point a2, Point b1, Point b2) 78 { 79 double c1 = Cross(a2-a1,b1-a1),c2 = Cross(a2-a1,b2-a1), 80 c3 = Cross(b2-b1,a1-b1),c4 = Cross(b2-b1,a2-b1); 81 return dcmp(c1)*dcmp(c2)<0 && dcmp(c3)*dcmp(c4)<0; 82 } 83 //允许在端点处相交，判断是否发生加入下面一段判断一个点是否在一条线段上的代码 84 bool OnSegment(Point P,Point a1, Point a2) 85 { 86 return dcmp(Cross(a1-P,a2-P)) == 0 && dcmp(Dot(a1-P,a2-P))<0; 87 } 88 //多边形的有向面积 89 double PolygonArea(Point* p, int n) 90 { 91 double area = 0; 92 for (int i = 1; i < n-1; ++i) 93 { 94 area +=Cross(p[i]-p[0],p[i+1]-p[0]); 95 } 96 return area/2; 97 } 98 const int maxn = 60; 99 int T,A,B; 100 Point P[maxn],Q[maxn]; 101 double Min,Max; 102 void update(Point P, Point A, Point B) 103 { 104 Min = min(Min, DistanceToSegment(P, A, B));//点P到线段AB的距离 105 Max = max(Max, Length(P-A)); 106 Max = max(Max, Length(P-B)); 107 } 108 int main(int argc, char const *argv[]) 109 { 110 scanf("%d",&T); 111 for (int kase = 1; kase <= T; kase++) 112 { 113 scanf("%d%d",&A, &B); 114 for (int i = 0; i < A; ++i) cin>>P[i].x>>P[i].y; 115 for (int i = 0; i < B; ++i) cin>>Q[i].x>>Q[i].y; 116 117 double LenA = 0, LenB = 0; 118 for (int i = 0; i < A-1; ++i) LenA += Length(P[i+1]-P[i]); 119 for (int i = 0; i < B-1; ++i) LenB += Length(Q[i+1]-Q[i]); 120 121 int Sa = 0, Sb = 0; 122 Point Pa = P[0], Pb = Q[0]; 123 Min = 1e9,Max = -1e9; 124 while(Sa < A-1 && Sb < B-1) 125 { 126 double La = Length(P[Sa+1]-Pa); 127 double Lb = Length(Q[Sb+1]-Pb); 128 double TT = min(La/LenA,Lb/LenB); 129 Vector Va = (P[Sa+1]-Pa)/La*TT*LenA; 130 Vector Vb = (Q[Sb+1]-Pb)/Lb*TT*LenB; 131 update(Pa,Pb,Pb+Vb-Va); 132 Pa = Pa + Va; 133 Pb = Pb + Vb; 134 if (Pa == P[Sa+1]) Sa++; 135 if (Pb == Q[Sb+1]) Sb++; 136 } 137 printf("Case %d: %.0f\n", kase, Max-Min); 138 } 139 return 0; 140 }
