Problem C
Fighting against a Polygonal Monster
Input: Standard Input
Output: Standard Output
You're fighting against a monster. You have a special weapon that can shoot a laser beam. The laser beam can be seen as a cylinder along the z-axis, So when look from above the XY plane, the laser beam is actually a circle centered at (0,0) (yes, you cannot change the position of center). The monster is a convex polygon (well, you may think of it as a very thin prism) with n vertices surrounding the origin (i.e. (0,0) is strictly inside the monster, not on his boundary).
You know that, when the common area of the laser beam and the monster is at least R, the monster dies. Since larger laser beam consumes more power, you're interested in the minimum radius of the laser beam.
Write a program to find the minimum radius. It is guaranteed that (0 £ R £A), whereA is the area of the monster.
Input
The input consists at most 10 cases. Each case starts with a single integer n (3 £ n £ 50) and a floating-point number R followed by n lines of two real numbers: the coordinates of the monster. The points are arranged in counter-clockwise order or clockwise order. The last test case is followed by a single zero, which should not be processed. The meaning of N and R are given in the problem statement.
Output
For each test case, print the case number and the minimum radius, to 2 decimal places. Inputs will be such that small precision errors will not change the visible output if you use double-precision floating-point numbers.
Sample Input Output for Sample Input
| 3 1.60 -1 -1 1 -1 0 1 0
| Case 1: 0.93 |
Problem setter: Rujia Liu
Alternative Solutions: Wenbin Tang, Yiming Li
Originally used in the Sichuan Team Selection Contest for NOI, 2006.
思路：二分半径我想大部分人都能想到。这个题目面积讨论起来其实很麻烦的。。。我的讨论方法是拿出凸多边形的一条边，看看这条边的两个端点在圆的什么地方，如果都在圆外，显然面积是扇形，如果两个端点都在圆内，显然是三角形，如果一个在外面，一个在里面，显然是一个扇形+三角形，如果一个在圆上，一个在圆外，是扇形，如果一个在圆上，一个在圆内，是三角形，如果两个都在外面，有两种情况，交点不在线段上面的话是扇形，如果都在线段上面是两个扇形+一个三角形。不可能一个点在线段上，一个不在的。。。但是不知道为什么我的代码会出现这种情况。。。。不过我改了个写法就过了，改成用OnSegment不会有问题了。慢慢写，把情况分清楚，你能做出来的。
代码：
#include<iostream>
#include<cstdio>
#include<string.h>
#include<math.h>
#include<string>
#include<cstring>
#include<map>
#include<algorithm>
#include<vector>
using namespace std;
#define eps 1e-10
const double inf = 1e16;
const double PI = acos(-1.0);
struct Point
{
Point(const Point&p) { x = p.x , y = p.y; }
Point (double xx=0,double yy=0) : x(xx) , y(yy) { }
double x;
double y;
};
typedef Point Vector;
Vector operator+(Vector v1,Vector v2) { return Vector(v1.x+v2.x,v1.y+v2.y); }
Vector operator-(Vector v1,Vector v2) { return Vector(v1.x-v2.x,v1.y-v2.y); }
Vector operator*(Vector v, double p) { return Vector(v.x*p,v.y*p); }
Vector operator/(Vector v,double p) { return Vector(v.x/p,v.y/p); }
bool operator < (Point a,Point b) { return a.x < b.x || (a.x==b.x && a.y > b.y); }
int dcmp(double x)
{
if (fabs(x) < eps) return 0;
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;
}
inline double toRad(double x) { return x * PI/180; }
inline double toDegreed(double rad) { return rad*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 fabs(Cross(b-a,c-a)); }
//旋转
Vector Rotate(Vector A,double rad)
{
return Vector(A.x*cos(rad)-A.y*sin(rad),A.x*sin(rad)+A.y*cos(rad));
}
//单位法线
Vector Normal(Vector A) { double L = Length(A); return Vector(-A.y/L,A.x/L); }
//点和直线
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;
}
double DistanceToLine(Point P,Point A,Point B)
{
Vector v1 = B-A , v2 = P-A;
return fabs(Cross(v1,v2))/Length(v1);
}
double DistanceToSegment(Point P,Point A,Point B)
{
if (A==B) return Length(P-A);
Vector v1 = B-A , v2 = P-A , v3 = P-B;
if (dcmp(Dot(v1,v2)) < 0) return Length(v2);
else if (dcmp(Dot(v1,v3)) > 0) return Length(v3);
else return fabs(Cross(v1,v2))/Length(v1);
}
//点在直线上的投影
Point GetLineProjection(Point P,Point A,Point B)
{
Vector v = B-A;
return A+v*(Dot(v,P-A)/Dot(v,v));
}
//线段相交(不包括端点）
bool SegmentProperIntersection(Point a1,Point a2,Point b1,Point b2)
{
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);
return dcmp(c1)*dcmp(c2) < 0 && dcmp(c3)*dcmp(c4)<0;
}
//点在线段上（不包括端点）
bool OnSegment(Point p,Point a,Point b)
{
return dcmp(Cross(a-p,b-p))==0 && dcmp(Dot(a-p,b-p)) < 0;
}
//--------------------------------------------------------------------------------------------
//直线和直线
struct Line
{
Point P; //直线上任意一点
Vector v; // 方向向量。它的左边就是对应的半平面
double ang; //极角
Line() { }
Line(Point P,Vector v) { this->P = P ; this->v = v; ang = atan2(v.y,v.x); }
bool operator < (const Line& L) const { return ang < L.ang; } //排序用的比较运算符
Point point(double t) { return v*t+P; }
};
//点p在有向直线L的左边（线上不算）
bool OnLeft(Line L , Point p) { return Cross(L.v,p-L.P) > 0; }
//二直线交点。假定交点唯一存在
Point GetIntersection(Line a,Line b)
{
Vector u = a.P-b.P;
double t = Cross(b.v,u) / Cross(a.v,b.v);
return a.P+a.v*t;
}
//半平面交的过程
int HalfplaneIntersection(Line* L,int n,Point* poly)
{
sort(L,L+n); //按极角排序
int first , last; //双端队列的第一个元素和最后一个元素的下表
Point *p = new Point[n]; //p[i]为q[i]和q[i+1]的交点
Line *q = new Line[n]; //双端队列
q[first=last=0] = L[0]; //双端队列初始化为只有一个半平面L[0]
for (int i = 1 ; i < n ; ++i) {
while (first < last && !OnLeft(L[i],p[last-1])) --last;
while (first < last && !OnLeft(L[i],p[first])) ++first;
q[++last] = L[i];
if (fabs(Cross(q[last].v,q[last-1].v)) < eps) {
--last;
if (OnLeft(q[last],L[i].P)) q[last] = L[i];
}
if (first < last) p[last-1] = GetIntersection(q[last-1],q[last]);
}
while (first < last && !OnLeft(q[first],p[last-1])) --last;
//删除无用平面（*）
if (last - first <=1 ) return 0; //空集(**)
p[last] = GetIntersection(q[last],q[first]); //计算首尾两个半平面的交点
//从deque复制到输出中
int m = 0;
for (int i = first ; i <= last ; ++i) poly[m++] = p[i];
return m ;
}
//--------------------------------------------
//与圆相关
struct Circle
{
Circle() { }
Point c;
double r;
Circle(Point c, double r) : c(c) , r(r) { }
Point point (double a) { return Point(c.x+cos(a)*r,c.y+sin(a)*r); }
};
int getLineCircleIntersection(Line L,Circle C,double &t1,double &t2,vector<Point>& sol)
{
double a = L.v.x , b = L.P.x-C.c.x , c= L.v.y, d = L.P.y-C.c.y;
double e = a*a+c*c , f = 2*(a*b+c*d) , g = b*b+d*d-C.r*C.r;
double delta = f*f-4*e*g; //判别式
if (dcmp(delta) < 0) return 0; //相离
if (dcmp(delta)==0) { //相切
t1 = t2 = -f/(2*e);
sol.push_back(L.point(t1));
return 1;
}
//相交
t1 = (-f-sqrt(delta)) / (2*e); sol.push_back(L.point(t1));
t2 = (-f+sqrt(delta)) / (2*e); sol.push_back(L.point(t2));
return 2;
}
double angle(Vector v) { return atan2(v.y,v.x); }
int getCircleCircleIntersection(Circle C1,Circle C2,vector<Point>& 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)); //向量C1C2的极角
//C1C2到C1P1的角
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;
}
//国电p到圆C的切线。v[i]是第i条切线的向量。返回切线条数
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;
}
}
int getTangents(Circle A,Circle B,Point* a, Point* b)
{
int cnt = 0;
if (A.r < B.r) { swap(A,B); swap(a,b); }
double d2 = Dot(A.c-B.c,A.c-B.c);
double rdiff = A.r-B.r;
double rsum = A.r+B.r;
if (dcmp(d2-rdiff*rdiff) < 0) return 0; //内含
double base = atan2(B.c.y-A.c.y,B.c.x-A.c.x);
if (d2==0 && dcmp(A.r-B.r)==0) return -1; //无限多条切线
if (dcmp(d2-rdiff*rdiff)==0) { //内切，1条切线
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(PI+base+ang); ++cnt;
a[cnt] = A.point(base-ang); b[cnt] = B.point(PI+base-ang); ++cnt;
}
return cnt;
}
//点p和圆的关系: 0:在圆上 1：在圆外 -1：在圆内
int PointCircleRelation(Point p,Circle c)
{
return dcmp(Dot(p-c.c,p-c.c)-c.r*c.r);
}
//A在B内
bool InCircle(Circle A,Circle B)
{
if (dcmp(A.r-B.r)>0) return false;
double d2 = Dot(A.c-B.c,A.c-B.c);
double rdiff = A.r-B.r;
double rsum = A.r+B.r;
if (dcmp(d2-rdiff*rdiff) <= 0) return true; //内含或内切或重合
return false;
}
//----------------------------------------------------------------------------
//与球相关的转化
//角度转换为弧度
double torad(double deg) { return deg/180*PI; }
//经纬度（角度）转化为空间坐标
void get_coord(double R,double lat,double lng,double& x,double& y,double& z)
{
lat = torad(lat);
lng = torad(lng);
x = R*cos(lat)*cos(lng);
y = R*cos(lat)*sin(lng);
z = R*sin(lat);
}
//-----------------------------------------------------------------------
const int maxn = 50+5;
int n;
double S;
Point p[maxn];
const Point center = Point(0,0);
bool Ok(double R)
{
Circle c(center,R);
double T = 0 , A = 0;
for (int i = 0 ; i < n ; ++i) {
Vector v1 = p[i] , v2 = p[(i+1)%n];
int k1 = PointCircleRelation(v1,c) , k2 = PointCircleRelation(v2,c);
if (k1==0 && k2==0) { T += Area2(center,v1,v2); continue; }
if (k1+k2==1) { A += Angle(v1,v2); continue; }
if (k1+k2<=-1) { T += Area2(center,v1,v2); continue; }
vector<Point> inter;
double t1 , t2;
getLineCircleIntersection(Line(v1,v2-v1),c,t1,t2,inter);
if (k1==-1 && k2==1) {
if (OnSegment(inter[0],v1,v2)) {
A += Angle(v2,inter[0]);
T += Area2(center,v1,inter[0]);
} else if (OnSegment(inter[1],v1,v2)) {
A += Angle(v2,inter[1]);
T += Area2(center,v1,inter[1]);
};
continue;
}
if (k1==1 && k2==-1) {
if (OnSegment(inter[0],v1,v2)) {
A += Angle(v1,inter[0]);
T += Area2(center,v2,inter[0]);
} else if (OnSegment(inter[1],v1,v2)) {
A += Angle(v1,inter[1]);
T += Area2(center,v2,inter[1]);
}
continue;
}
if (inter.size() < 2) A += Angle(v1,v2);
else if (inter.size()==2) {
if (dcmp(t1)>0 && dcmp(t1-1)<0) {
if (dcmp(t2)>0 && dcmp(t2-1)<0) {
A += Angle(v1,v2)-Angle(inter[0],inter[1]);
T += Area2(center,inter[0],inter[1]);
} else if (k1==1 && k2==1) while (true);
}
else A += Angle(v1,v2);
}
}
return dcmp((T+A*R*R)/2-S) >= 0;
}
int main()
{
int k = 0;
while (scanf("%d",&n)==1,n) {
++k;
scanf("%lf",&S);
double left = 0 , right = 0;
for (int i = 0 ; i < n ; ++i)
{
scanf("%lf%lf",&p[i].x,&p[i].y);
right = max(right,Length(p[i]));
}
double mid = left+(right-left)*0.5;
double ans = inf;
while (dcmp(left-right)<0) {
if (Ok(mid)) {
ans = min(ans,mid);
right = mid;
} else
left = mid;
mid = left+(right-left)*0.5;
}
printf("Case %d: %.2lf\n",k,ans+eps);
}
}
609
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
