题目链接:https://atcoder.jp/contests/abc157/tasks/abc157_f
题目大意
平面上有$N$个点,第$i$个点的坐标为$(x_i,~y_i)$,其权重为$c_i$。
现在要在平面上选一点$(X,~Y)$,设$A_i~=~c_i imes sqrt{left(X - x_i ight)^2 + left(Y-y_i ight)^2}$。
定义$T$为$A$数组中第$K$小的数。
问如何选点,能使$T$最小,输出最小值。
分析
二分$T$。
可以把$N$个点看成是$N$个圆盘,第$i$个圆盘圆心为$(x_i,~y_i)$,半径为$frac{T}{c_i}$。
然后二分的标准为:是否存在一块区域,叠了大于等于$K$个圆盘。
实际代码中,可以用点来代替区域。
暴力枚举所有可能的点即可。
如果两圆盘相交,那么交点就是可能的点。
如果两圆盘内含,那么可以选择圆心为代表点。
时间复杂度为$O(log(frac{c_{max}~~~*~~(x_{max}~~~+~~y_{max}~~)}{EPS})~*~N^3)$
代码如下
1 #include <bits/stdc++.h> 2 using namespace std; 3 4 /*-------------------Define Start-------------------*/ 5 typedef bool BL; // 布尔类型 6 typedef char SB; // 有符号1字节,8位 7 typedef unsigned char UB; // 无符号1字节,8位 8 typedef short SW; // 有符号短整型,16位 9 typedef unsigned short UW; // 无符号短整型,16位 10 typedef long SDW; // 有符号整型,32位 11 typedef unsigned long UDW; // 无符号整型,32位 12 typedef long long SLL; // 有符号长整型,64位 13 typedef unsigned long long ULL; // 无符号长整型,64位 14 typedef char CH; // 单个字符 15 typedef float R32; // 单精度浮点数 16 typedef double R64; // 双精度浮点数 17 18 #define Rep(i, n) for (register SDW i = 0; i < (n); ++i) 19 #define For(i, s, t) for (register SDW i = (s); i <= (t); ++i) 20 #define rFor(i, t, s) for (register SDW i = (t); i >= (s); --i) 21 #define foreach(i, c) for (__typeof(c.begin()) i = c.begin(); i != c.end(); ++i) 22 #define ms0(a) memset(a,0,sizeof(a)) 23 #define msI(a) memset(a,0x7f,sizeof(a)) 24 #define LOWBIT(x) ((x)&(-x)) 25 26 #define MP make_pair 27 #define PB push_back 28 #define ft first 29 #define sd second 30 #define ALL(x) x.begin(),x.end() 31 32 #define pr(x) cout << #x << " = " << x << " " 33 #define prln(x) cout << #x << " = " << x << endl 34 35 const ULL mod = 1e9 + 7; //常用模数(可根据题目需要修改) 36 const ULL inf = 0x7fffffff; //用来表示无限大 37 const ULL infLL = 0x7fffffffffffffffLL; //用来表示无限大 38 const R64 EPS = 1e-8; 39 40 // 重载<<操作符,用来打印vector 41 template < typename T > 42 ostream& operator<< (ostream& out, vector< T > vec) { 43 foreach(i, vec) { 44 if(i != vec.begin()) { 45 out << " "; 46 } 47 out << *i; 48 } 49 return out; 50 } 51 52 // 浮点数正负零判断 53 SDW sgn(R64 x) { 54 if(fabs(x) < EPS) { 55 return 0; 56 } 57 return x > 0 ? 1 : -1; 58 } 59 60 inline R64 pow2(R64 x) { 61 return x * x; 62 } 63 /*-------------------Define End-------------------*/ 64 65 template<typename T> 66 struct Point{ 67 SDW id; 68 T X,Y; 69 Point< T >(T x = 0,T y = 0) : X(x), Y(y) {} 70 inline T set(T x, T y) { return X = x, Y = y; } 71 inline Point<T> unitVec() { return Point< T >(X/abs(X), Y/abs(Y)); } 72 inline T centerDist2() { return X*X + Y*Y; } 73 inline double centerDist() { return sqrt(X*X + Y*Y); } 74 75 inline bool operator== (const Point<T>& x) const { return X == x.X && Y == x.Y; } 76 inline bool operator< (const Point<T>& x) const { 77 if(Y == x.Y)return X < x.X; 78 return Y < x.Y; 79 } 80 inline bool operator> (const Point<T>& x)const{return x < *this;} 81 82 bool between(const Point<T> &x, const Point<T> &y) { 83 if(x.X == y.X && X == x.X && Y > min(x.Y, y.Y) && Y < max(x.Y, y.Y)) return true; 84 if(x.Y == y.Y && Y == x.Y && X > min(x.X, y.X) && X < max(x.X, y.X)) return true; 85 return false; 86 } 87 88 inline Point<T> operator* (const T& k) { return Point<T>(k*X, k*Y); } 89 inline Point<T> operator/ (const T& k) { return Point<T>(X/k, Y/k); } 90 inline Point<T> operator+ (const Point<T>& x) { return Point<T>(X+x.X,Y+x.Y); } 91 inline Point<T> operator- (const Point<T>& x) { return Point<T>(X-x.X,Y-x.Y); } 92 93 T operator^(const Point<T> &x) const { return X*x.Y - Y*x.X; } 94 T operator*(const Point<T> &x) const { return X*x.X + Y*x.Y; } 95 }; 96 97 template<typename T> 98 istream &operator>> (istream &in, Point<T> &x) { 99 in >> x.X >> x.Y; 100 return in; 101 } 102 103 template<typename T> 104 ostream &operator<< (ostream &out, const Point<T> &x) { 105 out << "(" << x.X << ", " << x.Y << ")"; 106 return out; 107 } 108 109 template<typename T> 110 double dist(Point<T> a,Point<T> b) { 111 return sqrt((a-b)*(a-b)); 112 } 113 114 struct Circle{ 115 Point< R64 > p; // Circle center 116 R64 c, r; 117 }; 118 typedef Circle Disk; 119 120 istream& operator>> (istream& in, Disk &d) { 121 cin >> d.p; 122 in >> d.c; 123 return in; 124 } 125 126 /* 127 设两圆圆心分别为A,B 128 两圆半径分别为Ra,Rb 129 两圆交点分别为C(x1,y1),D(x2,y2),x1 <= x2 130 AB,CD交于点E(x0,y0) 131 AB的斜率为Kab 132 CD的斜率为Kcd 133 过点C作X轴的垂线,过点E作Y轴垂线,两直线交于点F 134 */ 135 BL getIntersectionOfTwoCircles(Circle &C1, Circle &C2, R64 &x1, R64 &y1, R64 &x2, R64 &y2) { 136 Point< R64 > A = C1.p; 137 Point< R64 > B = C2.p; 138 R64 Ra = C1.r, Rb = C2.r; 139 140 R64 AB = dist(A, B); // 圆心距 141 if(sgn(Ra + Rb - AB) < 0 || sgn(fabs(Ra - Rb) - AB) > 0) { // 两个圆没有交点的情况 142 return false; 143 } 144 145 R64 AE = (pow2(Ra) - pow2(Rb) + pow2(AB)) / (2 * AB); 146 R64 CE = sqrt(pow2(Ra) - pow2(AE)); 147 R64 x0, y0; 148 R64 Kab, Kcd; 149 R64 EF, CF; 150 151 if(sgn(A.X - B.X) == 0) { // 交点连线平行于X轴 152 x0 = A.X; 153 y0 = A.Y + (B.Y - A.Y) * (AE / AB); 154 155 x1 = x0 - CE; 156 y1 = y0; 157 x2 = x0 + CE; 158 y2 = y0; 159 } 160 else if(sgn(A.Y - B.Y) == 0) { // 交点连线平行于Y轴 161 x0 = A.X + (B.X - A.X) * (AE / AB); 162 y0 = A.Y; 163 164 x1 = x0; 165 y1 = y0 - CE; 166 x2 = x0; 167 y2 = y0 + CE; 168 } 169 else { 170 Kab = (A.Y - B.Y) / (A.X - B.X); 171 Kcd = -1 / Kab; 172 x0 = A.X + (B.X - A.X) * (AE / AB); 173 y0 = A.Y + Kab * (x0 - A.X); 174 EF = sqrt(pow2(CE) / (1 + pow2(Kcd))); 175 176 x1 = x0 - EF; 177 y1 = y0 + Kcd * (x1 - x0); 178 x2 = x0 + EF; 179 y2 = y0 + Kcd * (x2 - x0); 180 } 181 182 return true; 183 } 184 185 const UDW maxN = 6e1 + 7; 186 SDW N, K; 187 Disk disks[maxN]; 188 R64 ans; 189 190 void input(){ 191 cin >> N >> K; 192 Rep(i, N) { 193 cin >> disks[i]; 194 } 195 } 196 197 // 计数有多少disk覆盖点p 198 BL checkPoint(Point< R64 > p) { 199 R64 dis; 200 SDW cnt = 0; 201 202 Rep(i, N) { 203 dis = dist(p, disks[i].p); 204 205 if(sgn(disks[i].r - dis) >= 0) { 206 ++cnt; 207 } 208 } 209 return cnt >= K; 210 } 211 212 // 查找是否有一块区域是K个或K个以上Disk的交集 213 BL findAns(R64 T) { 214 Rep(i, N) { 215 disks[i].r = T / disks[i].c; 216 } 217 218 // 内含情况,枚举圆心 219 Rep(i, N) { 220 if(checkPoint(disks[i].p)) { 221 return true; 222 } 223 } 224 225 // 相交情况,枚举交点 226 Rep(i, N) { 227 For(j, i + 1, N - 1) { 228 R64 x1, y1, x2, y2; 229 230 if(!getIntersectionOfTwoCircles(disks[i], disks[j], x1, y1, x2, y2)) { 231 continue; 232 } 233 234 if(checkPoint(Point< R64 >(x1, y1)) || checkPoint(Point< R64 >(x2, y2))) { 235 return true; 236 } 237 } 238 } 239 240 return false; 241 } 242 243 void solve(){ 244 R64 L = 0; 245 R64 R = 2e5; 246 247 while(sgn(L - R) < 0) { 248 R64 mid = (L + R) / 2; 249 250 if(findAns(mid)) { 251 R = mid; 252 } 253 else { 254 L = mid + 1e-6; 255 } 256 } 257 258 ans = R; 259 } 260 261 void output(){ 262 printf("%f", ans); 263 } 264 265 int main() { 266 input(); 267 solve(); 268 output(); 269 return 0; 270 }