1. 点、线、凸边形
1 /******************************************************* 2 二维几何基础 3 【注意】数组下标从1开始。 4 *******************************************************/ 5 6 #include <iostream> 7 #include <cstdio> 8 #include <cstdlib> 9 #include <cmath> 10 #include <iomanip> 11 #include <algorithm> 12 #include <string> 13 #define re register 14 #define il inline 15 #define ll long long 16 #define ld long double 17 using namespace std; 18 const ll MAXN = 1e5+5; 19 const ll INF = 1e9; 20 const ld EPS = 1e-9; 21 22 //点坐标 23 struct POINT 24 { 25 ld x, y; 26 POINT() : x(0), y(0) {} 27 POINT(ld _x, ld _y) : x(_x), y(_y) {} 28 }; 29 //向量 30 typedef POINT VECTOR; 31 32 POINT xy[MAXN]; //顺时针多边形顶点存储 33 POINT xyB[MAXN]; //判断点存储 34 35 36 //符号函数 37 ll sgn(ld x) {return x < -EPS ? -1 : x > EPS;} 38 //向量+向量=向量，点+向量=点 39 VECTOR operator + (VECTOR a, VECTOR b) {return {a.x+b.x,a.y+b.y};} 40 //向量-向量=向量，点-点=向量 41 VECTOR operator - (VECTOR a, VECTOR b) {return {a.x-b.x,a.y-b.y};} 42 //向量*数=向量 43 VECTOR operator * (VECTOR a, ld k) {return {a.x*k,a.y*k};} 44 VECTOR operator * (ld k, VECTOR a) {return {a.x*k,a.y*k};} 45 //向量/数=向量 46 VECTOR operator / (VECTOR a, ld k) {return {a.x/k,a.y/k};} 47 //向量==向量，点==点 48 bool operator == (const VECTOR a, const VECTOR b) {return !sgn(a.x-b.x) && !sgn(a.y-b.y);} 49 //向量偏序关系（先x轴再y轴） 50 bool operator < (const VECTOR a, const VECTOR b) {return !sgn(a.x-b.x) ? a.y < b.y : a.x < b.x;} 51 //向量旋转（逆时针） 52 VECTOR rot(VECTOR a, ld sita) {return {a.x*cos(sita)-a.y*sin(sita),a.x*sin(sita)+a.y*cos(sita)};} 53 //取下端点 54 POINT min(POINT p1, POINT p2) {return p1.y<p2.y ? p1:p2;} 55 //取上端点 56 POINT max(POINT p1, POINT p2) {return p1.y<p2.y ? p2:p1;} 57 //点积 58 ld dot(POINT p1, POINT p2, POINT p) {return (p.x-p1.x)*(p2.x-p1.x)+(p.y-p1.y)*(p2.y-p1.y);} //（三点式：共p1端点） 59 ld dot(VECTOR a, VECTOR b) {return a.x*b.x+a.y*b.y;} //向量式 60 //叉积 61 ld cross(POINT p1, POINT p2, POINT p) {return (p.x-p1.x)*(p2.y-p1.y)-(p.y-p1.y)*(p2.x-p1.x);} //（三点式：共p1端点） 62 ld cross(VECTOR a, VECTOR b) {return a.x*b.y-a.y*b.x;} //向量式（aXb） 63 //向量长度 64 ld vlen(VECTOR a) {return sqrt(dot(a,a));} 65 //向量夹角余弦值 66 ld vcos(VECTOR a, VECTOR b) {return dot(a,b)/(vlen(a)*vlen(b));} 67 //向量夹角正弦值 68 ld vsin(VECTOR a, VECTOR b) {return fabs(cross(a,b))/(vlen(a)*vlen(b));} 69 //向量极角 70 ld vagl(VECTOR a, VECTOR b) {return acos(dot(a,b)/(vlen(a)*vlen(b)));} //向量间 71 ld vagl(VECTOR a) {return acos(a.x/vlen(a));} 72 //求直线斜率【注意】确保斜率存在 73 ld slope(VECTOR a) {return a.y/a.x;} //向量式 74 ld slope(POINT p, POINT q) {return (p.y-q.y)/(p.x-q.x);} //两点式 75 //单位向量【注意】确保非零向量 76 VECTOR unit(VECTOR a) {return a/vlen(a);} 77 //两直线交点 78 POINT intersectline(POINT p, VECTOR v, POINT q, VECTOR w) {return p+v*cross(w,p-q)/cross(v,w);} //（参数式：P=P0+t*v，P0为直线上某一点，v为直线的方向向量） 79 //点在直线上的投影 80 POINT proline(POINT a, POINT b, POINT p) {return a+(b-a)*(dot(b-a,p-a)/dot(b-a,b-a));} 81 //点关于直线的对称点 82 POINT refline(POINT a, POINT b, POINT p) {return proline(a,b,p)*2-p;} 83 //判断两直线是否平行 84 bool parallel(POINT p1, POINT p2, POINT q1, POINT q2) {return !sgn(cross(p2-p1,q2-q1)) && sgn(cross(p1-q1,p2-q1));} 85 //判断两直线重合 86 bool superposition(POINT p1, POINT p2, POINT q1, POINT q2) {return !sgn(cross(p2-p1,q2-q1)) && !sgn(cross(p1-q1,p2-q1));} 87 //点到直线距离 88 ld disline(POINT a, POINT b, POINT p) {return fabs(cross(b-a,p-a))/vlen(b-a);} //不取绝对值得到的是有向距离 89 //点到线段距离（两种情况：点的投影在线段上，则为垂直距离；点的投影不在线段上，则为到两端距离的最小值） 90 ld disseg(POINT a, POINT b, POINT p) {return a==b ? vlen(p-a):dot(b-a,p-a)<0 ? vlen(p-a):dot(b-a,p-b)>0 ? vlen(p-b):fabs(cross(b-a,p-a))/vlen(b-a);} 91 //线段相交判断（严格） 92 bool intersectseg(POINT p1, POINT p2, POINT q1, POINT q2) {return cross(p1-q1,q2-q1)*cross(p2-q1,q2-q1)<0 && cross(q1-p1,p2-p1)*cross(q2-p1,p2-p1)<0;} 93 //判断点p(x,y)是否在线段p1p2上（两种情况：1.包括端点；2.不包含端点） 94 bool onseg(POINT p1, POINT p2, POINT p) {return !sgn(cross(p1,p2,p)) && sgn(dot(p,p1,p2))<=0;} //包含端点 95 bool onseg_strict(POINT p1, POINT p2, POINT p) {return !sgn(cross(p1,p2,p)) && sgn(dot(p,p1,p2))<0;} //不包含端点 96 97 //点射法判断点是否在多边形内部（边界也算在内部） 98 //复杂度O(n)（不论凹凸） 99 bool inpolygon_dot(ld x, ld y, ll n) 100 { 101 ll cnt = 0; 102 for(re ll i = 0; i < n; ++i) 103 { 104 POINT p1 = min(xy[i+1], xy[(i+1)%n+1]); //取下端点 105 POINT p2 = max(xy[i+1], xy[(i+1)%n+1]); //取上端点 106 //特判点在线段上 107 if(onseg(p1,p2,{x,y})) return true; 108 //从点(x,y)向x反方向引出射线 109 //计算点射到的多边形的边数边数 110 if(sgn(p1.y-y)<0 && sgn(y-p2.y) <=0 && sgn(cross(p1,p2,{x,y}))>0) ++cnt; 111 } 112 if(cnt%2) return true; 113 else return false; 114 } 115 116 //二分法判断点是否在多边形内部（不包括边界） 117 //复杂度O(logn)（要求凸多边形） 118 bool inpolygon_bis(ld x, ld y, ll n) 119 { 120 POINT p = {x,y}; 121 ll l = 2, r = n; 122 //判断是否在幅角最小和最大的两条边上 123 if(onseg(xy[1],xy[l],p) || onseg(xy[1],xy[n],p)) return false; 124 //二分法判断是否在多边形内部 125 while(l < r) 126 { 127 ll mid = (l+r)>>1; //【注意】如此取中点最终：r==l或r==l-1（只有此两种情况） 128 ll d = sgn(cross(xy[mid]-xy[1],p-xy[1])); 129 if(d < 0) l = mid+1; 130 else if(d > 0) r = mid-1; 131 else 132 { 133 if(onseg_strict(xy[1],xy[mid],p)) 134 { 135 return true; 136 } 137 else 138 { 139 return false; 140 } 141 } 142 } 143 //判断在最终二分得到的对角线两端的边是否都满足严格在内的条件 144 if(l >= r && (sgn(cross(xy[l]-xy[l-1],p-xy[l-1]))>=0 || sgn(cross(xy[l%n+1]-xy[l],p-xy[l]))>=0)) 145 { 146 return false; 147 } 148 return true; 149 } 150 151 //计算多边形面积（凹凸均可） 152 ld polygonarea(ll n) 153 { 154 ld aera = 0; 155 for(re ll i = 0; i < n; ++i) 156 { 157 aera += cross(xy[i+1], xy[(i+1)%n+1]); 158 } 159 //计算出来的aera为有向面积 160 return fabs(aera)/2; 161 } 162 163 int main() 164 { 165 std::ios::sync_with_stdio(false); 166 //判断直线平行/重合/相交 167 ll n; 168 cin >> n; 169 cout << "INTERSECTING LINES OUTPUT" << endl; 170 for(re ll i = 1; i <= n; ++i) 171 { 172 POINT p1, p2, p3, p4; 173 cin >> p1.x >> p1.y >> p2.x >> p2.y; 174 cin >> p3.x >> p3.y >> p4.x >> p4.y; 175 if(parallel(p1,p2,p3,p4)) cout << "NONE" << endl; 176 else if(superposition(p1,p2,p3,p4)) cout << "LINE" << endl; 177 else 178 { 179 POINT q = intersectline(p1,p2-p1,p3,p4-p3); 180 cout << "POINT " << fixed << setprecision(2) << q.x << " " << q.y << endl; 181 } 182 } 183 cout << "END OF OUTPUT" << endl; 184 return 0; 185 /* 186 //计算多边形的面积 187 ll n; 188 while(true) 189 { 190 std::cin >> n; 191 if(!n) break; 192 for(re ll i = 1; i <= n; ++i) 193 { 194 std::cin >> xy[i].x; 195 std::cin >> xy[i].y; 196 } 197 std::cout << std::fixed << std::setprecision(1) << polygonarea(n) << std::endl; 198 } 199 return 0; 200 */ 201 /* 202 //判断点是否在多边形内（不包括边界） 203 ll n; 204 std::cin >> n; 205 for(re ll i = 1; i <= n; ++i) 206 { 207 std::cin >> xy[i].x; 208 std::cin >> xy[i].y; 209 } 210 ll m; 211 cin >> m; 212 for(re ll i = 1; i <= m; ++i) 213 { 214 cin >> xyB[i].x; 215 cin >> xyB[i].y; 216 } 217 for(re ll i = 1; i <= m; ++i) 218 { 219 if(!inpolygon_bis(xyB[i].x, xyB[i].y, n)) 220 { 221 printf("NO "); 222 return 0; 223 } 224 } 225 printf("YES "); 226 return 0; 227 */ 228 /* 229 //判断点是否在多边形内（包括边界） 230 ll N, M, k = 0; 231 while(!(std::cin >> N).eof()) 232 { 233 if(!N) break; 234 std::cin >> M; 235 if(k) printf(" "); 236 ++k; 237 for(re ll i = 1; i <= N; ++i) 238 { 239 std::cin >> xy[i].x; 240 std::cin >> xy[i].y; 241 //printf("(%f,%f) ", xy[i].x, xy[i].y); 242 } 243 printf("Problem %lld: ", k); 244 for(re ll i = 1; i <= M; ++i) 245 { 246 ld x, y; 247 std::cin >> x; 248 std::cin >> y; 249 //printf("(%f,%f) ", x, y); 250 if(inpolygon(x, y, N)) 251 { 252 printf("Within "); 253 } 254 else 255 { 256 printf("Outside "); 257 } 258 } 259 } 260 return 0; 261 */ 262 }
2. 圆
1 /***************************************************************** 2 圆 3 【注意】三角与反三角最多用一次（损失精度大）。 4 *****************************************************************/ 5 6 #include <bits/stdc++.h> 7 #include <iomanip> 8 #include <stack> 9 #include <cstdio> 10 #include <cmath> 11 #include <algorithm> 12 #define re register 13 #define il inline 14 #define ll long long 15 #define ld long double 16 using namespace std; 17 const ld PI = 3.14159265358979323846; 18 const ll MAXN = 5e5+5; 19 const ld INF = 1e9; 20 const ld EPS = 1e-10; 21 22 //点坐标 23 struct POINT 24 { 25 ld x, y; 26 POINT() : x(0), y(0) {} 27 POINT(ld _x, ld _y) : x(_x), y(_y) {} 28 }; 29 //向量 30 typedef POINT VECTOR; 31 32 //符号函数 33 ll sgn(ld x) {return x < -EPS ? -1 : x > EPS;} 34 //向量+向量=向量，点+向量=点 35 VECTOR operator + (VECTOR a, VECTOR b) {return {a.x+b.x,a.y+b.y};} 36 //向量-向量=向量，点-点=向量 37 VECTOR operator - (VECTOR a, VECTOR b) {return {a.x-b.x,a.y-b.y};} 38 //向量*数=向量 39 VECTOR operator * (VECTOR a, ld k) {return {a.x*k,a.y*k};} 40 VECTOR operator * (ld k, VECTOR a) {return {a.x*k,a.y*k};} 41 //向量/数=向量 42 VECTOR operator / (VECTOR a, ld k) {return {a.x/k,a.y/k};} 43 //向量==向量，点==点 44 bool operator == (const VECTOR a, const VECTOR b) {return !sgn(a.x-b.x) && !sgn(a.y-b.y);} 45 //向量偏序关系（先x轴再y轴） 46 bool operator < (const VECTOR a, const VECTOR b) {return !sgn(a.x-b.x) ? a.y < b.y : a.x < b.x;} 47 //取下端点 48 POINT min(POINT p1, POINT p2) {return p1.y<p2.y ? p1:p2;} 49 //取上端点 50 POINT max(POINT p1, POINT p2) {return p1.y<p2.y ? p2:p1;} 51 //点积 52 ld dot(POINT p1, POINT p2, POINT p) {return (p.x-p1.x)*(p2.x-p1.x)+(p.y-p1.y)*(p2.y-p1.y);} //（三点式：共p1端点） 53 ld dot(VECTOR a, VECTOR b) {return a.x*b.x+a.y*b.y;} //向量式 54 //叉积 55 ld cross(POINT p1, POINT p2, POINT p) {return (p.x-p1.x)*(p2.y-p1.y)-(p.y-p1.y)*(p2.x-p1.x);} //（三点式：共p1端点） 56 ld cross(VECTOR a, VECTOR b) {return a.x*b.y-a.y*b.x;} //向量式（aXb） 57 //向量长度 58 ld vlen(VECTOR a) {return sqrt(dot(a,a));} 59 //向量旋转（逆时针） 60 VECTOR rot(VECTOR a, ld sita) {return {a.x*cos(sita)-a.y*sin(sita),a.x*sin(sita)+a.y*cos(sita)};} 61 //取单位向量 62 VECTOR unit(VECTOR a) {return a/vlen(a);} 63 //取正交单位向量 64 VECTOR norm(VECTOR a) {return VECTOR(-a.y,a.x)/vlen(a);} 65 //向量夹角余弦值 66 ld vcos(VECTOR a, VECTOR b) {return dot(a,b)/(vlen(a)*vlen(b));} 67 //向量夹角正弦值 68 ld vsin(VECTOR a, VECTOR b) {return fabs(cross(a,b))/(vlen(a)*vlen(b));} 69 //向量极角 70 //【注意】使用系统反三角函数返回nan是因为超出了定义域（大多数由于精度造成） 71 //【注意】故使用系统反三角函数需处理边界 72 ld vagl(VECTOR a, VECTOR b) 73 { 74 ld d = dot(a,b)/(vlen(a)*vlen(b)); 75 if(sgn(d-1) == 0) return 0; 76 else if(sgn(d+1) == 0) return PI; 77 else return acos(d); 78 } //向量间[0,PI] 79 ld vagl(VECTOR a) {return atan2(a.y, a.x);} //单向量[-PI,PI] 80 81 //直线 82 struct LINE 83 { 84 POINT p; //直线上定点 85 VECTOR v; //直线方向向量 86 LINE() {} 87 //点向式 88 LINE(POINT _p, VECTOR _v):p(_p), v(_v) {} 89 //点角式 90 LINE(POINT _p, ld a):p(_p), v(VECTOR(cos(a), sin(a))) {} 91 //通过参数t唯一确定直线上点坐标（点向式：P = p + t*v） 92 POINT point(ld t) {return p + t*v;} 93 //平移向量r 94 LINE trans(VECTOR r) {return LINE(p+r,v);} 95 }; 96 97 //点到直线距离 98 ld PLdist(POINT p, LINE L) {return fabs(sin(vagl(p-L.p,L.v))*vlen(p-L.p));} 99 100 //圆 101 struct CIRCLE 102 { 103 POINT c; //圆心坐标 104 ld r; //半径 105 CIRCLE() {} 106 CIRCLE(POINT _c, ld _r):c(_c),r(_r) {} 107 //通过圆心角唯一确定圆上点坐标 108 POINT point(ld a) //a为圆心角 109 { 110 return POINT(c.x+cos(a)*r, c.y+sin(a)*r); 111 } 112 }; 113 114 //两点中垂线 115 //返回中垂线的个数 116 ll PerpendicularBisector(POINT p1, POINT p2, LINE &L) 117 { 118 if(p1 == p2) return 0; 119 L.p = (p1+p2)/2; 120 L.v = rot(p2-p1, PI/2); 121 return 1; 122 } 123 124 //两直线交点 125 POINT LLitesct; 126 127 //直线与直线的交点 128 //返回交点个数 129 ll LineLineIntersection(LINE L1, LINE L2, POINT &p) 130 { 131 if(sgn(cross(L1.v, L2.v)) == 0) 132 { 133 if(sgn(cross(L2.p-L1.p,L2.v) == 0)) return -1; //两直线重合 134 else return 0; //两直线平行 135 } 136 //以L1为准 137 ld t1 = -cross(L2.p-L1.p, L2.v)/cross(L2.v, L1.v); 138 p = L1.p + t1*L1.v; 139 //以L2为准 140 //ld t2 = -cross(L1.p-L2.p, L1.v)/cross(L1.v, L2.v); 141 //p = L2.p + t2*L2.v; 142 return 1; 143 } 144 145 //直线与圆交点数组 146 vector <POINT> LCitesct; 147 148 //直线与圆的交点 149 //解方程：(at+b)^2+(ct+d)^2 = r^2; 150 //化简得：et^2+ft+g = 0 151 //返回交点个数 152 ll LineCircleIntersection(LINE L, CIRCLE C, ld &t1, ld &t2, vector <POINT>& sol) 153 { 154 ld a = L.v.x, b = L.p.x-C.c.x, c = L.v.y, d = L.p.y-C.c.y; 155 ld e = a*a + c*c, f = 2*(a*b + c*d), g = b*b + d*d - C.r*C.r; 156 ld delta = f*f - 4*e*g; 157 if(sgn(delta) < 0) //相离 158 return 0; 159 else if(sgn(delta) == 0) //相切 160 { 161 t1 = t2 = -f/(2*e); 162 sol.push_back(L.point(t1)); 163 return 1; 164 } 165 else 166 { 167 t1 = (-f - sqrt(delta))/(2*e); 168 t2 = (-f + sqrt(delta))/(2*e); 169 sol.push_back(L.point(t1)); 170 sol.push_back(L.point(t2)); 171 return 2; 172 } 173 } 174 ll LineCircleIntersection(LINE L, CIRCLE C, POINT *p) 175 { 176 ld t1, t2; 177 ld a = L.v.x, b = L.p.x-C.c.x, c = L.v.y, d = L.p.y-C.c.y; 178 ld e = a*a + c*c, f = 2*(a*b + c*d), g = b*b + d*d - C.r*C.r; 179 ld delta = f*f - 4*e*g; 180 if(sgn(delta) < 0) //相离 181 return 0; 182 else if(sgn(delta) == 0) //相切 183 { 184 t1 = t2 = -f/(2*e); 185 p[0] = L.point(t1); 186 return 1; 187 } 188 else 189 { 190 t1 = (-f - sqrt(delta))/(2*e); 191 t2 = (-f + sqrt(delta))/(2*e); 192 p[0] = L.point(t1); 193 p[1] = L.point(t2); 194 return 2; 195 } 196 } 197 198 //圆与圆的交点数组 199 vector <POINT> CCitesct; 200 201 //两圆交点 202 ll CircleCircleIntersection(CIRCLE C1, CIRCLE C2, vector <POINT>& sol) 203 { 204 ld d = vlen(C1.c-C2.c); 205 if(sgn(d) == 0) 206 { 207 if(sgn(C1.r-C2.r) == 0) return -1; //两圆重合 208 else return 0; //同心相含 209 } 210 if(sgn(C1.r+C2.r-d) < 0) return 0; //相离 211 if(sgn(fabs(C1.r-C2.r)-d) > 0) return 0; //异心相含 212 213 //设两圆交点为p1,p2 214 ld a = vagl(C2.c - C1.c); //向量C1C2的极角 215 ld b = acos((C1.r*C1.r + d*d - C2.r*C2.r)/(2*C1.r*d)); //向量C1C2到C1P1（或C1P2）的夹角 216 POINT p1 = C1.point(a-b); 217 POINT p2 = C1.point(a+b); 218 sol.push_back(p1); 219 if(p1 == p2) return 1; //相切 220 sol.push_back(p2); //相交 221 return 2; 222 } 223 224 //圆切线方向向量 225 VECTOR LCtgt[2]; 226 227 //过定点作圆的切线 228 ll Tangents(POINT p, CIRCLE C, VECTOR *v) 229 { 230 VECTOR u = C.c - p; 231 ld dist = vlen(u); 232 if(dist < C.r) return 0; //点在圆内，没有切线 233 else if(sgn(dist-C.r) == 0) //点在圆上，一条切线 234 { 235 v[0] = rot(u, PI/2); 236 return 1; 237 } 238 else //点在圆外，两条切线 239 { 240 ld a = asin(C.r/dist); 241 v[0] = rot(u, -a); 242 v[1] = rot(u, +a); 243 return 2; 244 } 245 } 246 247 //两圆公切点（分别在两圆上） 248 POINT CCtgt1[4], CCtgt2[4]; 249 250 //两圆公切线 251 ll Tangents(CIRCLE C1, CIRCLE C2, POINT *a, POINT *b) 252 { 253 ll cnt = 0; 254 if(C1.r < C2.r) {swap(C1, C2); swap(a, b);} 255 ld d = dot(C1.c-C2.c, C1.c-C2.c); 256 ld rdiff = C1.r - C2.r; 257 ld rsum = C1.r + C2.r; 258 if(d < rdiff*rsum) return 0; //相含 259 ld base = atan2(C2.c.y-C1.c.y, C2.c.x-C1.c.x); //向量C1C2极角[-PI,PI] 260 if(sgn(d) == 0 && sgn(C1.r-C2.r) == 0) return -1; //重合 261 if(sgn(d-rdiff*rdiff) == 0) //内切 262 { 263 a[cnt] = C1.point(base); 264 b[cnt++] = C2.point(base); 265 return 1; 266 } 267 ld ang = acos((C1.r-C2.r)/sqrt(d)); //夹角[0,PI] 268 a[cnt] = C1.point(base+ang); 269 b[cnt++] = C2.point(base+ang); 270 a[cnt] = C1.point(base-ang); 271 b[cnt++] = C2.point(base-ang); 272 if(sgn(d-rsum*rsum)) //外切 273 { 274 a[cnt] = C1.point(base); 275 b[cnt++] = C2.point(base+PI); 276 } 277 else if(d > rsum*rsum) //相离 278 { 279 ld ang = acos((C1.r+C2.r)/sqrt(d)); 280 a[cnt] = C1.point(base+ang); 281 b[cnt++] = C2.point(base+ang); 282 a[cnt] = C1.point(base-ang); 283 b[cnt++] = C2.point(base-ang); 284 } 285 return cnt; //切点个数 286 } 287 288 //三点外接圆 289 //返回外接圆的个数 290 ll CircumscribedCircle(POINT a, POINT b, POINT c, CIRCLE &C) 291 { 292 if(sgn(cross(a-b, b-c)) == 0) return 0; //三点共线 293 LINE L1, L2; 294 PerpendicularBisector(a, b, L1); 295 PerpendicularBisector(b, c, L2); 296 POINT p; 297 LineLineIntersection(L1, L2, p); 298 C.c = p; 299 C.r = vlen(p-a); 300 return 1; 301 } 302 303 //三点内接圆 304 //返回内接圆的个数 305 ll InscribedCircle(POINT a, POINT b, POINT c, CIRCLE &C) 306 { 307 if(sgn(cross(a-b, b-c)) == 0) return 0; //三点共线 308 LINE L1 = LINE(a, (vagl(b-a)+vagl(c-a))/2), L2 = LINE(b, (vagl(a-b)+vagl(c-b))/2); 309 POINT p; 310 LineLineIntersection(L1, L2, p); 311 C.c = p; 312 C.r = PLdist(p,LINE(a,b-a)); 313 return 1; 314 } 315 316 317 int main() 318 { 319 //UVA-12304 320 //ios::sync_with_stdio(false); 321 string str; 322 while(cin >> str) 323 { 324 if(str == "CircumscribedCircle") 325 { 326 POINT a, b, c; 327 cin >> a.x >> a.y >> b.x >> b.y >> c.x >> c.y; 328 CIRCLE C; 329 CircumscribedCircle(a,b,c,C); 330 cout << fixed << "(" << C.c.x << "," << C.c.y << "," << C.r << ")" << endl; 331 } 332 else if(str == "InscribedCircle") 333 { 334 POINT a, b, c; 335 cin >> a.x >> a.y >> b.x >> b.y >> c.x >> c.y; 336 CIRCLE C; 337 InscribedCircle(a,b,c,C); 338 cout << fixed << "(" << C.c.x << "," << C.c.y << "," << C.r << ")" << endl; 339 } 340 //求已知圆的切线的极角 341 else if(str == "TangentLineThroughPoint") 342 { 343 CIRCLE C; 344 POINT p; 345 cin >> C.c.x >> C.c.y >> C.r >> p.x >> p.y; 346 VECTOR v[2]; 347 ll cnt = Tangents(p,C,v); 348 cout << "["; 349 if(cnt == 1) 350 { 351 ld r = vagl(v[0]); 352 ld a = (sgn(r) < 0 ? PI+r : sgn(r-PI) == 0 ? 0 : r)*180/PI; 353 cout << fixed << a; 354 } 355 else if(cnt == 2) 356 { 357 ld mxr = vagl(v[0]); 358 ld mir = vagl(v[1]); 359 ld mxa = (sgn(mxr) < 0 ? PI+mxr : sgn(mxr-PI) == 0 ? 0 : mxr)*180/PI; 360 ld mia = (sgn(mir) < 0 ? PI+mir : sgn(mir-PI) == 0 ? 0 : mir)*180/PI; 361 if(mxa < mia) swap(mxa,mia); 362 cout << fixed << mia << "," << mxa; 363 } 364 cout << "]" << endl; 365 } 366 //求过定点且与已知直线相切的圆的圆心 367 else if(str == "CircleThroughAPointAndTangentToALineWithRadius") 368 { 369 CIRCLE P; 370 POINT p1, p2; 371 cin >> P.c.x >> P.c.y >> p1.x >> p1.y >> p2.x >> p2.y >> P.r; 372 VECTOR v = rot(p2-p1,PI/2); 373 ld t = P.r/vlen(v); 374 LINE L1 = LINE(p1+t*v,p2-p1), L2 = LINE(p1-t*v,p2-p1); 375 ll cnt = 0; 376 ld t1, t2; 377 vector <POINT> C; 378 cnt += LineCircleIntersection(L1,P,t1,t2,C); 379 cnt += LineCircleIntersection(L2,P,t1,t2,C); 380 sort(C.begin(),C.end()); 381 cout << "["; 382 if(cnt == 1) cout << fixed << "(" << C[0].x << "," << C[0].y << ")"; 383 else if(cnt == 2) cout << fixed << "(" << C[0].x << "," << C[0].y << ")" << "," << "(" << C[1].x << "," << C[1].y << ")"; 384 cout << "]" << endl; 385 } 386 //求给定半径且与两条不相交直线相切的圆的圆心 387 else if(str == "CircleTangentToTwoLinesWithRadius") 388 { 389 POINT p1, p2, p3, p4; 390 ld r; 391 cin >> p1.x >> p1.y >> p2.x >> p2.y >> p3.x >> p3.y >> p4.x >> p4.y >> r; 392 LINE L1 = LINE(p1,p2-p1), L2 = LINE(p3,p4-p3); 393 VECTOR v1 = norm(p2-p1), v2 = norm(p4-p3); 394 LINE L11 = L1.trans(r*v1), L12 = L1.trans(-r*v1); 395 LINE L21 = L2.trans(r*v2), L22 = L2.trans(-r*v2); 396 vector <POINT> C; 397 POINT p; 398 LineLineIntersection(L11,L21,p); 399 C.push_back(p); 400 LineLineIntersection(L11,L22,p); 401 C.push_back(p); 402 LineLineIntersection(L12,L21,p); 403 C.push_back(p); 404 LineLineIntersection(L12,L22,p); 405 C.push_back(p); 406 sort(C.begin(),C.end()); 407 cout << "["; 408 cout << fixed << "(" << C[0].x << "," << C[0].y << ")" << "," << "(" << C[1].x << "," << C[1].y << ")" << ","; 409 cout << fixed << "(" << C[2].x << "," << C[2].y << ")" << "," << "(" << C[3].x << "," << C[3].y << ")"; 410 cout << "]" << endl; 411 412 } 413 else if(str == "CircleTangentToTwoDisjointCirclesWithRadius") 414 { 415 CIRCLE C1, C2; 416 ld r; 417 cin >> C1.c.x >> C1.c.y >> C1.r >> C2.c.x >> C2.c.y >> C2.r >> r; 418 C1.r += r, C2.r += r; 419 vector <POINT> p; 420 ll cnt = CircleCircleIntersection(C1,C2,p); 421 sort(p.begin(),p.end()); 422 cout << "["; 423 if(cnt == 1) 424 { 425 cout << fixed << "(" << p[0].x << "," << p[0].y << ")"; 426 } 427 else if(cnt == 2) 428 { 429 cout << fixed << "(" << p[0].x << "," << p[0].y << ")" << "," << "(" << p[1].x << "," << p[1].y << ")"; 430 } 431 cout << "]" << endl; 432 } 433 } 434 return 0; 435 }
3. 半平面
1 /***************************************************************** 2 半平面交 3 【注意】三角与反三角最多用一次（损失精度大）。 4 【注意】数组下标一律从1开始。 5 *****************************************************************/ 6 7 //#include <bits/stdc++.h> 8 #include <iostream> 9 #include <cstring> 10 #include <cstdlib> 11 #include <iomanip> 12 #include <stack> 13 #include <cstdio> 14 #include <cmath> 15 #include <algorithm> 16 #define re register 17 #define il inline 18 #define ll int 19 #define ld double 20 using namespace std; 21 const ld PI = 3.14159265358979323846; 22 const ll MAXN = 5e5+5; 23 const ld INF = 1e9; 24 const ld EPS = 1e-8; 25 26 //点坐标 27 struct POINT 28 { 29 ld x, y; 30 POINT() : x(0), y(0) {} 31 POINT(ld _x, ld _y) : x(_x), y(_y) {} 32 }; 33 //向量 34 typedef POINT VECTOR; 35 36 //符号函数 37 ll sgn(ld x) {return x < -EPS ? -1 : x > EPS;} 38 //向量+向量=向量，点+向量=点 39 VECTOR operator + (VECTOR a, VECTOR b) {return {a.x+b.x,a.y+b.y};} 40 //向量-向量=向量，点-点=向量 41 VECTOR operator - (VECTOR a, VECTOR b) {return {a.x-b.x,a.y-b.y};} 42 //向量*数=向量 43 VECTOR operator * (VECTOR a, ld k) {return {a.x*k,a.y*k};} 44 VECTOR operator * (ld k, VECTOR a) {return {a.x*k,a.y*k};} 45 //向量/数=向量 46 VECTOR operator / (VECTOR a, ld k) {return {a.x/k,a.y/k};} 47 //向量==向量，点==点 48 bool operator == (const VECTOR a, const VECTOR b) {return !sgn(a.x-b.x) && !sgn(a.y-b.y);} 49 //向量偏序关系（先x轴再y轴） 50 bool operator < (const VECTOR a, const VECTOR b) {return !sgn(a.x-b.x) ? a.y < b.y : a.x < b.x;} 51 //取下端点 52 POINT min(POINT p1, POINT p2) {return p1.y<p2.y ? p1:p2;} 53 //取上端点 54 POINT max(POINT p1, POINT p2) {return p1.y<p2.y ? p2:p1;} 55 //点积 56 ld dot(POINT p1, POINT p2, POINT p) {return (p.x-p1.x)*(p2.x-p1.x)+(p.y-p1.y)*(p2.y-p1.y);} //（三点式：共p1端点） 57 ld dot(VECTOR a, VECTOR b) {return a.x*b.x+a.y*b.y;} //向量式 58 //叉积 59 ld cross(POINT p1, POINT p2, POINT p) {return (p.x-p1.x)*(p2.y-p1.y)-(p.y-p1.y)*(p2.x-p1.x);} //（三点式：共p1端点） 60 ld cross(VECTOR a, VECTOR b) {return a.x*b.y-a.y*b.x;} //向量式（aXb） 61 //向量长度 62 ld vlen(VECTOR a) {return sqrt(dot(a,a));} 63 //向量旋转（逆时针） 64 VECTOR rot(VECTOR a, ld sita) {return {a.x*cos(sita)-a.y*sin(sita),a.x*sin(sita)+a.y*cos(sita)};} 65 //取单位向量 66 VECTOR unit(VECTOR a) {return a/vlen(a);} 67 //取正交单位向量 68 VECTOR norm(VECTOR a) {return VECTOR(-a.y,a.x)/vlen(a);} 69 //向量夹角余弦值 70 ld vcos(VECTOR a, VECTOR b) {return dot(a,b)/(vlen(a)*vlen(b));} 71 //向量夹角正弦值 72 ld vsin(VECTOR a, VECTOR b) {return fabs(cross(a,b))/(vlen(a)*vlen(b));} 73 //向量极角 74 //【注意】使用系统反三角函数返回nan是因为超出了定义域（大多数由于精度造成） 75 //【注意】故使用系统反三角函数需处理边界 76 ld vagl(VECTOR a, VECTOR b) 77 { 78 ld d = dot(a,b)/(vlen(a)*vlen(b)); 79 if(sgn(d-1) == 0) return 0; 80 else if(sgn(d+1) == 0) return PI; 81 else return acos(d); 82 } //向量间[0,PI] 83 ld vagl(VECTOR a) {return atan2(a.y, a.x);} //单向量[-PI,PI] 84 85 //凸包 86 ll vexcnt = 0; //凸包顶点数 87 POINT xy[MAXN], xyt[MAXN]; //多边形顶点存储，临时多边形顶点存储 88 POINT convex[MAXN]; //凸包顶点数组 89 90 //求凸包周长 91 ld convexperimeter(POINT *poly, ll n) 92 { 93 ld ans = 0; 94 for(re ll i = 1; i <= n; ++i) 95 { 96 ans += vlen(poly[i]-poly[i%n+1]); 97 } 98 return ans; 99 } 100 //计算多边形面积 101 ld polygonarea(POINT *poly, ll n) 102 { 103 ld aera = 0; 104 for(re ll i = 0; i < n; ++i) 105 { 106 aera += cross(poly[i+1], poly[(i+1)%n+1]); 107 } 108 return fabs(aera)/2; //不加绝对值为有向面积 109 } 110 111 //andrew算法求凸包（凸包数组正序为逆时针） 112 //1.严格凸包（没有重复点和三点共线）2.非严格凸包（允许重复点和三点共线） 113 void andrew(ll n) 114 { 115 vexcnt = 0; //初始化数据 116 memcpy(xyt,xy,sizeof(xy)); //备份原来顶点集 117 sort(xyt+1,xyt+n+1); //排序后xyt[1]和xyt[n]一定在凸包中 118 //求解凸包顶点集 119 //正向扫描（下凸包） 120 convex[++vexcnt] = xyt[1]; //xyt[1]一定在凸包中 121 ll j = 2; 122 while(j <= n) 123 { 124 if(vexcnt <= 1) //因为xyt[2]不一定在凸包集中 125 { 126 convex[++vexcnt] = xyt[j++]; 127 } 128 else 129 { 130 //取前面两个点 131 //严格凸包：sgn(cross(convex[vexcnt]-convex[vexcnt-1],xyt[j]-convex[vexcnt]))>0 132 //非严格凸包：sgn(cross(convex[vexcnt]-convex[vexcnt-1],xyt[j]-convex[vexcnt]))>=0 133 if(sgn(cross(convex[vexcnt]-convex[vexcnt-1],xyt[j]-convex[vexcnt]))>=0) convex[++vexcnt] = xyt[j++]; 134 else --vexcnt; 135 } 136 } 137 //反向扫描（上凸包） 138 //至少包含了xyt[1]和xyt[n] 139 ll record = vexcnt; //记录下凸包的结果 140 ll k = n-1; 141 while(k >= 1) 142 { 143 if(vexcnt <= record) 144 { 145 convex[++vexcnt] = xyt[k--]; 146 } 147 else 148 { 149 //取前面两个点 150 //严格凸包：sgn(cross(convex[vexcnt]-convex[vexcnt-1],xyt[k]-convex[vexcnt]))>0 151 //非严格凸包：sgn(cross(convex[vexcnt]-convex[vexcnt-1],xyt[k]-convex[vexcnt]))>=0 152 if(sgn(cross(convex[vexcnt]-convex[vexcnt-1],xyt[k]-convex[vexcnt]))>0) convex[++vexcnt] = xyt[k--]; 153 else --vexcnt; 154 } 155 } 156 while(convex[--vexcnt]==convex[1]); //因为和xyt[1]相等的点都在下凸包中了 157 //for(re ll i = 1; i <= vexcnt; ++i) cout << "(" << convex[i].x << "," << convex[i].y << ")" << endl; 158 return; 159 } 160 161 //直线 162 struct LINE 163 { 164 POINT p; //直线上定点 165 VECTOR v; //直线方向向量 166 LINE() {} 167 //点向式 168 LINE(POINT _p, VECTOR _v):p(_p), v(_v) {} 169 //点角式 170 LINE(POINT _p, ld a):p(_p), v(VECTOR(cos(a), sin(a))) {} 171 //通过参数t唯一确定直线上点坐标（点向式：P = p + t*v） 172 POINT point(ld t) {return p + t*v;} 173 //平移向量r 174 LINE trans(VECTOR r) {return LINE(p+r,v);} 175 }; 176 177 //有向直线 178 struct DIRLINE 179 { 180 POINT p; //定点 181 VECTOR v; //方向向量 182 ld a; //极角 183 DIRLINE() {} 184 DIRLINE(POINT _p, VECTOR _v):p(_p),v(_v),a(atan2(v.y,v.x)) {} 185 bool operator < (const DIRLINE &DL) const {return a < DL.a;} 186 }; 187 188 //点到直线距离 189 ld PLdist(POINT p, LINE L) {return fabs(sin(vagl(p-L.p,L.v))*vlen(p-L.p));} 190 //点在有向直线左边 191 bool OnLeft(DIRLINE DL, POINT p) {return sgn(cross(DL.v,p-DL.p)) > 0;} 192 //点在有向直线右边 193 bool OnRight(DIRLINE DL, POINT p) {return sgn(cross(DL.v,p-DL.p)) < 0;} 194 195 //两点中垂线 196 //返回中垂线的个数 197 ll PerpendicularBisector(POINT p1, POINT p2, LINE &L) 198 { 199 if(p1 == p2) return 0; 200 L.p = (p1+p2)/2; 201 L.v = rot(p2-p1, PI/2); 202 return 1; 203 } 204 205 //两直线交点 206 POINT LLitesct; 207 208 //直线与直线的交点 209 //返回交点个数 210 ll LineLineIntersection(LINE L1, LINE L2, POINT &p) 211 { 212 if(sgn(cross(L1.v, L2.v)) == 0) 213 { 214 if(sgn(cross(L2.p-L1.p,L2.v)) == 0) return -1; //两直线重合 215 else return 0; //两直线平行 216 } 217 //以L1为准 218 ld t1 = -cross(L2.p-L1.p, L2.v)/cross(L2.v, L1.v); 219 p = L1.p + t1*L1.v; 220 //以L2为准 221 //ld t2 = -cross(L1.p-L2.p, L1.v)/cross(L1.v, L2.v); 222 //p = L2.p + t2*L2.v; 223 return 1; 224 } 225 ll LineLineIntersection(DIRLINE DL1, DIRLINE DL2, POINT &p) 226 { 227 if(sgn(cross(DL1.v, DL2.v)) == 0) 228 { 229 if(sgn(cross(DL2.p-DL1.p,DL2.v)) == 0) return -1; //两直线重合 230 else return 0; //两直线平行 231 } 232 //以L1为准 233 ld t1 = -cross(DL2.p-DL1.p, DL2.v)/cross(DL2.v, DL1.v); 234 p = DL1.p + t1*DL1.v; 235 //以L2为准 236 //ld t2 = -cross(L1.p-L2.p, L1.v)/cross(L1.v, L2.v); 237 //p = L2.p + t2*L2.v; 238 return 1; 239 } 240 241 //半平面交 242 //返回凸包顶点数 243 //【注意】若半平面交可能是无界区域，则需在外面加一个坐标很大的包围框 244 ll HalfPlaneIntersection(DIRLINE *DL, ll n, POINT *poly) 245 { 246 sort(DL+1,DL+n+1); 247 ll first, last; //双端队列的首尾元素下标 248 DIRLINE *q = new DIRLINE[n+1]; //双端队列 249 POINT *p = new POINT[n+1]; //p[i]为q[i]和q[i+1]的交点 250 first = last = 1; 251 q[1] = DL[1]; //初始化只有一个半平面L[0] 252 for(re ll i = 2; i <= n; ++i) 253 { 254 while(first < last && !OnLeft(DL[i], p[last-1])) --last; 255 while(first < last && !OnLeft(DL[i], p[first])) ++first; 256 q[++last] = DL[i]; 257 if(sgn(cross(q[last].v, q[last-1].v)) == 0) //两向量平行 258 { 259 --last; 260 if(OnLeft(q[last], DL[i].p)) q[last] = DL[i]; //同向取内测的一个 261 } 262 if(first < last) LineLineIntersection(q[last-1],q[last],p[last-1]); 263 } 264 while(first < last && !OnLeft(q[first],p[last-1])) --last; 265 //删除无用平面 266 if(last-first <= 1) return 0; //空集 267 LineLineIntersection(q[last],q[first],p[last]); //计算首尾两个半平面交的交点 268 ll cnt = 0; 269 for(re ll i = first; i <= last; ++i) poly[++cnt] = p[i]; 270 return cnt; 271 } 272 273 int main() 274 { 275 //UVALive-3890 276 ios::sync_with_stdio(false); 277 ll n; 278 //scanf("%lld", &n); 279 while(cin >> n) 280 //while(~scanf("%d", &n)) 281 { 282 if(n == 0) break; 283 POINT *p = new POINT[n+1]; 284 for(re ll i = 1; i <= n; ++i) 285 { 286 cin >> p[i].x >> p[i].y; 287 //scanf("%lf%lf%lf%lf", &p1.x, &p1.y, &p2.x, &p2.y); 288 } 289 VECTOR *v = new VECTOR[n+1]; 290 VECTOR *r = new VECTOR[n+1]; 291 for(re ll i = 1; i <= n; ++i) 292 { 293 v[i] = p[i%n+1]-p[i]; 294 r[i] = norm(v[i]); //单位正交向量 295 } 296 //二分答案 297 ld left = 0, right = 20000; 298 DIRLINE *DL = new DIRLINE[n+1]; 299 POINT *poly = new POINT[n+1]; 300 while(right-left > 1e-6) 301 { 302 ld mid = (left+right)/2; 303 for(re ll i = 1; i <= n; ++i) 304 { 305 DL[i] = DIRLINE(p[i]+r[i]*mid,v[i]); 306 } 307 ll cnt = HalfPlaneIntersection(DL,n,poly); 308 if(!cnt) right = mid; 309 else left = mid; 310 } 311 cout << fixed << setprecision(6) << left << endl; 312 } 313 return 0; 314 }