本文共 8274 字，大约阅读时间需要 27 分钟。
#include #include #include #include #include #include #include #include #include #include #define MM(a) memset(a,0,sizeof(a)) typedef long long ll; typedef unsigned long long ULL; const double eps = 1e-14; const int inf = 0x3f3f3f3f; const double pi=acos(-1); using namespace std; struct Point { double x, y; Point() {} Point(double x, double y) { this->x = x; this->y = y; } void read() { scanf("%lf%lf", &x, &y); } }; bool cmp(Point a,Point b) { if(a.x==b.x) return a.y c = c; this->r = r; } Point point(double a) { return Point(c.x + cos(a) * r, c.y + sin(a) * r); } }; typedef Point Vector; Vector operator + (Vector A, Vector B) { return Vector(A.x + B.x, A.y + B.y); } Vector operator - (Vector A, Vector B) { return Vector(A.x - B.x, A.y - B.y); } Vector operator * (Vector A, double p) { return Vector(A.x * p, A.y * p); } Vector operator / (Vector A, double p) { return Vector(A.x / p, A.y / p); } const double PI = acos(-1.0); int dcmp(double x) { if (fabs(x) < eps) return 0; else return x < 0 ? -1 : 1; } bool operator == (const Point& a, const Point& b) { return dcmp(a.x - b.x) == 0 && dcmp(a.y - b.y) == 0; } bool operator < (const Point& a, const Point& b) { return a.x < b.x || (a.x == b.x && a.y < b.y); } double torad(double ang) { return ang/180*pi; } double Dot(Vector A, Vector B) {return A.x * B.x + A.y * B.y;} //点积 double Length(Vector A) {return sqrt(Dot(A, A));} //向量的模 double Angle(Vector A, Vector B) {return acos(Dot(A, B) / Length(A) / Length(B));} //向量夹角 double Cross(Vector A, Vector B) {return A.x * B.y - A.y * B.x;} //叉积 double Area2(Point A, Point B, Point C) {return Cross(B - A, C - A);} //有向面积 double angle(Vector v) {return atan2(v.y, v.x);} Point GetLineIntersection(Point P, Vector v, Point Q, Vector w) { Vector u = P - Q; double t = Cross(w, u) / Cross(v, w); return P + v * t; } Vector Rotate(Vector A, double rad) { return Vector(A.x * cos(rad) - A.y * sin(rad), A.x * sin(rad) + A.y * cos(rad)); } double DistanceToLine(Point P, Point A, Point B) { Vector v1 = B - A, v2 = P - A; return fabs(Cross(v1, v2)) / Length(v1); } Vector AngleBisector(Point p, Vector v1, Vector v2){//给定两个向量，求角平分线 double rad = Angle(v1, v2); return Rotate(v1, dcmp(Cross(v1, v2)) * 0.5 * rad); } //求线与x轴的真实角(0<=X<180) double RealAngleWithX(Vector a){ Vector b(1, 0); if (dcmp(Cross(a, b)) == 0) return 0.0; else if (dcmp(Dot(a, b) == 0)) return 90.0; double rad = Angle(a, b); rad = (rad / PI) * 180.0; if (dcmp(a.y) < 0) rad = 180.0 - rad; return rad; } //求直线与圆的交点 int getLineCircleIntersection(Point p, Vector v, Circle c, vector &sol) { double a1 = v.x, b1 = p.x - c.c.x, c1 = v.y, d1 = p.y - c.c.y; double e1 = a1 * a1 + c1 * c1, f1 = 2 * (a1 * b1 + c1 * d1), g1 = b1 * b1 + d1 * d1 - c.r * c.r; double delta = f1 * f1 - 4 * e1 * g1, t; if(dcmp(delta) < 0) return 0; else if(dcmp(delta) == 0){ t = (-f1) / (2 * e1); sol.push_back(p + v * t); return 1; } else{ t = (-f1 + sqrt(delta)) / (2 * e1); sol.push_back(p + v * t); t = (-f1 - sqrt(delta)) / (2 * e1); sol.push_back(p + v * t); return 2; } } //两圆相交 int getCircleCircleIntersection(Circle C1, Circle C2, vector &sol) { double d = Length(C1.c - C2.c); if (dcmp(d) == 0) { if (dcmp(C1.r - C2.r) == 0) return -1; // 重合 return 0; } if (dcmp(C1.r + C2.r - d) < 0) return 0; if (dcmp(fabs(C1.r - C2.r) - d) > 0) return 0; double a = angle(C2.c - C1.c); double da = acos((C1.r * C1.r + d * d - C2.r * C2.r) / (2 * C1.r * d)); Point p1 = C1.point(a - da), p2 = C1.point(a + da); sol.push_back(p1); if(p1 == p2) return 1; sol.push_back(p2); return 2; } //点到圆的切线 int getTangents(Point p, Circle C, Vector *v) { Vector u = C.c - p; double dist = Length(u); if (dist < C.r) return 0; else if (dcmp(dist - C.r) == 0) { v[0] = Rotate(u, PI / 2); return 1; } else { double ang = asin(C.r / dist); v[0] = Rotate(u, -ang); v[1] = Rotate(u, +ang); return 2; } } //两圆公切线 //a[i], b[i]分别是第i条切线在圆A和圆B上的切点 int getCircleTangents(Circle A, Circle B, Point *a, Point *b) { int cnt = 0; if (A.r < B.r) { swap(A, B); swap(a, b); } //圆心距的平方 double 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); double rdiff = A.r - B.r; double rsum = A.r + B.r; double base = angle(B.c - A.c); //重合有无限多条 if (d2 == 0 && dcmp(A.r - B.r) == 0) return -1; //内切 if (dcmp(d2 - rdiff * rdiff) == 0) { a[cnt] = A.point(base); b[cnt] = B.point(base); cnt++; return 1; } //有外公切线 double ang = acos((A.r - B.r) / sqrt(d2)); a[cnt] = A.point(base + ang); b[cnt] = B.point(base + ang); cnt++; a[cnt] = A.point(base - ang); b[cnt] = B.point(base - ang); cnt++; //一条内切线，两条内切线 if (dcmp(d2 - rsum*rsum) == 0) { a[cnt] = A.point(base); b[cnt] = B.point(PI + base); cnt++; } else if (dcmp(d2 - rsum*rsum) > 0) { double ang = acos((A.r + B.r) / sqrt(d2)); a[cnt] = A.point(base + ang); b[cnt] = B.point(base + ang); cnt++; a[cnt] = A.point(base - ang); b[cnt] = B.point(base - ang); cnt++; } return cnt; } //三角形外切圆 Circle CircumscribedCircle(Point p1, Point p2, Point p3) { double Bx = p2.x - p1.x, By = p2.y - p1.y; double Cx = p3.x - p1.x, Cy = p3.y - p1.y; double D = 2 * (Bx * Cy - By * Cx); double cx = (Cy * (Bx * Bx + By * By) - By * (Cx * Cx + Cy * Cy)) / D + p1.x; double cy = (Bx * (Cx * Cx + Cy * Cy) - Cx * (Bx * Bx + By * By)) / D + p1.y; Point p = Point(cx, cy); return Circle(p, Length(p1 - p)); } //三角形内切圆 Circle InscribedCircle(Point p1, Point p2, Point p3) { double a = Length(p2 - p3); double b = Length(p3 - p1); double c = Length(p1 - p2); Point p = (p1 * a + p2 * b + p3 * c) / (a + b + c); return Circle(p, DistanceToLine(p, p1, p2)); } //求经过点p1,与直线(p2, w)相切,半径为r的一组圆 int CircleThroughAPointAndTangentToALineWithRadius(Point p1, Point p2, Vector w, double r, vector &sol) { Circle c1 = Circle(p1, r); double t = r / Length(w); Vector u = Vector(-w.y, w.x); Point p4 = p2 + u * t; int tot = getLineCircleIntersection(p4, w, c1, sol); u = Vector(w.y, -w.x); p4 = p2 + u * t; tot += getLineCircleIntersection(p4, w, c1, sol); return tot; } //给定两个向量，求两向量方向内夹着的圆的圆心。圆与两线均相切，圆的半径已给定 Point Centre_CircleTangentTwoNonParallelLineWithRadius(Point p1, Vector v1, Point p2, Vector v2, double r){ Point p0 = GetLineIntersection(p1, v1, p2, v2); Vector u = AngleBisector(p0, v1, v2); double rad = 0.5 * Angle(v1, v2); double l = r / sin(rad); double t = l / Length(u); return p0 + u * t; } //求与两条不平行的直线都相切的4个圆，圆的半径已给定 int CircleThroughAPointAndTangentALineWithRadius(Point p1, Vector v1, Point p2, Vector v2, double r, Point *sol) { int ans = 0; sol[ans++] = Centre_CircleTangentTwoNonParallelLineWithRadius(p1, v1, p2, v2, r); sol[ans++] = Centre_CircleTangentTwoNonParallelLineWithRadius(p1, v1 * -1, p2, v2, r); sol[ans++] = Centre_CircleTangentTwoNonParallelLineWithRadius(p1, v1, p2, v2 * -1, r); sol[ans++] = Centre_CircleTangentTwoNonParallelLineWithRadius(p1, v1 * -1, p2, v2 * -1, r); return ans; } //求与两个相离的圆均外切的一组圆，三种情况 int CircleTangentToTwoDisjointCirclesWithRadius(Circle c1, Circle c2, double r, Point *sol){ double dis1 = c1.r + r + r + c2.r; double dis2= Length(c1.c - c2.c); if(dcmp(dis1 - dis2) < 0) return 0; Vector u = c2.c - c1.c; double t = (r + c1.r) / Length(u); if(dcmp(dis1 - dis2)==0){ Point p0 = c1.c + u * t; sol[0] = p0; return 1; } double aa = Length(c1.c - c2.c); double bb = r + c1.r, cc = r + c2.r; double rad = acos((aa * aa + bb * bb - cc * cc) / (2 * aa * bb)); Vector w = Rotate(u, rad); Point p0 = c1.c + w * t; sol[0] = p0; w = Rotate(u, -rad); p0 = c1.c + w * t; sol[1] = p0; return 2; } //判断点与圆的位置关系（自己写的） int pointincircle(Point a,Circle o) { double l=Length(o.c-a); if(dcmp(l-o.r)>0) return 1; else if(dcmp(l-o.r)==0) return 0; else if(dcmp(l-o.r)<0) return -1; } int ConvexHull(Point *p, int n, Point* ch) //求凸包 { sort(p, p + n,cmp);//先按照x，再按照y int m = 0; for(int i = 0; i < n; i++) { while(m > 1 && Cross(ch[m-1] - ch[m-2], p[i] - ch[m-2]) <= 0) m--; ch[m++] = p[i]; } int k = m; for(int i = n-2; i >= 0; i--) { while(m > k && Cross(ch[m-1] - ch[m-2], p[i] - ch[m-2]) <= 0) m--; ch[m++] = p[i]; } if(n > 1) m--; return m; } double Polygonarea(Point *p,int n) { double area=0; for(int i=1;i
转载地址：http://xtgsi.baihongyu.com/