首页 > 代码库 > 几何模板总结——算法竞赛入门经典(第二版)
几何模板总结——算法竞赛入门经典(第二版)
1 #include <iostream> 2 #include <string> 3 #include <cstdio> 4 #include <cstring> 5 #include <cmath> 6 #include <vector> 7 #include <algorithm> 8 9 using namespace std; 10 11 const double PI = acos(-1); 12 13 struct Point { 14 double x, y; 15 Point(double x = 0, double y = 0) : x(x), y(y) {} 16 }; 17 18 typedef Point Vector; 19 Vector operator + (const Vector &A, const Vector &B) { return Vector(A.x + B.x, A.y + B.y); } 20 Vector operator - (const Point &A, const Point &B) { return Vector(A.x - B.x, A.y - B.y); } 21 Vector operator * (const Vector &A, const double &p) { return Vector(A.x * p, A.y * p); } 22 Vector operator / (const Vector &A, const double &p) { return Vector(A.x / p, A.y / p); } 23 bool operator < (const Point &a, const Point &b) { return a.x < b.x || (a.x == b.x && a.y < b.y); } 24 25 const double eps = 1e-10; 26 int dcmp(double x) { 27 if(fabs(x) < eps) return 0; 28 else return x < 0 ? -1 : 1; 29 } 30 31 bool operator == (const Point &a, const Point &b) { return dcmp(a.x - b.x) == 0 && dcmp(a.y - b.y) == 0; } 32 33 double Dot(const Vector &A, const Vector &B) { return A.x * B.x + A.y * B.y; } 34 double Length(const Vector &A) { return sqrt(Dot(A, A)); } 35 double Angle(const Vector &A, const Vector &B) { return acos(Dot(A, B) / Length(A) / Length(B)); } 36 37 double Cross(const Vector &A, const Vector B) { return A.x * B.y - A.y * B.x; } 38 double Area2(const Point &A, const Point &B, const Point &C) { return Cross(B - A, C - A); } 39 40 Vector Rotate(const Vector &A, double rad) { return Vector(A.x * cos(rad) - A.y * sin(rad), A.x * sin(rad) + A.y * cos(rad)); } 41 42 Vector Normal(const Vector &A) { 43 double L = Length(A); 44 return Vector(-A.y / L, A.x / L); 45 } 46 47 Point GetLineIntersection(const Point &P, const Vector &v, const Point &Q, const Vector &w) { 48 Vector u = P - Q; 49 double t = Cross(w, u) / Cross(v, w); 50 return P + v * t; 51 } 52 53 double DistanceToLine(const Point &P, const Point &A, const Point &B) { 54 Vector v1 = B - A, v2 = P - A; 55 return fabs(Cross(v1, v2)) / Length(v1); 56 } 57 double DistanceToLine_(const Point &P, const Point &A, const Point &B) { 58 Vector v1 = B - A, v2 = P - A; 59 return Cross(v1, v2) / Length(v1); 60 } 61 double DistanceToSegment(const Point &P, const Point &A, const Point &B) { 62 if(A == B) return Length(P - A); 63 Vector v1 = B - A, v2 = P - A, v3 = P - B; 64 if(dcmp(Dot(v1, v2)) < 0) return Length(v2); 65 else if(dcmp(Dot(v1, v3)) > 0) return Length(v3); 66 else return fabs(Cross(v1, v2)) / Length(v1); 67 } 68 69 Point GetLineProjection(const Point &P, const Point &A, const Point &B) { 70 Vector v = B - A; 71 return A + v * (Dot(v, P - A) / Dot(v, v)); 72 } 73 74 bool SegmentProperIntersection(const Point &a1, const Point &a2, const Point &b1, const Point &b2) { 75 double c1 = Cross(a2 - a1, b1 - a1), c2 = Cross(a2 - a1, b2 - a1), c3 = Cross(b2 - b1, a1 - b1), c4 = Cross(b2 - b1, a2 - b1); 76 return dcmp(c1) * dcmp(c2) < 0 && dcmp(c3) * dcmp(c4) < 0; 77 } 78 79 bool OnSegment(const Point &p, const Point &a1, const Point &a2) { 80 return dcmp(Cross(a1 - p, a2 - p)) == 0 && dcmp(Dot(a1 - p, a2 - p)) < 0; 81 } 82 83 double ConvexPolygonArea(Point *p, int n) { 84 double area = 0; 85 for(int i = 1; i < n - 1; ++i) area += Cross(p[i] - p[0], p[i + 1] - p[0]); 86 return area / 2; 87 } 88 89 double PolygonArea(Point *p, int n) { 90 double area = 0; 91 92 for(int i = 1; i < n - 1; ++i) area += Cross(p[i] - p[0], p[i + 1] - p[0]); 93 return area / 2; 94 } 95 96 struct Line { 97 Point p; 98 Vector v; 99 Line(Point p, Vector v):p(p),v(v) { }100 Point point(double t) {101 return p + v*t;102 }103 Line move(double d) {104 return Line(p + Normal(v)*d, v);105 }106 };107 108 Line getLine(double x1, double y1, double x2, double y2) {109 Point p1(x1,y1);110 Point p2(x2,y2);111 return Line(p1, p2-p1);112 }113 114 struct Circle {115 Point c;116 double r;117 118 Circle(Point c, double r) : c(c), r(r) {}119 120 Point point(const double &a) const {121 return Point(c.x + cos(a) * r, c.y + sin(a) * r);122 }123 };124 125 int getLineCircleIntersection(Line L, Circle C, double &t1, double &t2, vector<Point> &sol) {126 double a = L.v.x, b = L.p.x - C.c.x, c = L.v.y, d = L.p.y - C.c.y;127 double e = a * a + c * c, f = 2 * (a * b + c * d), g = b * b + d * d - C.r * C.r;128 double delta = f * f - 4 * e * g;129 if(dcmp(delta) < 0) return 0;130 if(dcmp(delta) == 0) {131 t1 = t2 = -f / (2 * e);132 sol.push_back(L.point(t1));133 return 1;134 }135 t1 = (-f - sqrt(delta)) / (2 * e);136 sol.push_back(L.point(t1));137 t2 = (-f + sqrt(delta)) / (2 * e);138 sol.push_back(L.point(t2));139 return 2;140 }141 142 double angle(const Vector &v) { return atan2(v.y, v.x); }143 144 int getCircleCircleIntersection(Circle C1, Circle C2, vector<Point> &sol) {145 double d = Length(C1.c - C2.c);146 if(dcmp(d) == 0) {147 if(dcmp(C1.r - C2.r) == 0) return -1;148 return 0;149 }150 if(dcmp(C1.r + C2.r - d) < 0) return 0;151 if(dcmp(fabs(C1.r - C2.r) - d) > 0) return 0;152 153 double a = angle(C2.c - C1.c);154 double da = acos((C1.r * C1.r + d * d - C2.r * C2.r) / (2 * C1.r * d));155 156 Point p1 = C1.point(a - da), p2 = C1.point(a + da);157 158 sol.push_back(p1);159 if(p1 == p2) return 1;160 sol.push_back(p2);161 return 2;162 }163 164 int getTangents(Point p, Circle C, Vector *v) {165 Vector u = C.c - p;166 double dist = Length(u);167 if(dist < C.r) return 0;168 else if(dcmp(dist - C.r) == 0) {169 // v[0] = Rotate(u, PI / 2);170 return 1;171 }172 else {173 double ang = asin(C.r / dist);174 v[0] = Rotate(u, -ang);175 v[1] = Rotate(u, +ang);176 return 2;177 }178 }179 180 int getTangents(Circle A, Circle B, Point *a, Point *b) {181 int cnt = 0;182 if(A.r < B.r) { swap(A, B); swap(a, b); }183 int d2 = (A.c.x - B.c.x) * (A.c.x - B.c.x) + (A.c.y - B.c.y) * (A.c.y - B.c.y);184 int rdiff = A.r - B.r;185 int rsum = A.r + B.r;186 if(d2 < rdiff * rdiff) return 0;187 188 double base = atan2(B.c.y - A.c.y, B.c.x - A.c.x);189 if(d2 == 0 && A.r == B.r) return -1;190 if(d2 == rdiff * rdiff) {191 a[cnt] = A.point(base);192 b[cnt] = B.point(base);193 ++cnt;194 return 1;195 }196 double ang = acos((A.r - B.r) / sqrt(d2));197 a[cnt] = A.point(base + ang);198 b[cnt] = B.point(base + ang);199 ++cnt;200 a[cnt] = A.point(base - ang);201 b[cnt] = B.point(base - ang);202 ++cnt;203 if(d2 == rsum * rsum) {204 a[cnt] = A.point(base);205 b[cnt] = B.point(PI + base);206 ++cnt;207 }208 else if(d2 > rsum * rsum) {209 double ang = acos((A.r + B.r) / sqrt(d2));210 a[cnt] = A.point(base + ang);211 b[cnt] = B.point(PI + base + ang);212 ++cnt;213 a[cnt] = A.point(base - ang);214 b[cnt] = B.point(PI + base - ang);215 ++cnt;216 }217 return cnt;218 }219 220 double torad(const double °) {221 return deg / 180 * PI;222 }223 224 void get_coord(const double &R, double lat, double lng, double &x, double &y, double &z) {225 lat = torad(lat);226 lng = torad(lng);227 x = R * cos(lat) * cos(lng);228 y = R * cos(lat) * sin(lng);229 z = R * sin(lat);230 }231 232 int ConvexHull(Point *p, int n, Point *ch) {233 sort(p, p + n);234 int m = 0;235 for(int i = 0; i < n; ++i) {236 while(m > 1 && Cross(ch[m - 1] - ch[m - 2], p[i] - ch[m - 2]) <= 0) --m;237 ch[m++] = p[i];238 }239 int k = m;240 for(int i = n - 2; i >= 0; --i) {241 while(m > k && Cross(ch[m - 1] - ch[m - 2], p[i] - ch[m - 2]) <= 0) --m;242 ch[m++] = p[i];243 }244 if(n > 1) --m;245 return m;246 }247 248 int isPointInPolygon(Point p, Point *poly, int n) {249 int wn = 0;250 for(int i = 0;i < n; ++i) {251 if(OnSegment(p, poly[i], poly[(i+1)%n])) return -1;252 int k=dcmp(Cross(poly[(i+1)%n]-poly[i], p-poly[i]));253 int d1=dcmp(poly[i].y-p.y);254 int d2=dcmp(poly[(i+1)%n].y-p.y);255 if(k>0&&d1<=0&&d2>0) wn++;256 if(k<0&&d2<=0&&d1>0) wn--;257 }258 if(wn!=0) return 1;259 else return 0;260 }261 262 Point read_point() {263 Point P;264 scanf("%lf%lf",&P.x,&P.y);265 return P;266 }267 268 int main() {269 270 }
几何模板总结——算法竞赛入门经典(第二版)
声明:以上内容来自用户投稿及互联网公开渠道收集整理发布,本网站不拥有所有权,未作人工编辑处理,也不承担相关法律责任,若内容有误或涉及侵权可进行投诉: 投诉/举报 工作人员会在5个工作日内联系你,一经查实,本站将立刻删除涉嫌侵权内容。