[SCOI2015]小凸想跑步
先推一波式子,设\(p(x,y)\),那么我们尝试写出点\((x_i,y_i),(x_{i+1},y_{i+1})\)和\(p(x,y)\)形成的三角形面积,就是用叉积写一波
\[2S=(x_i-x)(y_{i+1}-y)-(x_{i+1}-x)(y_i-y)
\]
大力拆开式子,发现\(2S=(y_{i}-y_{i+1})x+(x_{i+1}-x_i)y+x_iy_{i+1}-x_{i+1}y_i\)
不难发现这是一个\(Ax+By+C\)的形式,肯定跟直线有点关系,我们设\(A_i=y_{i}-y_{i+1},B_i=x_{i+1}-x_i,C_i=x_iy_{i+1}-x_{i+1}y_i\)
我们需要使得\(\forall i\in[1,n]\)有\(A_0x+B_0+C_0<A_ix+B_iy+C_i\)
即\((A_0-A_i)x+(B_0-B_i)y+(C_0-C_i)<0\)
不难发现这样的\((x,y)\)都在一个半平面内
确切的讲当\(B_0-B_i>0\)的时候,这样的\((x,y)\)都在直线\((A_0-A_i)x+(B_0-B_i)y+(C_0-C_i)=0\)的下方;当\(B_0-B_i<0\)时,就在直线的上方
于是我们对得到的这\(n\)条限制以及凸包上原来的\(n\)条边求一个半平面交就好了,答案就是半平面交的面积除以原来凸包的面积
代码
#include<bits/stdc++.h>
#define re register
const double eps=1e-6;
const int maxn=2e5+5;
int n,h,t,N,tot;
struct Pt{double x,y;}p[maxn];
struct Line{Pt s,t;double det;}L[maxn],q[maxn];
inline int dcmp(double a,double b) {return a+eps>b&&a-eps<b;}
inline Pt operator^(double t,const Pt &A) {return (Pt){t*A.x,t*A.y};}
inline double operator*(const Pt &A,const Pt &B) {return A.x*B.y-A.y*B.x;}
inline Pt operator+(const Pt &A,const Pt &B) {return (Pt){A.x+B.x,A.y+B.y};}
inline Pt operator-(const Pt &A,const Pt &B) {return (Pt){A.x-B.x,A.y-B.y};}
inline int OnRight(const Pt &c,const Line &A) {return (A.t-c)*(A.s-c)>0;}
inline Pt sec(const Line &A,const Line &B) {
Pt v1=A.t-A.s,v2=B.t-B.s,v3=A.s-B.s;
double t=(v3*v2)/(v2*v1);
return A.s+(t^v1);
}
inline int cmp(const Line &A,const Line &B) {
return dcmp(A.det,B.det)?OnRight(A.s,B):(A.det<B.det);
}
inline void Half() {
for(re int i=1;i<=N;i++)
L[i].det=atan2((L[i].t-L[i].s).y,(L[i].t-L[i].s).x);
std::sort(L+1,L+N+1,cmp);tot=1;h=1,t=0;
for(re int i=2;i<=N;i++) {
if(!dcmp(L[i].det,L[i-1].det)) ++tot;
L[tot]=L[i];
}
for(re int i=1;i<=tot;i++) {
while(h<t&&OnRight(sec(q[t],q[t-1]),L[i])) --t;
while(h<t&&OnRight(sec(q[h],q[h+1]),L[i])) ++h;
q[++t]=L[i];
}
while(h<t&&OnRight(sec(q[t],q[t-1]),q[h])) --t;
while(h<t&&OnRight(sec(q[h],q[h+1]),q[t])) ++h;
}
inline void getABC(double &A,double &B,double &C,int i) {
C=p[i]*p[i+1];A=p[i].y-p[i+1].y,B=p[i+1].x-p[i].x;
}
int main() {
scanf("%d",&n);double Sum=0,ans=0;
for(re int i=1;i<=n;i++)scanf("%lf%lf",&p[i].x,&p[i].y);p[n+1]=p[1];
for(re int i=1;i<=n;i++) Sum+=fabs((p[i]-p[1])*(p[i+1]-p[1]));
for(re int i=1;i<=n;i++) L[++N].s=p[i],L[N].t=p[i+1];
double A1,B1,C1,A,B,C;getABC(A1,B1,C1,1);Line nw;
for(re int i=2;i<=n;i++) {
getABC(A,B,C,i);C=C1-C,B=B1-B,A=A1-A;
if(!dcmp(A,0)) nw.s=(Pt){-C/A,0},nw.t=(Pt){-C/A-B,A};
else if(!dcmp(B,0)) nw.s=(Pt){0,-C/B},nw.t=(Pt){-B,-C/B+A};
L[++N]=nw;
}
Half();q[t+1]=q[h];int cnt=0;
for(re int i=h;i<=t;i++) p[++cnt]=sec(q[i],q[i+1]);
p[cnt+1]=p[1];
for(re int i=1;i<=cnt;i++) ans+=fabs((p[i]-p[1])*(p[i+1]-p[1]));
printf("%.4lf\n",ans/Sum);return 0;
}