题目连接:http://acm.csu.edu.cn:20080/csuoj/problemset/problem?pid=2300
题目大意:给定一个凸多边形和一个圆,移动圆的位置使得和多边形的交最大,求交的最大面积。
此题显然是个凸函数,所以可以三分套三分,也可模拟退火,感觉还是三分套三分稳一点。
关于如何求多边形和圆的交,这里有一篇博客:http://www.cnblogs.com/lxglbk/archive/2012/08/12/2634192.html
代码里也写了点注释,可以参考一下。
AC代码:
#include<bits/stdc++.h>
using namespace std;
const int N=15;
const double pi=acos(-1.0);
const double eps=1e-9;
const double inf=1e9;
struct Point{
double x,y;
Point(){};
Point(double x,double y):x(x),y(y){}
Point operator-(const Point p)const{
return Point(x-p.x,y-p.y);
}
Point operator+(const Point p)const{
return Point(x+p.x,y+p.y);
}
double operator*(const Point p)const{
return x*p.x+y*p.y;
}
double operator^(const Point p)const{
return x*p.y-y*p.x;
}
}a[N],p[N];
inline double dis(Point a,Point b){
return sqrt((b-a)*(b-a));
}
inline double cross(Point a,Point b,Point c){//叉积,关于a
return (b-a)^(c-a);
}
inline double dot(Point a,Point b,Point c){//点积,关于a
return (b-a)*(c-a);
}
inline bool online(Point a,Point b,Point c){//c是否在线段ab上
return cross(a,b,c)==0&&c.x>=min(a.x,b.x)&&\
c.y>=min(a.y,b.y)&&c.x<=max(a.x,b.x)&&c.y<=max(a.y,b.y);
}
double area(Point b,Point c,double r){//圆和一个顶点为圆心的三角形的交
Point a(0.0,0.0);
if(dis(b,c)<eps)
return 0.0;
double h=fabs(cross(a,b,c))/dis(b,c);
if(dis(a,b)>r-eps&&dis(a,c)>r-eps){//两个端点都在圆的外面
double angle=acos(dot(a,b,c)/dis(a,b)/dis(a,c));
if(h>r-eps){
return 0.5*r*r*angle;
}
else if(dot(b,a,c)>0&&dot(c,a,b)){
double angle1=2*acos(h/r);
return 0.5*r*r*fabs(angle-angle1)+0.5*r*r*sin(angle1);
}
else{
return 0.5*r*r*angle;
}
}
else if(dis(a,b)<r+eps&&dis(a,c)<r+eps){ //两个端点都在圆内的情况
return 0.5*fabs(cross(a,b,c));
}
else{//一个端点在圆外,一个端点在圆内的情况
if(dis(a,b)>dis(a,b)){//默认b在圆内
swap(b,c);
}
if(fabs(dis(a,b)<eps)){//ab距离为零直接返回
return 0.0;
}
if(dot(b,a,c)<eps){//角abc为钝角
double angle1=acos(h/dis(a,b));
double angle2=acos(h/r)-angle1;
double angle3=acos(h/dis(a,c))-angle2;
return 0.5*r*dis(a,b)*sin(angle1+angle2)+0.5*r*r*angle3;
}
else{ //角abc为锐角
double angle1=acos(h/dis(a,b));
double angle2=acos(h/r);
double angle3=acos(h/dis(a,c))-angle2;
return 0.5*r*dis(a,b)*sin(angle1+angle2)+0.5*r*r*angle3;
}
}
}
int n;
double r;
double gets(Point o){//圆和多边形的交
for(int i=0;i<n;i++)
p[i]=a[i]-o;
o=Point(0.0,0.0);
double sum=0;
for(int i=0;i<n;i++){
int nxt=(i+1)%n;
double s=area(p[i],p[nxt],r);
if(cross(o,p[i],p[nxt])>0)
sum+=s;
else
sum-=s;
}
return fabs(sum);
}
double cal(double x){
double ly=inf,ry=-inf;
for(int i=0,j=1;i<n;i++,j=(j+1)%n){
if((a[i].x-x)*(a[j].x-x)<=0){//仅当边和X=x相交时才要更新ly,ry
double tmp=a[i].y+(a[j].y-a[i].y)/(a[j].x-a[i].x)*(x-a[i].x);
ly=min(ly,tmp);
ry=max(ry,tmp);
}
}
double s1,s2;
while(ry-ly>eps){
double m1=(ly+ly+ry)/3;
double m2=(ly+ry+ry)/3;
s1=gets(Point(x,m1)),s2=gets(Point(x,m2));
if(s1>s2)
ry=m2;
else
ly=m1;
}
return s1;
}
int main(){
double lx,rx;
scanf("%d%lf",&n,&r);
lx=inf,rx=-inf;
for(int i=0;i<n;i++){
scanf("%lf%lf",&a[i].x,&a[i].y);
lx=min(lx,a[i].x);
rx=max(rx,a[i].x);
}
double s1,s2;
while(rx-lx>eps){
double m1=(lx+lx+rx)/3;
double m2=(lx+rx+rx)/3;
s1=cal(m1),s2=cal(m2);
if(s1>s2)
rx=m2;
else
lx=m1;
}
printf("%.6f\n",s1);
return 0;
}