BZOJ4445 SCOI2015小凸想跑步(半平面交)
考虑怎样的点满足条件。设其为(xp,yp),则要满足(x0-xp,y0-yp)×(x1-xp,y1-yp)<=(xi-xp,yi-yp)×(xi+1-xp,yi+1-yp)对任意i成立。拆开式子,有(x0-xp)*(y1-yp)-(y0-yp)*(x1-xp)<=(xi-xp)*(yi+1-yp)-(yi-yp)*(xi+1-xp),也即x0y1-x0yp-xpy1-y0x1+y0xp+ypx1<=xiyi+1-xiyp-xpyi+1-yixi+1+yixp+ypxi+1。移项,得(y0-y1+yi+1-yi)xp+(x1-x0+xi-xi+1)yp<=xiyi+1-yixi+1-x0y1+y0x1。一长串乱七八糟的限制都是对(xp,yp)这个二元组的,半平面交即可。当然再套上一个凸包自身的限制。
#include<iostream> #include<cstdio> #include<cmath> #include<cstdlib> #include<cstring> #include<algorithm> using namespace std; #define ll long long #define N 100010 #define vector point #define double long double char getc(){char c=getchar();while ((c<'A'||c>'Z')&&(c<'a'||c>'z')&&(c<'0'||c>'9')) c=getchar();return c;} int gcd(int n,int m){return m==0?n:gcd(m,n%m);} int read() { int x=0,f=1;char c=getchar(); while (c<'0'||c>'9') {if (c=='-') f=-1;c=getchar();} while (c>='0'&&c<='9') x=(x<<1)+(x<<3)+(c^48),c=getchar(); return x*f; } const double eps=1E-8; struct point { double x,y; vector operator +(const vector&a) const { return (vector){x+a.x,y+a.y}; } vector operator -(const vector&a) const { return (vector){x-a.x,y-a.y}; } double operator *(const vector&a) const { return x*a.y-y*a.x; } vector operator *(const double&a) const { return (vector){x*a,y*a}; } }a[N]; struct line { point a;vector p; bool operator <(const line&a) const { return atan2(p.x,p.y)>atan2(a.p.x,a.p.y); } }b[N<<1]; int n,m,head,tail; point p[N<<1]; line q[N<<1]; double area; bool onright(point a,line p) { return (a-p.a)*p.p>=0; } point cross(line a,line b) { return a.a+a.p*((b.p*(b.a-a.a))/(b.p*a.p)); } int main() { #ifndef ONLINE_JUDGE freopen("bzoj4445.in","r",stdin); freopen("bzoj4445.out","w",stdout); const char LL[]="%I64d\n"; #else const char LL[]="%lld\n"; #endif n=read(); for (int i=0;i<n;i++) a[i].x=read(),a[i].y=read(); for (int i=0;i<n;i++) area+=a[i]*a[(i+1)%n]; for (int i=0;i<n;i++) m++,b[m].a=a[i],b[m].p=a[(i+1)%n]-a[i]; for (int i=1;i<n;i++) { double A=(a[0].y-a[1].y+a[(i+1)%n].y-a[i].y),B=(a[1].x-a[0].x+a[i].x-a[(i+1)%n].x),C=(a[i]*a[(i+1)%n]-a[0]*a[1]); m++;b[m].p=(vector){-B,A}; if (fabs(B)>0.5) b[m].a=(point){0,C/B}; else b[m].a=(point){C/A,0}; } sort(b+1,b+m+1); head=tail=1;q[1]=b[1]; for (int i=2;i<=m;i++) { while (head<tail&&onright(p[tail],b[i])) tail--; while (head<tail&&onright(p[head+1],b[i])) head++; q[++tail]=b[i]; if (fabs(q[tail].p*q[tail-1].p)<eps) { tail--; if (onright(q[tail].a,b[i])) q[tail]=b[i]; } if (head<tail) p[tail]=cross(q[tail],q[tail-1]); } while (head<tail&&onright(p[tail],q[head])) tail--; p[head]=cross(q[head],q[tail]); double area2=0; for (int i=head;i<tail;i++) area2+=p[i]*p[i+1];area2+=p[tail]*p[head]; #undef double printf("%.4f",(double)(area2/area)); return 0; }