计算几何
谢Achen提醒,差点忘了有计算几何这东西。
poj2187 平面最远点对,凸包+旋转卡壳模板
//Serene #include<algorithm> #include<iostream> #include<cstring> #include<cstdlib> #include<cstdio> #include<cmath> using namespace std; #define ll long long #define For(i,a,b) for(int i=(a);i<=(b);++i) #define Rep(i,a,b) for(int i=(a);i>=(b);--i) const int maxn=5e4+7; int n; char cc; ll ff; template<typename T>void read(T& aa) { aa=0;cc=getchar();ff=1; while((cc<'0'||cc>'9')&&cc!='-') cc=getchar(); if(cc=='-') ff=-1,cc=getchar(); while(cc>='0'&&cc<='9') aa=aa*10+cc-'0',cc=getchar(); aa*=ff; } const double eps=1e-10; double pf(double x){return x*x;} inline int dcmp(const long double &x){return abs(x)<eps? 0:(x>0? 1:-1);} struct Xl{ double x,y; Xl(double x=0.,double y=0.):x(x),y(y){} bool operator == (const Xl &b) const{return (!dcmp(x-b.x))&&(!dcmp(y-b.y));} bool operator < (const Xl &b) const{return !dcmp(x-b.x)? y<b.y:x<b.x;} Xl operator + (const Xl &b) const{return Xl(x+b.x,y+b.y);} Xl operator - (const Xl &b) const{return Xl(x-b.x,y-b.y);} Xl operator * (const double &b) const{return Xl(x*b,y*b);} Xl operator / (const double &b) const{return Xl(x/b,y/b);} double len() const{return pf(x)+pf(y);} }node[maxn],zz[maxn]; int t; double D_(const Xl& a,const Xl& b){return a.x*b.x+a.y*b.y;} double X_(const Xl& a,const Xl& b){return a.x*b.y-a.y*b.x;} bool cmp(const Xl& a,const Xl& b){ return dcmp(X_(a-node[1],b-node[1]))==0? (a-node[1]).len()<(b-node[1]).len() : dcmp(X_(a-node[1],b-node[1]))>0; } void graham() { For(i,2,n) if(node[i]<node[1]) swap(node[i],node[1]); sort(node+2,node+n+1,cmp); zz[++t]=node[1]; zz[++t]=node[2]; For(i,3,n) { while(t>1&&dcmp(X_(zz[t]-zz[t-1],node[i]-zz[t]))<=0) --t;// "<" -> WA zz[++t]=node[i]; } } double RC() { zz[t+1]=zz[1]; int now=2;double rs=0; For(i,1,t) { while(X_(zz[i+1]-zz[i],zz[now]-zz[i])<X_(zz[i+1]-zz[i],zz[now+1]-zz[i])) if((++now)>t) now=1; rs=max(rs,(zz[now]-zz[i]).len()); } return rs; } int main() { read(n); For(i,1,n) { read(node[i].x); read(node[i].y); } graham(); printf("%lld",(ll)RC()); return 0; } /* 4 0 0 0 1 1 3 1 0 */
hdu1007 平面最近点对模板
分治,时间复杂度是nlog(n)的,由鸽巢原理可知,在那个for套for的循环里面枚举实际上最多只会枚举6个点(不会证)
//Serene #include<algorithm> #include<iostream> #include<cstring> #include<cstdlib> #include<cstdio> #include<cmath> using namespace std; #define ll long long #define db double #define For(i,a,b) for(int i=(a);i<=(b);++i) #define Rep(i,a,b) for(int i=(a);i>=(b);--i) const int maxn=1e5+7; const db INF=1e18; int n; db ans; char cc; ll ff; template<typename T>void read(T& aa) { aa=0;cc=getchar();ff=1; while((cc<'0'||cc>'9')&&cc!='-') cc=getchar(); if(cc=='-') ff=-1,cc=getchar(); while(cc>='0'&&cc<='9') aa=aa*10+cc-'0',cc=getchar(); aa*=ff; } const db eps=1e-10; db pf(db x){return x*x;} struct Xl{ db x,y; Xl(db x=0.,db y=0.):x(x),y(y){} Xl operator + (const Xl& b)const{return Xl(x+b.x,y+b.y);} Xl operator - (const Xl& b)const{return Xl(x-b.x,y-b.y);} db len()const{return sqrt(pf(x)+pf(y));} }node[maxn],zz[maxn]; bool cmp1(const Xl& a,const Xl& b) { return a.x==b.x? a.y<b.y : a.x<b.x; } bool cmp2(const Xl& a,const Xl& b) { return a.y<b.y; } void solve(int l,int r) { if(l>=r) return; if(l+1==r) { ans=min(ans,(node[l]-node[r]).len()); return; } int mid=(l+r)>>1,t=0; solve(l,mid); solve(mid+1,r); For(i,l,r) if(fabs(node[i].x-node[mid].x)<ans) zz[++t]=node[i]; sort(zz+1,zz+t+1,cmp2); For(i,1,t) for(int j=i+1;j<=t&&zz[j].y-zz[i].y<ans;++j) ans=min(ans,(zz[j]-zz[i]).len()); } int main() { read(n);db x,y; while(n) { For(i,1,n) { scanf("%lf%lf",&x,&y); node[i]=Xl(x,y); } sort(node+1,node+n+1,cmp1); ans=INF; solve(1,n); printf("%.2f\n",ans/2); read(n); } return 0; }
poj3335 判断多边形内核存在 ,半平面交模板
注意onleft中要写>=而不是>这样一个点的情况才能把边保留下来
最后不能也判面积否则就把一个点的情况判掉了。
为此debug了一下午。
//Serene #include<algorithm> #include<iostream> #include<cstring> #include<cstdlib> #include<cstdio> #include<cmath> using namespace std; #define ll long long #define For(i,a,b) for(int i=(a);i<=(b);++i) const int maxn=300+7; int T,n; char cc; ll ff; template<typename T>void read(T& aa) { aa=0;cc=getchar();ff=1; while((cc<'0'||cc>'9')&&cc!='-') cc=getchar(); if(cc=='-') ff=-1,cc=getchar(); while(cc>='0'&&cc<='9') aa=aa*10+cc-'0',cc=getchar(); aa*=ff; } const double eps=1e-8; inline int dcmp(const long double &x){return fabs(x)<eps? 0:(x>0? 1:-1);} struct Xl{ double x,y; Xl(double x=0.,double y=0.):x(x),y(y){} Xl operator + (const Xl &b) const{return Xl(x+b.x,y+b.y);} Xl operator - (const Xl &b) const{return Xl(x-b.x,y-b.y);} Xl operator * (const double &b) const{return Xl(x*b,y*b);} double angle() const{return atan2(y,x);} }node[maxn],p[maxn]; double D_(const Xl& a,const Xl& b){return a.x*b.x+a.y*b.y;} double X_(const Xl& a,const Xl& b){return a.x*b.y-a.y*b.x;} struct Li{ Xl u,v; double rad; Li(Xl u=Xl(0.,0.),Xl v=Xl(0.,0.)):u(u),v(v){rad=(v-u).angle();} bool operator < (const Li &b) const{return rad<b.rad;} }li[maxn],zz[maxn]; Xl lijd(const Li &a,const Li& b) { Xl x=a.v-a.u,y=b.v-b.u,z=a.u-b.u; long double c=X_(y,z); long double t=c/(c+X_(x+z,y)); return a.u+x*t; } bool onleft(const Xl &a,const Li &l) {return dcmp(X_(l.v-l.u,a-l.u))>=0;}//">" -> WA bool work() { For(i,1,n-1) li[i]=Li(node[i],node[i+1]); li[n]=Li(node[n],node[1]); sort(li+1,li+n+1); int s=1,t=0; zz[++t]=li[1]; For(i,2,n) { while(s<t&&!onleft(p[t-1],li[i])) t--; while(s<t&&!onleft(p[s],li[i])) s++; zz[++t]=li[i]; if(!dcmp(X_(li[i].v-li[i].u,zz[t-1].v-zz[t-1].u))) if(onleft(li[i].u,zz[--t])) zz[t]=li[i]; if(s<t) p[t-1]=lijd(zz[t],zz[t-1]); } while(s<t&&!onleft(p[t-1],li[s])) --t; return t-s>1;// } int main() { read(T); ll x,y; while(T--) { read(n); For(i,1,n) { read(x); read(y); node[n-i+1]=Xl(x,y); } if(work()) printf("YES\n"); else printf("NO\n"); } return 0; } /* 1 11 29 3 88 16 89 33 95 59 97 70 56 31 87 97 31 44 28 21 21 60 7 52 */
bzoj2168 半平面交模板 onleft的>=改成>就过了
//Serene #include<algorithm> #include<iostream> #include<cstring> #include<cstdlib> #include<cstdio> #include<cmath> using namespace std; #define ll long long #define db double #define For(i,a,b) for(int i=(a);i<=(b);++i) #define Rep(i,a,b) for(int i=(a);i>=(b);--i) const int maxn=500+7; int n,m; char cc; ll ff; template<typename T>void read(T& aa) { aa=0;cc=getchar();ff=1; while((cc<'0'||cc>'9')&&cc!='-') cc=getchar(); if(cc=='-') ff=-1,cc=getchar(); while(cc>='0'&&cc<='9') aa=aa*10+cc-'0',cc=getchar(); aa*=ff; } const db eps=1e-8; inline int dcmp(const db &x){return fabs(x)<eps? 0 : (x>0? 1:-1);} struct Xl{ db x,y; Xl(db x=0.,db y=0.):x(x),y(y){} Xl operator + (const Xl& b)const{return Xl(x+b.x,y+b.y);} Xl operator - (const Xl& b)const{return Xl(x-b.x,y-b.y);} Xl operator * (const db& b)const{return Xl(b*x,b*y);} db angle()const{return atan2(y,x);} }node[maxn],p[maxn]; db D_(const Xl& a,const Xl& b){return a.x*b.x+a.y*b.y;} db X_(const Xl& a,const Xl& b){return a.x*b.y-a.y*b.x;} struct Li{ Xl u,v; db rad; Li(Xl u=Xl(0.,0.),Xl v=Xl(0.,0.)):u(u),v(v){rad=(v-u).angle();} bool operator < (const Li& b)const{return rad<b.rad;} }li[maxn],zz[maxn]; int tot; Xl lijd(const Li& a,const Li& b) { Xl x=a.v-a.u,y=b.v-b.u,z=a.u-b.u; db c=X_(z,y); db t=c/(c+X_(y,x+z)); return a.u+x*t; } bool onleft(const Xl& a,const Li& l){return X_(l.v-l.u,a-l.u)>0;} db HPI() { sort(li+1,li+tot+1); int s=1,t=0; zz[++t]=li[1]; For(i,2,tot) { while(s<t&&!onleft(p[t-1],li[i])) t--; while(s<t&&!onleft(p[s],li[i])) s++; zz[++t]=li[i]; if(!dcmp(X_(li[i].v-li[i].u,zz[t-1].v-zz[t-1].u))) if(onleft(li[i].u,zz[--t])) zz[t]=li[i]; if(s<t) p[t-1]=lijd(zz[t],zz[t-1]); } while(s<t&&!onleft(p[t-1],li[s])) t--; if(t-s<2) return 0.000; db rs=0; p[t]=lijd(zz[t],zz[s]); For(i,s+1,t-1) rs+=X_(p[i]-p[s],p[i+1]-p[s]); return fabs(rs)/2.0; } int main() { read(n); ll x,y; For(i,1,n) { read(m); For(j,1,m) { read(x); read(y); node[j]=Xl(x,y); } For(j,1,m-1) li[++tot]=Li(node[j],node[j+1]); li[++tot]=Li(node[m],node[1]); } printf("%.3f",HPI()); return 0; }
bzoj4445 半平面交
简单推推式子
bzoj上精度坑人,最好开long double
//Serene #include<algorithm> #include<iostream> #include<cstring> #include<cstdlib> #include<cstdio> #include<cmath> using namespace std; #define ll long long #define db long double #define For(i,a,b) for(int i=(a);i<=(b);++i) #define Rep(i,a,b) for(int i=(a);i>=(b);--i) const int maxn=2e5+7; const db PI=acos(-1),eps=1e-18; int n; char cc; ll ff; template<typename T>void read(T& aa) { aa=0;cc=getchar();ff=1; while((cc<'0'||cc>'9')&&cc!='-') cc=getchar(); if(cc=='-') ff=-1,cc=getchar(); while(cc>='0'&&cc<='9') aa=aa*10+cc-'0',cc=getchar(); aa*=ff; } db pf(const db& x){return x*x;} inline int dcmp(const long double &x){return abs(x)<eps? 0:(x>0? 1:-1);} struct Xl{ db x,y; Xl(db x=0.,db y=0.):x(x),y(y){} Xl operator + (const Xl& b)const{return Xl(x+b.x,y+b.y);} Xl operator - (const Xl& b)const{return Xl(x-b.x,y-b.y);} Xl operator * (const db& b)const{return Xl(x*b,y*b);} db len()const{return sqrt(pf(x)+pf(y));} double angle() const{return atan2(y,x);} }node[maxn],p[maxn],R; db sr; double D_(const Xl& a,const Xl& b){return a.x*b.x+a.y*b.y;} double X_(const Xl& a,const Xl& b){return a.x*b.y-a.y*b.x;} struct Li{ Xl u,v;db rad; Li(Xl u=Xl(0.,0.),Xl v=Xl(0.,0.)):u(u),v(v){rad=(v-u).angle();} bool operator < (const Li &b) const{return rad<b.rad;} }li[maxn],zz[maxn]; int tot; bool onleft(const Xl& a,const Li& l) { return dcmp(X_(l.v-l.u,a-l.u))>0;} bool get_line(Xl R,db c) { db a=-R.y,b=R.x,x,y; if(dcmp(a)==0&&dcmp(b)==0) { if(dcmp(c)>=0) return 1; return 0; } Xl u,v; if(dcmp(a)==0) { y=-c/b; u=Xl(0,y); v=Xl(1,y); if(dcmp(b)<0) swap(u,v); } else if(dcmp(b)==0) { x=-c/a; u=Xl(x,1); v=Xl(x,0); if(dcmp(a)<0) swap(u,v); } else { x=0; y=-c/b; u=Xl(x,y); x=1; y=(-a-c)/b; v=Xl(x,y); if(dcmp(b)<0) swap(u,v); } li[++tot]=Li(u,v); return 1; } Xl lijd(const Li& a,const Li& b) { Xl x=a.v-a.u,y=b.v-b.u,z=a.u-b.u; db S=X_(z,y); db t=S/(S+X_(y,x+z)); return a.u+x*t; } db HPL() { sort(li+1,li+tot+1); int s=1,t=0; zz[++t]=li[1]; For(i,2,tot) { while(s<t&&!onleft(p[t-1],li[i])) t--; while(s<t&&!onleft(p[s],li[i])) s++; zz[++t]=li[i]; if(dcmp(X_(li[i].v-li[i].u,zz[t-1].v-zz[t-1].u))==0) if(onleft(li[i].u,zz[--t])) zz[t]=li[i]; if(s<t) p[t-1]=lijd(zz[t],zz[t-1]); } while(s<t&&!onleft(p[t-1],li[1])) t--; if(s+1>t-1) return 0.0; db rs=0; p[t]=lijd(zz[t],zz[s]); For(i,s+1,t-1) rs+=X_(p[i]-p[s],p[i+1]-p[s]); return dcmp(rs)==0? 0:fabs(rs)/2.0; } db work() { For(i,2,n-1) { R=node[1]-node[2]-node[i]+node[i+1]; sr=X_(node[i],node[i+1])-X_(node[1],node[2]); if(!get_line(R,sr)) return 0.0; } R=node[1]-node[2]-node[n]+node[1]; sr=X_(node[n],node[1])-X_(node[1],node[2]); if(!get_line(R,sr)) return 0.0; return HPL(); } int main() { read(n); ll x,y; For(i,1,n) { read(x); read(y); node[i]=Xl(x,y); } db S=0; For(i,2,n-1) S+=X_(node[i]-node[1],node[i+1]-node[1]); S=fabs(S)/2.0; For(i,1,n-1) li[i]=Li(node[i],node[i+1]); li[n]=Li(node[n],node[1]); tot=n; printf("%.4Lf",work()/S); return 0; } /* 3 0 0 1 0 0 1 */
hdu3007 最小圆覆盖模板
数学一本通上的问题分析要把人绕晕,看代码还是要直观些
这个做法复杂度是O(n)的,当加入圆的顺序随机时,因为三点定一圆,所以不在圆内概率是3/i
所以可以用random_shuffle()这个函数让点随机(要用到algorithm这个头文件)
//Serene #include<algorithm> #include<iostream> #include<cstring> #include<cstdlib> #include<cstdio> #include<cmath> using namespace std; #define ll long long #define db double #define For(i,a,b) for(int i=(a);i<=(b);++i) #define Rep(i,a,b) for(int i=(a);i>=(b);--i) const int maxn=500+7; int n; char cc; ll ff; template<typename T>void read(T& aa) { aa=0;cc=getchar();ff=1; while((cc<'0'||cc>'9')&&cc!='-') cc=getchar(); if(cc=='-') ff=-1,cc=getchar(); while(cc>='0'&&cc<='9') aa=aa*10+cc-'0',cc=getchar(); aa*=ff; } const db eps=1e-10; db pf(db x){return x*x;} inline int dcmp(const long db &x){return abs(x)<eps? 0:(x>0? 1:-1);} struct Xl{ db x,y; Xl(db x=0.,db y=0.):x(x),y(y){} Xl operator + (const Xl& b)const{return Xl(x+b.x,y+b.y);} Xl operator - (const Xl& b)const{return Xl(x-b.x,y-b.y);} Xl operator * (const db& b)const{return Xl(x*b,y*b);} Xl operator / (const db& b)const{return Xl(x/b,y/b);} db len()const{return sqrt(pf(x)+pf(y));} }node[maxn],c; db D_(const Xl& a,const Xl& b){return a.x*b.x+a.y*b.y;} db X_(const Xl& a,const Xl& b){return a.x*b.y-a.y*b.x;} db r; Xl cir_enter(const Xl& a,const Xl& b,const Xl& c) { db bx=b.x-a.x,by=b.y-a.y,bpf=pf(bx)+pf(by); db cx=c.x-a.x,cy=c.y-a.y,cpf=pf(cx)+pf(cy); db area=2*(bx*cy-by*cx); db x=(cy*bpf-by*cpf)/area; db y=(bx*cpf-cx*bpf)/area; return Xl(x+a.x,y+a.y); } void mincir() { random_shuffle(node+1,node+n+1); c=node[1]; r=0; For(i,2,n) if(dcmp((node[i]-c).len()-r)>0) { c=node[i]; r=0; For(j,1,i-1) if(dcmp((node[j]-c).len()-r)>0) { c=(node[i]+node[j])/2; r=(node[i]-c).len(); For(k,1,j-1) if(dcmp((node[k]-c).len()-r)>0) { c=cir_enter(node[i],node[j],node[k]); r=(node[i]-c).len(); } } } } int main() { read(n);db x,y; while(n) { For(i,1,n) { scanf("%lf%lf",&x,&y); node[i]=Xl(x,y); } mincir(); printf("%.2f %.2f %.2f\n",c.x,c.y,r); read(n); } return 0; } /* 3 1.00 1.00 3.00 0.00 3.00 2.00 0 */
bzoj4750 妖怪
上凸包
这道题我的思维方式有点奇怪
对于每个妖怪我们看做一个点$(x,y)$
那么它的战斗力为$\frac{a+b}{a} \times x + \frac{a+b}{b} \times y$
我们考虑对于一个固定的$(a,b)$,战斗力最大的是上凸壳上的一个点
现在反过来考虑对于一个点$(x_0,y_0)$,什么时候它是战斗力最大的
令$r=\frac{b}{a}$
那么有$r \times x_0 + y_0 > r \times x_1 + y_1 $
1) $x_0 > x_1$:$r \geq - \frac{y_0-y_1}{x_0-x_1}$
2) $x_0 < x_1$:$r \leq - \frac{y_1-y_0}{x_1-x_0}$
我们发现,对于上凸壳的一个点a,假如它左边的点是b,右边的点是c
那么$r$就是在$-angle(a-b)$到$-angle(c-a)$之间,并且$r \geq 0$(因为$a,b \geq 0$)
对于所有点中一个战斗力最大的那个点$(x,y)$,什么样的$(a,b)$可以使得它的战斗力尽量小呢
$\frac{a+b}{a} \times x = \frac{a+b}{b} \times y$
$\frac{b}{a} = \frac{y}{x}$
//Serene #include<algorithm> #include<iostream> #include<cstring> #include<cstdlib> #include<cstdio> #include<cmath> using namespace std; #define ll long long #define db double #define For(i,a,b) for(int i=(a);i<=(b);++i) #define Rep(i,a,b) for(int i=(a);i>=(b);--i) const int maxn=1e6+7; const db INF=1e16; int n; db ans=INF; char cc; ll ff; template<typename T>void read(T& aa) { aa=0;cc=getchar();ff=1; while((cc<'0'||cc>'9')&&cc!='-') cc=getchar(); if(cc=='-') ff=-1,cc=getchar(); while(cc>='0'&&cc<='9') aa=aa*10+cc-'0',cc=getchar(); aa*=ff; } const db eps=1e-12; int dcmp(const db &x){return fabs(x)<=eps? 0:(x<0? -1:1);} struct Xl{ db x,y; Xl(db x=0.,db y=0.):x(x),y(y){} Xl operator - (const Xl& b) const{return Xl(x-b.x,y-b.y);} bool operator < (const Xl& b) const{return dcmp(x-b.x)==0? y<b.y:x<b.x;} db angle() const{return y/x;} }node[maxn],zz[maxn]; db X_(const Xl& a,const Xl& b) {return a.x*b.y-a.y*b.x;} bool onleft(const Xl& p,const Xl& u,const Xl& v) { return dcmp(X_(v-u,p-u))>=0; } db p[maxn]; db get_ans(const Xl& p,db l,db r) { l=-l; r=-r; l=max(l,0.0); r=max(r,0.0); db t=p.y/p.x; if(t>r) t=r; else if(t<l) t=l; if(dcmp(t)==0) return INF; db a=100,b=a*t; return (a+b)/a*p.x+(a+b)/b*p.y; } void solve() { sort(node+1,node+n+1); int t=0; For(i,1,n) { while(t&&node[i].x==zz[t].x) t--; while(t>1&&onleft(node[i],zz[t-1],zz[t])) t--; zz[++t]=node[i]; } For(i,1,t-1) p[i]=(zz[i+1]-zz[i]).angle(); p[t]=-INF; For(i,1,t) ans=min(ans,get_ans(zz[i],p[i-1],p[i])); } int main() { read(n); db x,y; For(i,1,n) { read(x); read(y); node[i]=Xl(x,y); } solve(); printf("%.4f",ans); return 0; }
bzoj1069 最大土地面积
没有什么脑子的模板题,一直以为这道题很难,然后今天看了看数据范围……
//Serene #include<algorithm> #include<iostream> #include<cstring> #include<cstdlib> #include<cstdio> #include<cmath> using namespace std; #define ll long long #define db double #define For(i,a,b) for(int i=(a);i<=(b);++i) #define Rep(i,a,b) for(int i=(a);i>=(b);--i) const int maxn=4000+7; int n; db ans[maxn][maxn]; char cc; ll ff; template<typename T>void read(T& aa) { aa=0;cc=getchar();ff=1; while((cc<'0'||cc>'9')&&cc!='-') cc=getchar(); if(cc=='-') ff=-1,cc=getchar(); while(cc>='0'&&cc<='9') aa=aa*10+cc-'0',cc=getchar(); aa*=ff; } const db eps=1e-10; db pf(db x){return x*x;} int dcmp(db x){return fabs(x)<=eps? 0:(x>0? 1:-1);} struct Xl{ db x,y; Xl(db x=0.,db y=0.):x(x),y(y){} bool operator < (const Xl& b) const{return dcmp(x-b.x)==0? y<b.y:x<b.x;} Xl operator + (const Xl& b) const{return Xl(x+b.x,y+b.y);} Xl operator - (const Xl& b) const{return Xl(x-b.x,y-b.y);} db len() const{return pf(x)+pf(y);} }node[maxn],zz[maxn]; int t; db D_(const Xl& a,const Xl& b) {return a.x*b.x+a.y*b.y;} db X_(const Xl& a,const Xl& b) {return a.x*b.y-a.y*b.x;} bool cmp(const Xl& a,const Xl& b) { return dcmp(X_(a-node[1],b-node[1]))==0? (a-node[1]).len()<(b-node[1]).len() : dcmp(X_(a-node[1],b-node[1]))>0; } void graham() { For(i,2,n) if(node[i]<node[1]) swap(node[i],node[1]); sort(node+2,node+n+1,cmp); zz[++t]=node[1]; zz[++t]=node[2]; For(i,3,n) { while(t>1&&dcmp(X_(zz[t]-zz[t-1],node[i]-zz[t])<=0)) --t; zz[++t]=node[i]; } } void RC() { For(i,1,t) zz[i+t]=zz[i]; zz[0]=zz[t]; int now; Xl p; For(i,1,t) { now=(i+1)%t; For(j,i+1,i+t-1) { p=zz[j]-zz[i]; while(X_(p,zz[now]-zz[i])<X_(p,zz[now+1]-zz[i])) now=(now+1)%t; ans[i][j]=X_(p,zz[now]-zz[i]); } } } int main() { read(n); db x,y; For(i,1,n) { scanf("%lf%lf",&x,&y); node[i]=Xl(x,y); } graham(); RC(); db totans=0; For(i,1,t) For(j,i+1,t) totans=max(totans,ans[i][j]+ans[j][i+t]); printf("%.3lf",totans/2); return 0; }
bzoj2388 旅行规划
分块+凸包+神奇的三分
//Serene #include<algorithm> #include<iostream> #include<cstring> #include<cstdlib> #include<cstdio> #include<cmath> using namespace std; #define ll long long #define db double #define For(i,a,b) for(int i=(a);i<=(b);++i) #define Rep(i,a,b) for(int i=(a);i>=(b);--i) const int maxn=2e5+7,maxt=500+7; const ll INF=1e16; int n,m,sz,tot; ll s[maxn],A[maxn],B[maxn]; char cc; ll ff; template<typename T>void read(T& aa) { aa=0;cc=getchar();ff=1; while((cc<'0'||cc>'9')&&cc!='-') cc=getchar(); if(cc=='-') ff=-1,cc=getchar(); while(cc>='0'&&cc<='9') aa=aa*10+cc-'0',cc=getchar(); aa*=ff; } const db eps=1e-8; int dcmp(db x) {return fabs(x)<=eps? 0:(x<0? -1:1);} struct Xl{ db x,y; Xl(db x=0.,db y=0.):x(x),y(y){} Xl operator - (const Xl& b) const{return Xl(x-b.x,y-b.y);} bool operator < (const Xl& b) const{return dcmp(x-b.x==0)? y<b.y : x<b.x;} }node[maxn],p[maxt][maxt]; int top[maxt]; db X_(const Xl& a,const Xl& b) {return a.x*b.y-a.y*b.x;} db D_(const Xl& a,const Xl& b) {return a.x*b.x+a.y*b.y;} void get_ham(int pos) { int l=max(pos*sz,1),r=min((pos+1)*sz-1,n); For(i,l,r) node[i]=Xl(i,s[i]); int t=0; p[pos][++t]=node[l]; For(i,l+1,r) { while(t>1&&dcmp(X_(p[pos][t]-p[pos][t-1],node[i]-p[pos][t-1]))>=0) t--; p[pos][++t]=node[i]; } top[pos]=t; } void pd(int pos) { if(A[pos]==0&&B[pos]==0) return; int l=max(pos*sz,1),r=min((pos+1)*sz-1,n); For(i,l,r) s[i]+=A[pos]+(i-l)*B[pos]; A[pos]=B[pos]=0; } void chge(int l,int r,ll x) { pd(l/sz); pd(r/sz); if(l/sz==r/sz) { For(i,l,r) s[i]+=x*(i-l+1); get_ham(l/sz); For(i,r+1,min((r/sz+1)*sz-1,n)) s[i]+=x*(r-l+1); For(i,l/sz+1,tot) A[i]+=(r-l+1)*x; return; } For(i,l,min((l/sz+1)*sz-1,n)) s[i]+=x*(i-l+1); get_ham(l/sz); For(i,max((r/sz)*sz,1),r) s[i]+=x*(i-l+1); For(i,r+1,min((r/sz+1)*sz-1,n)) s[i]+=x*(r-l+1); get_ham(r/sz); For(i,l/sz+1,r/sz-1) A[i]+=(i*sz-l+1)*x,B[i]+=x; For(i,r/sz+1,tot) A[i]+=(r-l+1)*x; } ll cal(int pos) {return A[pos/sz]+B[pos/sz]*(pos-(pos/sz)*sz)+s[pos];} ll get_ans(int pos) { int l=1,r=top[pos],mid; ll a,b,c; while(l<=r) { if(l+1>=r) return max(cal(p[pos][l].x),cal(p[pos][r].x)); mid=(l+r)>>1; a=cal(p[pos][mid-1].x); b=cal(p[pos][mid].x); c=cal(p[pos][mid+1].x); if(a<b&&b<c) l=mid+1; else if(a>b&&b>c) r=mid-1; else return b; } return 0; } ll Yth(int l,int r) { pd(l/sz); pd(r/sz); ll rs=-INF; if(l/sz==r/sz) { For(i,l,r) rs=max(rs,s[i]); return rs; } For(i,l,min((l/sz+1)*sz-1,n)) rs=max(rs,s[i]); For(i,max((r/sz)*sz,1),r) rs=max(rs,s[i]); For(i,l/sz+1,r/sz-1) rs=max(rs,get_ans(i)); return rs; } int main() { read(n); sz=min((int)sqrt(n),200)+1; tot=n/sz; ll op,l,r,x; For(i,1,n) read(s[i]),s[i]+=s[i-1]; For(i,0,tot) get_ham(i); read(m); For(i,1,m) { read(op); read(l); read(r); if(op) printf("%lld\n",Yth(l,r)); else { read(x); chge(l,r,x); } } return 0; }
bzoj1038瞭望塔
半平面交+线段最小距离 板子
//Serene #include<algorithm> #include<iostream> #include<cstring> #include<cstdlib> #include<cstdio> #include<cmath> using namespace std; #define ll long long #define db double #define For(i,a,b) for(int i=(a);i<=(b);++i) #define Rep(i,a,b) for(int i=(a);i>=(b);--i) const int maxn=1000+7; int n; char cc; ll ff; template<typename T>void read(T& aa) { aa=0;cc=getchar();ff=1; while((cc<'0'||cc>'9')&&cc!='-') cc=getchar(); if(cc=='-') ff=-1,cc=getchar(); while(cc>='0'&&cc<='9') aa=aa*10+cc-'0',cc=getchar(); aa*=ff; } const db eps=1e-14,INF=1e30; int dcmp(const db x){return fabs(x)<=eps? 0:(x>0? 1:-1);} db pf(const db x){return x*x;} struct Xl{ db x,y; Xl(db x=0.,db y=0.):x(x),y(y){} Xl operator + (const Xl& b) const{return Xl(x+b.x,y+b.y);} Xl operator - (const Xl& b) const{return Xl(x-b.x,y-b.y);} Xl operator * (const db& b) const{return Xl(x*b,y*b);} Xl operator / (const db& b) const{return Xl(x/b,y/b);} db angle() const{return atan2(y,x);} db len() const{return sqrt(pf(x)+pf(y));} }node[maxn],p[maxn]; db X_(const Xl& a,const Xl& b) {return a.x*b.y-a.y*b.x;} db D_(const Xl& a,const Xl& b) {return a.x*b.x+a.y*b.y;} struct Li{ Xl u,v; db rad; Li(){} Li(Xl u,Xl v):u(u),v(v){rad=(v-u).angle();} bool operator < (const Li& b) const{return rad<b.rad;} }li[maxn],zz[maxn]; bool onleft(const Xl& a,const Li& l) {return dcmp(X_(l.v-l.u,a-l.u))>0;} Xl lijd(const Li& a,const Li& b) { Xl x=a.v-a.u,y=b.v-b.u,z=a.u-b.u; db t=X_(z,y); t=t/(t+X_(y,x+z)); return a.u+x*t; } db xdjl(Li a,Li b) { return min((a.u-b.u).len(),(a.v-b.v).len()); } db get_ans(db l,db r,Li a,Li b) { if(l>r) return INF; Xl x=a.v-a.u,y=b.v-b.u; a.u=a.u+x*((l-a.u.x)/(a.v.x-a.u.x)); b.u=b.u+y*((l-b.u.x)/(b.v.x-b.u.x)); x=a.v-a.u; y=b.v-b.u; a.v=a.u+x*((r-a.u.x)/(a.v.x-a.u.x)); b.v=b.u+y*((r-b.u.x)/(b.v.x-b.u.x)); return xdjl(a,b); } db HPI() { sort(li+1,li+n); int s=1,t=0; zz[++t]=li[1]; For(i,2,n-1) { while(s<t&&!onleft(p[t-1],li[i])) t--; while(s<t&&!onleft(p[s],li[i])) s++; zz[++t]=li[i]; if(s<t&&dcmp(li[i].rad-zz[t-1].rad)==0) if(onleft(li[i].u,zz[--t])) zz[t]=li[i]; if(s<t) p[t-1]=lijd(zz[t],zz[t-1]); } while(s<t&&!onleft(p[t-1],li[s])) t--; p[s-1].x=-INF; p[t+1].x=INF; db rs=INF,l,r; int pos=s; For(i,1,n-1) { while(pos<t&&p[pos].x<node[i].x) ++pos; l=max(node[i].x,p[pos-1].x); r=min(node[i+1].x,p[pos].x); rs=min(rs,get_ans(l,r,zz[pos],Li(node[i],node[i+1]))); while(pos<t&&p[pos].x<=node[i+1].x) { ++pos; l=max(node[i].x,p[pos-1].x); r=min(node[i+1].x,p[pos].x); rs=min(rs,get_ans(l,r,zz[pos],Li(node[i],node[i+1]))); } } return rs; } int main() { read(n); For(i,1,n) scanf("%lf",&node[i].x); For(i,1,n) scanf("%lf",&node[i].y); For(i,1,n-1) li[i]=Li(node[i],node[i+1]); printf("%.3f",HPI()); return 0; }