标签:
题意:求某凸多边形内部离边界最远的点到边界的距离
首先介绍半平面、半平面交的概念:
半平面:对于一条有向直线,它的方向的左手侧就是它所划定的半平面范围。如图所示:
半平面交:多个半平面的交集。有点类似二元函数的线性规划。如图
求半平面交:用的kuangbin模板= =
sol:二分答案
二分距离值,按这个值把边界向内缩,再求新的半平面交。如图:
绿色的是原图形,蓝色是按距离值向里面缩进去之后得到的新图形。对这个新图做半平面交即可。
若半平面交存在,说明与边界的距离是该值的点存在(半平面交里面的点都是)
那么如何把图形向里面缩呢?其实在纸上画画就容易多了:
因为半平面是一条有向直线,所以只要有直线上的任意一点和方向就可以了,这个点在不在对应的那条边上都无所谓。
如图,棕色是一开始的边,黑色是它对应的半平面。
设当前二分答案值是mid,对这条边做一个长度为mid的法向量(图中绿色)
棕色边起点+法向量即得到新半平面上的一个点(图中蓝色)
新半平面的方向和原来一样,直接平移过来就行了。
在做法向量的时候采用向量的概念能简化问题。这里在kuangbin模板的基础上加上了向量类Vector:
1 struct Vector:public point 2 { 3 Vector(){} 4 Vector(double a,double b) 5 { 6 x=a; y=b; 7 } 8 Vector(point _a,point _b) //a->b 9 { 10 double dx=_b.x-_a.x; 11 double dy=_b.y-_a.y; 12 x=dx; y=dy; 13 } 14 Vector(line v) 15 { 16 double dx=v.b.x-v.a.x; 17 double dy=v.b.y-v.a.y; 18 x=dx; y=dy; 19 } 20 double length() 21 { 22 return (sqrt(x*x+y*y)); 23 } 24 Vector Normal() //返回this的单位长度的法向量 25 { 26 double L=sqrt(x*x+y*y); 27 Vector Vans=Vector(-y/L,x/L); 28 return Vans; 29 } 30 }; 31
重载了三种构造函数:给出向量的坐标值(x,y)、给出向量首尾两点、用一条线段作为向量
----------------------------------------------
AC Code:
1 #include<vector> 2 #include<list> 3 #include<map> 4 #include<set> 5 #include<deque> 6 #include<queue> 7 #include<stack> 8 #include<bitset> 9 #include<algorithm> 10 #include<functional> 11 #include<numeric> 12 #include<utility> 13 #include<iostream> 14 #include<sstream> 15 #include<iomanip> 16 #include<cstdio> 17 #include<cmath> 18 #include<cstdlib> 19 #include<cctype> 20 #include<string> 21 #include<cstring> 22 #include<cstdio> 23 #include<cmath> 24 #include<cstdlib> 25 #include<ctime> 26 #include<climits> 27 #include<complex> 28 #define mp make_pair 29 #define pb push_back 30 using namespace std; 31 const double eps=1e-6; 32 const double pi=acos(-1.0); 33 const double inf=1e20; 34 const int maxp=1111; 35 int dblcmp(double d) 36 { 37 if (fabs(d)<eps)return 0; 38 return d>eps?1:-1; 39 } 40 inline double sqr(double x){return x*x;} 41 struct point 42 { 43 double x,y; 44 point(){} 45 point(double _x,double _y): 46 x(_x),y(_y){}; 47 void input() 48 { 49 scanf("%lf%lf",&x,&y); 50 } 51 void output() 52 { 53 printf("%.2f %.2f\n",x,y); 54 } 55 bool operator==(point a)const 56 { 57 return dblcmp(a.x-x)==0&&dblcmp(a.y-y)==0; 58 } 59 bool operator<(point a)const 60 { 61 return dblcmp(a.x-x)==0?dblcmp(y-a.y)<0:x<a.x; 62 } 63 double len() 64 { 65 return hypot(x,y); 66 } 67 double len2() 68 { 69 return x*x+y*y; 70 } 71 double distance(point p) 72 { 73 return hypot(x-p.x,y-p.y); 74 } 75 point add(point p) 76 { 77 return point(x+p.x,y+p.y); 78 } 79 point sub(point p) 80 { 81 return point(x-p.x,y-p.y); 82 } 83 point mul(double b) 84 { 85 return point(x*b,y*b); 86 } 87 point div(double b) 88 { 89 return point(x/b,y/b); 90 } 91 double dot(point p) 92 { 93 return x*p.x+y*p.y; 94 } 95 double det(point p) 96 { 97 return x*p.y-y*p.x; 98 } 99 double rad(point a,point b) 100 { 101 point p=*this; 102 return fabs(atan2(fabs(a.sub(p).det(b.sub(p))),a.sub(p).dot(b.sub(p)))); 103 } 104 point trunc(double r) 105 { 106 double l=len(); 107 if (!dblcmp(l))return *this; 108 r/=l; 109 return point(x*r,y*r); 110 } 111 point rotleft() 112 { 113 return point(-y,x); 114 } 115 point rotright() 116 { 117 return point(y,-x); 118 } 119 point rotate(point p,double angle)//绕点p逆时针旋转angle角度 120 { 121 point v=this->sub(p); 122 double c=cos(angle),s=sin(angle); 123 return point(p.x+v.x*c-v.y*s,p.y+v.x*s+v.y*c); 124 } 125 }; 126 struct line 127 { 128 point a,b; 129 line(){} 130 line(point _a,point _b) 131 { 132 a=_a; 133 b=_b; 134 } 135 bool operator==(line v) 136 { 137 return (a==v.a)&&(b==v.b); 138 } 139 //倾斜角angle 140 line(point p,double angle) 141 { 142 a=p; 143 if (dblcmp(angle-pi/2)==0) 144 { 145 b=a.add(point(0,1)); 146 } 147 else 148 { 149 b=a.add(point(1,tan(angle))); 150 } 151 } 152 //ax+by+c=0 153 line(double _a,double _b,double _c) 154 { 155 if (dblcmp(_a)==0) 156 { 157 a=point(0,-_c/_b); 158 b=point(1,-_c/_b); 159 } 160 else if (dblcmp(_b)==0) 161 { 162 a=point(-_c/_a,0); 163 b=point(-_c/_a,1); 164 } 165 else 166 { 167 a=point(0,-_c/_b); 168 b=point(1,(-_c-_a)/_b); 169 } 170 } 171 void input() 172 { 173 a.input(); 174 b.input(); 175 } 176 void adjust() 177 { 178 if (b<a)swap(a,b); 179 } 180 double length() 181 { 182 return a.distance(b); 183 } 184 double angle()//直线倾斜角 0<=angle<180 185 { 186 double k=atan2(b.y-a.y,b.x-a.x); 187 if (dblcmp(k)<0)k+=pi; 188 if (dblcmp(k-pi)==0)k-=pi; 189 return k; 190 } 191 //点和线段关系 192 //1 在逆时针 193 //2 在顺时针 194 //3 平行 195 int relation(point p) 196 { 197 int c=dblcmp(p.sub(a).det(b.sub(a))); 198 if (c<0)return 1; 199 if (c>0)return 2; 200 return 3; 201 } 202 bool pointonseg(point p) 203 { 204 return dblcmp(p.sub(a).det(b.sub(a)))==0&&dblcmp(p.sub(a).dot(p.sub(b)))<=0; 205 } 206 bool parallel(line v) 207 { 208 return dblcmp(b.sub(a).det(v.b.sub(v.a)))==0; 209 } 210 //2 规范相交 211 //1 非规范相交 212 //0 不相交 213 int segcrossseg(line v) 214 { 215 int d1=dblcmp(b.sub(a).det(v.a.sub(a))); 216 int d2=dblcmp(b.sub(a).det(v.b.sub(a))); 217 int d3=dblcmp(v.b.sub(v.a).det(a.sub(v.a))); 218 int d4=dblcmp(v.b.sub(v.a).det(b.sub(v.a))); 219 if ((d1^d2)==-2&&(d3^d4)==-2)return 2; 220 return (d1==0&&dblcmp(v.a.sub(a).dot(v.a.sub(b)))<=0|| 221 d2==0&&dblcmp(v.b.sub(a).dot(v.b.sub(b)))<=0|| 222 d3==0&&dblcmp(a.sub(v.a).dot(a.sub(v.b)))<=0|| 223 d4==0&&dblcmp(b.sub(v.a).dot(b.sub(v.b)))<=0); 224 } 225 int linecrossseg(line v)//*this seg v line 226 { 227 int d1=dblcmp(b.sub(a).det(v.a.sub(a))); 228 int d2=dblcmp(b.sub(a).det(v.b.sub(a))); 229 if ((d1^d2)==-2)return 2; 230 return (d1==0||d2==0); 231 } 232 //0 平行 233 //1 重合 234 //2 相交 235 int linecrossline(line v) 236 { 237 if ((*this).parallel(v)) 238 { 239 return v.relation(a)==3; 240 } 241 return 2; 242 } 243 point crosspoint(line v) 244 { 245 double a1=v.b.sub(v.a).det(a.sub(v.a)); 246 double a2=v.b.sub(v.a).det(b.sub(v.a)); 247 return point((a.x*a2-b.x*a1)/(a2-a1),(a.y*a2-b.y*a1)/(a2-a1)); 248 } 249 double dispointtoline(point p) 250 { 251 return fabs(p.sub(a).det(b.sub(a)))/length(); 252 } 253 double dispointtoseg(point p) 254 { 255 if (dblcmp(p.sub(b).dot(a.sub(b)))<0||dblcmp(p.sub(a).dot(b.sub(a)))<0) 256 { 257 return min(p.distance(a),p.distance(b)); 258 } 259 return dispointtoline(p); 260 } 261 point lineprog(point p) 262 { 263 return a.add(b.sub(a).mul(b.sub(a).dot(p.sub(a))/b.sub(a).len2())); 264 } 265 point symmetrypoint(point p) 266 { 267 point q=lineprog(p); 268 return point(2*q.x-p.x,2*q.y-p.y); 269 } 270 }; 271 272 struct Vector:public point 273 { 274 Vector(){} 275 Vector(double a,double b) 276 { 277 x=a; y=b; 278 } 279 Vector(point _a,point _b) //a->b 280 { 281 double dx=_b.x-_a.x; 282 double dy=_b.y-_a.y; 283 x=dx; y=dy; 284 } 285 Vector(line v) 286 { 287 double dx=v.b.x-v.a.x; 288 double dy=v.b.y-v.a.y; 289 x=dx; y=dy; 290 } 291 double length() 292 { 293 return (sqrt(x*x+y*y)); 294 } 295 Vector Normal() 296 { 297 double L=sqrt(x*x+y*y); 298 Vector Vans=Vector(-y/L,x/L); 299 return Vans; 300 } 301 }; 302 303 struct halfplane:public line //半平面 304 { 305 double angle; 306 halfplane(){} 307 //表示向量 a->b逆时针(左侧)的半平面 308 halfplane(point _a,point _b) 309 { 310 a=_a; 311 b=_b; 312 } 313 halfplane(line v) 314 { 315 a=v.a; 316 b=v.b; 317 } 318 void calcangle() 319 { 320 angle=atan2(b.y-a.y,b.x-a.x); 321 } 322 bool operator<(const halfplane &b)const 323 { 324 return angle<b.angle; 325 } 326 }; 327 struct halfplanes //半平面交 328 { 329 int n; 330 halfplane hp[maxp]; 331 point p[maxp]; 332 int que[maxp]; 333 int st,ed; 334 void push(halfplane tmp) 335 { 336 hp[n++]=tmp; 337 } 338 void unique() 339 { 340 int m=1,i; 341 for (i=1;i<n;i++) 342 { 343 if (dblcmp(hp[i].angle-hp[i-1].angle))hp[m++]=hp[i]; 344 else if (dblcmp(hp[m-1].b.sub(hp[m-1].a).det(hp[i].a.sub(hp[m-1].a))>0))hp[m-1]=hp[i]; 345 } 346 n=m; 347 } 348 bool halfplaneinsert() 349 { 350 int i; 351 for (i=0;i<n;i++)hp[i].calcangle(); 352 sort(hp,hp+n); 353 unique(); 354 que[st=0]=0; 355 que[ed=1]=1; 356 p[1]=hp[0].crosspoint(hp[1]); 357 for (i=2;i<n;i++) 358 { 359 while (st<ed&&dblcmp((hp[i].b.sub(hp[i].a).det(p[ed].sub(hp[i].a))))<0)ed--; 360 while (st<ed&&dblcmp((hp[i].b.sub(hp[i].a).det(p[st+1].sub(hp[i].a))))<0)st++; 361 que[++ed]=i; 362 if (hp[i].parallel(hp[que[ed-1]]))return false; 363 p[ed]=hp[i].crosspoint(hp[que[ed-1]]); 364 } 365 while (st<ed&&dblcmp(hp[que[st]].b.sub(hp[que[st]].a).det(p[ed].sub(hp[que[st]].a)))<0)ed--; 366 while (st<ed&&dblcmp(hp[que[ed]].b.sub(hp[que[ed]].a).det(p[st+1].sub(hp[que[ed]].a)))<0)st++; 367 if (st+1>=ed)return false; 368 return true; 369 } 370 /* 371 void getconvex(polygon &con) 372 { 373 p[st]=hp[que[st]].crosspoint(hp[que[ed]]); 374 con.n=ed-st+1; 375 int j=st,i=0; 376 for (;j<=ed;i++,j++) 377 { 378 con.p[i]=p[j]; 379 } 380 }*/ 381 }; 382 383 point p[1000]; 384 Vector V[1000],Vt[1000]; 385 halfplanes TH; 386 int n; 387 388 int main() 389 { 390 //freopen("in.txt","r",stdin); 391 392 while (cin>>n) 393 { 394 if (n==0) break; 395 for (int i=0;i<n;i++) //n points:[0..n-1] 396 p[i].input(); 397 398 for (int i=0;i<n;i++) //v[i]:p[i]->p[i+1] 399 { 400 V[i]=Vector(p[i],p[(i+1)%n]); 401 //printf("vector: %.6lf %.6lf\n",V[i].x,V[i].y); 402 Vt[i]=V[i].Normal(); 403 } 404 405 double l=0,r=5001; 406 while (r-l>eps) 407 { 408 double mid=(l+r)/2; 409 //double mid=l+(r-l)/2; 410 TH.n=0; 411 //halfplanes TH; 412 413 for (int i=0;i<n;i++) 414 { 415 point t1=p[i].add(Vt[i].mul(mid)); 416 point t2=t1.add(V[i]); 417 418 line tmp=line(t1,t2); 419 TH.push(halfplane(tmp)); 420 } 421 422 //printf("%.6lf %d",mid,TH.n); 423 424 if (TH.halfplaneinsert()) 425 { 426 //cout<<"OK"<<endl; 427 l=mid; //l=mid+0.00001; 428 } 429 else 430 { 431 //cout<<"NO"<<endl; 432 r=mid; //r=mid-0.00001; 433 } 434 } 435 printf("%.6lf\n",l); 436 } 437 return 0; 438 }
推荐iPad上一个免费又好用的绘图板:Sketches
标签:
原文地址:http://www.cnblogs.com/pdev/p/4275911.html