[模板]旋转卡壳
旋转卡壳就是对凸包用两条平行线夹着，然后绕着凸包的定点和边旋转，在旋转的过程中我们可以得到凸包的宽，直径等特征值，当然还可以计算其他东西。
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 |
#include <bits/stdc++.h> using namespace std; const int maxn=100+10; const double eps=1e-8; int dcmp(double x) { if (fabs(x)<eps) return 0; if (x>0) return 1; return -1; } struct point{ double x,y; point(){} point(double x,double y):x(x),y(y){} point operator + (const point & o) const{ return point(x+o.x,y+o.y); } point operator - (const point & o) const{ return point(x-o.x,y-o.y); } double operator ^ (const point & o) const { return x*o.y-y*o.x; } bool operator < (const point & o) const{ if (y==o.y) return x<o.x; return y>o.y; } void read() { scanf("%lf%lf",&x,&y); } double dis(const point & o) { return sqrt((x-o.x)*(x-o.x)+(y-o.y)*(y-o.y)); } } p[maxn]; double dis(point a,point b,point c) { double res= ((a-b)^(a-c))/b.dis(c); cerr<<res<<endl; return res; } bool cmp_x(const point &p,const point &q) { if (p.x!=q.x) return p.x<q.x; return p.y<q.y; } vector<point> convex_hull(point * ps,int n) { sort(ps,ps+n,cmp_x); int k=0; vector<point> qs(n*2); for (int i=0;i<n;i++) { while(k>1&&dcmp((qs[k-1]-qs[k-2])^(ps[i]-qs[k-1]))<=0) k--; qs[k++]=ps[i]; } int t=k; for (int i=n-2;i>=0;i--) { while(k>t&&dcmp((qs[k-1]-qs[k-2])^(ps[i]-qs[k-1]))<=0) k--; qs[k++]=ps[i]; } qs.resize(k-1); return qs; } int main() { int n; while(scanf("%d",&n)!=EOF) { for (int i=0;i<n;i++) { p[i].read(); } vector<point> pt=convex_hull(p,n); int l,r; int maxy,miny; maxy=-100005; miny=100005; for (int i=0;i<pt.size();i++) { if (pt[i].y>maxy) { maxy=pt[i].y; r=i; } if (pt[i].y<miny) { miny=pt[i].y; l=i; } } n=pt.size(); double ans=4*100000.0; int sl=l,sr=r; while(l!=sr||r!=sl) { int nl=(sl+1)%n; int nr=(sr+1)%n; if (dcmp((pt[nl]-pt[sl])^(pt[nr]-pt[sr]))<0) { ans=min(ans,dis(pt[sr],pt[sl],pt[nl])); sl=nl; } else { ans=min(ans,dis(pt[sl],pt[sr],pt[nr])); sr=nr; } } printf("%.10lf\n",ans); } return 0; } |
Read More »
越人歌
越人歌 今夕何夕兮搴洲中流 今日何日兮得与王子同舟 蒙羞被好兮不訾诟耻 心几烦而不绝兮得知王子 山有木兮木有枝 心说君兮君不知
Read More »滕王阁序
滕王阁序 唐·王勃 豫章故郡，洪都新府。星分翼轸，地接衡庐。襟三江而带五湖，控蛮荆而引瓯越。物华天宝，龙光射牛斗之墟；人杰地灵，徐孺下陈蕃之榻。雄州雾列，俊采星驰。台隍枕夷夏之交，宾主尽东南之美。都督阎公之雅望，棨戟遥临；宇文新州之懿范，襜帷暂驻。十旬休假，胜友如云；千里逢迎，高朋满座。腾蛟起凤，孟学士之词宗；紫电青霜，王将军之武库。家君作宰，路出名区；童子何知，躬逢胜饯。 时维九月，序属三秋。潦…
Read More »[模板]ntt和fft
FFT：
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 |
struct Complex { double real, image; Complex( double _real, double _image) { real = _real; image = _image; } Complex(){} }; Complex operator + (const Complex &c1, const Complex &c2) { return Complex(c1.real + c2.real, c1.image + c2.image); } Complex operator - (const Complex &c1, const Complex &c2) { return Complex(c1.real - c2.real, c1.image - c2.image); } Complex operator * (const Complex &c1, const Complex &c2) { return Complex(c1.real*c2.real - c1.image*c2.image, c1.real*c2.image + c1.image*c2.real); } int rev(int id, int len) { int ret = 0; for(int i = 0; (1 << i) < len; i++) { ret <<= 1; if(id & (1 << i)) ret |= 1; } return ret; } Complex A[1 << 19]; void FFT(Complex *a, int len, int DFT) { for(int i = 0; i < len; i++) A[rev(i, len)] = a[i]; for(int s = 1; (1 << s) <= len; s++) { int m = (1 << s); Complex wm = Complex(cos(DFT*2*PI/m), sin(DFT*2*PI/m)); for(int k = 0; k < len; k += m) { Complex w = Complex(1, 0); for(int j = 0; j < (m >> 1); j++) { Complex t = w*A[k + j + (m >> 1)]; Complex u = A[k + j]; A[k + j] = u + t; A[k + j + (m >> 1)] = u - t; w = w*wm; } } } if(DFT == -1) for(int i = 0; i < len; i++) A[i].real /= len, A[i].image /= len; for(int i = 0; i < len; i++) a[i] = A[i]; return; } |
NTT
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 |
#include <iostream> #include <cstdio> #include <cstring> #include <algorithm> #include <cmath> using namespace std; const int maxn=300020; const int P=880803841; const int G=26; int wn[25]; long long mul(long long x,long long y) { return (x*y-(long long)(x/(long double)P*y+1e-3)*P+P) % P; } int qpow(int x,int k,int p) { int ret=1; while(k) { if (k&1) ret=1LL*ret*x%p; k>>=1; x=1LL*x*x%p; } return ret; } void getwn() { for (int i=1;i<=21;i++) { int t=1<<i; wn[i]=qpow(G,(P-1)/t,P); } } void change(int *y,int len) { for (int i=1,j=len/2;i<len-1;i++) { if (i<j) swap(y[i],y[j]); int k=len/2; while(j>=k) { j-=k; k/=2; } j+=k; } } void NTT(int *y,int len,int on) { change(y,len); int id=0; for (int h=2;h<=len;h<<=1) { id++; for (int j=0;j<len;j+=h) { int w=1; for (int k=j;k<j+h/2;k++) { int u=y[k]; int t=1LL*y[k+h/2]*w%P; y[k]=u+t; if (y[k]>=P) y[k]-=P; y[k+h/2]=u-t+P; if (y[k+h/2]>=P) y[k+h/2]-=P; w=1LL*w*wn[id]%P; } } } if (on==-1) { for (int i=1;i<len/2;i++) { swap(y[i],y[len-i]); } int inv=qpow(len,P-2,P); for (int i=1;i<len;i++) { y[i]=1LL*y[i]*inv%P; } } } int main() { getwn(); return 0; } |
[模板]半平面交
给定一系列半平面，求其交集的面积 下面代码p为给定点集，用于得到半平面l
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 |
#include <iostream> #include <cstdio> #include <cstring> #include <algorithm> #include <cmath> using namespace std; const int maxn=1500+10; const double eps=1e-8; int n,pn,dq[maxn],top,bot; struct point{ double x,y; point(){} point(double x,double y):x(x),y(y){} void read() { scanf("%lf%lf",&x,&y); } void print() { printf("(%lf,%lf)\n",x,y); } point operator +(const point &o)const { return point(x+o.x,y+o.y); } point operator -(const point &o)const { return point(x-o.x,y-o.y); } point operator *(double k) const { return point(x*k,y*k); } double operator *(const point &o)const { return x*o.x+y*o.y; } double operator ^(const point &o)const { return x*o.y-y*o.x; } } p[maxn]; struct line{ point a,b; double angle; line(){} line(point oa,point ob) { a=oa; b=ob; angle=atan2(ob.y-oa.y,ob.x-oa.x); } line & operator =(line l) { a.x=l.a.x; a.y=l.a.y; b.x=l.b.x; b.y=l.b.y; angle=l.angle; return *this; } } l[maxn]; int dblcmp(double k) { if (fabs(k)<eps) return 0; return k>0?1:-1; } double multi(point p0,point p1,point p2) { return (p1-p0)^(p2-p0); } bool cmp(const line &l1,const line &l2) { int d=dblcmp(l1.angle-l2.angle); if (!d) { return dblcmp(multi(l1.a,l2.a,l2.b))<0; } return d<0; } void getIntersect(line l1,line l2,point &p) { double A1=l1.b.y-l1.a.y; double B1=l1.a.x-l1.b.x; double C1=(l1.b.x-l1.a.x)*l1.a.y-(l1.b.y-l1.a.y)*l1.a.x; double A2=l2.b.y-l2.a.y; double B2=l2.a.x-l2.b.x; double C2=(l2.b.x-l2.a.x)*l2.a.y-(l2.b.y-l2.a.y)*l2.a.x; p.x=(C2*B1-C1*B2)/(A1*B2-A2*B1); p.y=(C1*A2-C2*A1)/(A1*B2-A2*B1); } bool judge(line l0,line l1,line l2) { point p; getIntersect(l1,l2,p); return dblcmp(multi(p,l0.a,l0.b))>0; } void HalfPlaneIntersect() { int i,j; sort(l,l+n,cmp); for (i=0,j=0;i<n;i++) { if (dblcmp(l[i].angle-l[j].angle)>0) { l[++j]=l[i]; } } n=j+1; dq[0]=0; dq[1]=1; top=1; bot=0; for (int i=2;i<n;i++) { while(top>bot&&judge(l[i],l[dq[top]],l[dq[top-1]])) top--; while(top>bot&&judge(l[i],l[dq[bot]],l[dq[bot+1]])) bot++; dq[++top]=i; } while(top>bot&&judge(l[dq[bot]],l[dq[top]],l[dq[top-1]])) top--; while(top>bot&&judge(l[dq[top]],l[dq[bot]],l[dq[bot+1]])) bot++; dq[++top]=dq[bot]; for (pn=0,i=bot;i<top;i++,pn++) { getIntersect(l[dq[i+1]],l[dq[i]],p[pn]); } } double getArea() { if (pn<3) return 0; double area=0; for (int i=1;i<pn-1;i++) { //p[i].print(); area+=multi(p[0],p[i],p[i+1]); } if (area<0) area=-area; return area/2; } int main() { int t,i; scanf("%d",&t); while(t--) { scanf("%d",&n); for (i=0;i<n;i++) { p[i].read(); } for (i=0;i<n-1;i++) { l[i]=line(p[i],p[i+1]); } l[i]=line(p[i],p[0]); HalfPlaneIntersect(); printf("%.2lf\n",getArea()); } return 0; } |
Read More »
[模板]最大空凸包
最大空凸包指的是以给定点集的子集为定点的凸包，且该凸包内部不包含其他点，求其中面积最大的那个。 下面代码中的dot数组表示输入的点集，n为点集大小。
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 |
#include <iostream> #include <cstdio> #include <cstring> #include <cmath> #include <algorithm> using namespace std; const int maxn=55; const double zero=1e-8; struct Vector{ double x,y; }; inline Vector operator - (Vector a, Vector b) { Vector c; c.x=a.x-b.x; c.y=a.y-b.y; return c; } inline double Sqr(double a) { return a*a; } inline int Sign(double a) { if (fabs(a)<=zero) return 0; return a<0?-1:1; } inline bool operator < (Vector a,Vector b) { return Sign(b.y-a.y)>0||Sign(b.y-a.y)==0&&Sign(b.x-a.x)>0; } inline double Max(double a,double b) { return a>b ? a : b; } inline double Length(Vector a) { return sqrt(Sqr(a.x)+Sqr(a.y)); } inline double Cross(Vector a,Vector b) { return a.x*b.y-a.y*b.x; } Vector dot[maxn],List[maxn]; double opt[maxn][maxn]; int seq[maxn]; int n,len; double ans; bool Compare(Vector a,Vector b) { int temp=Sign(Cross(a,b)); if(temp!=0) return temp>0; temp=Sign(Length(b)-Length(a)); return temp>0; } void solve(int vv) { int t,i,j,_len; for (i=len=0;i<n;i++) { if (dot[vv]<dot[i]) List[len++]=dot[i]-dot[vv]; } for (i=0;i<len;i++) { for (j=0;j<len;j++) { opt[i][j]=0; } } sort(List,List+len,Compare); double v; for (t=1;t<len;t++) { _len=0; for (i=t-1;i>=0&&Sign(Cross(List[t],List[i]))==0;i--); while(i>=0) { v=Cross(List[i],List[t])/2; seq[_len++]=i; for (j=i-1;j>=0&&Sign(Cross(List[i]-List[t],List[j]-List[t]))>0;j--); if (j>=0) v+=opt[i][j]; ans=Max(ans,v); opt[t][i]=v; i=j; } for (i=_len-2;i>=0;i--) { opt[t][seq[i]]=Max(opt[t][seq[i]],opt[t][seq[i+1]]); } } } int main() { int T; scanf("%d",&T); while(T--) { scanf("%d",&n); for (int i=0;i<n;i++) { scanf("%lf%lf",&dot[i].x,&dot[i].y); } ans=0; for (int i=0;i<n;i++) { solve(i); } printf("%.1lf\n",ans); } return 0; } |
Read More »
插值？插值！
这是一个神奇的东西，mark！
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 |
#include <cstdio> #include <cstring> #include <cmath> #include <algorithm> #include <vector> #include <string> #include <map> #include <set> #include <cassert> using namespace std; #define rep(i,a,n) for (int i=a;i<n;i++) #define per(i,a,n) for (int i=n-1;i>=a;i--) #define pb push_back #define mp make_pair #define all(x) (x).begin(),(x).end() #define fi first #define se second #define SZ(x) ((int)(x).size()) typedef vector<int> VI; typedef long long ll; typedef pair<int,int> PII; const ll mod=1000000007; ll powmod(ll a,ll b) {ll res=1;a%=mod; assert(b>=0); for(;b;b>>=1){if(b&1)res=res*a%mod;a=a*a%mod;}return res;} // head int _; ll n; namespace linear_seq { const int N=10010; ll res[N],base[N],_c[N],_md[N]; vector<ll> Md; void mul(ll *a,ll *b,ll k) { rep(i,0,k+k) _c[i]=0; rep(i,0,k) if (a[i]) rep(j,0,k) _c[i+j]=(_c[i+j]+a[i]*b[j])%mod; for (int i=k+k-1;i>=k;i--) if (_c[i]) rep(j,0,SZ(Md)) _c[i-k+Md[j]]=(_c[i-k+Md[j]]-_c[i]*_md[Md[j]])%mod; rep(i,0,k) a[i]=_c[i]; } int solve(ll n,VI a,VI b) { ll ans=0,pnt=0; ll k=SZ(a); assert(SZ(a)==SZ(b)); rep(i,0,k) _md[k-1-i]=-a[i];_md[k]=1; Md.clear(); rep(i,0,k) if (_md[i]!=0) Md.push_back(i); rep(i,0,k) res[i]=base[i]=0; res[0]=1; while ((1ll<<pnt)<=n) pnt++; for (int p=pnt;p>=0;p--) { mul(res,res,k); if ((n>>p)&1) { for (int i=k-1;i>=0;i--) res[i+1]=res[i];res[0]=0; rep(j,0,SZ(Md)) res[Md[j]]=(res[Md[j]]-res[k]*_md[Md[j]])%mod; } } rep(i,0,k) ans=(ans+res[i]*b[i])%mod; if (ans<0) ans+=mod; return ans; } VI BM(VI s) { VI C(1,1),B(1,1); int L=0,m=1,b=1; rep(n,0,SZ(s)) { ll d=0; rep(i,0,L+1) d=(d+(ll)C[i]*s[n-i])%mod; if (d==0) ++m; else if (2*L<=n) { VI T=C; ll c=mod-d*powmod(b,mod-2)%mod; while (SZ(C)<SZ(B)+m) C.pb(0); rep(i,0,SZ(B)) C[i+m]=(C[i+m]+c*B[i])%mod; L=n+1-L; B=T; b=d; m=1; } else { ll c=mod-d*powmod(b,mod-2)%mod; while (SZ(C)<SZ(B)+m) C.pb(0); rep(i,0,SZ(B)) C[i+m]=(C[i+m]+c*B[i])%mod; ++m; } } return C; } int gao(VI a,ll n) { VI c=BM(a); c.erase(c.begin()); rep(i,0,SZ(c)) c[i]=(mod-c[i])%mod; return solve(n,c,VI(a.begin(),a.begin()+SZ(c))); } }; int main() { while(~scanf("%lld",&n)){ printf("%d\n",linear_seq::gao(VI{1,1, 5, 11, 36, 95, 281, 781, 2245, 6336, 18061, 51205, 145601, 413351, 1174500, 3335651, 9475901, 26915305, 76455961, 217172736},n)); } } |
Read More »
后赤壁赋 [苏轼]
后赤壁赋 [宋] 苏轼 是岁十月之望， 步自雪堂， 将归于临皋。 二客从予过黄泥之坂。 霜露既降， 木叶尽脱， 人影在地， 仰见明月， 顾而乐之， 行歌相答。 已而叹曰： “有客无酒， 有酒无肴， 月白风清， 如此良夜何！” 客曰： “今者薄暮， 举网得鱼， 巨口细鳞， 状如松江之鲈。 顾安所得酒乎？” 归而谋诸妇。 妇曰： “我有斗酒， 藏之久矣， 以待子不时之需。” 于是携酒与鱼， 复游于赤壁…
Read More »[HDU4441] Queue Sequence
Queue Sequence Time Limit: 8000/4000 MS (Java/Others) Memory Limit: 32768/32768 K (Java/Others) Problem Description There’s a queue obeying the first in first out rule. Each time you…
Read More »