Назад к лучшим решениям Status: AC, problem ZSQRT2, contest ZEL07. By bobber (Mike Obuhov), 2007-03-03 22:15:08.
  1. #include <math.h>
  2. #include <stdlib.h>
  3. #include <stdio.h>
  4. #include <memory.h>
  5. #define max(a,b)(((a)>(b))?(a):(b))
  6. #define PI 3.14159265358979
  7. #define MX 1048576
  8. #define T double
  9. #define I int
  10. #define V void
  11. #define K const
  12. #define R return
  13. #define MS memset
  14. #define MC memcpy
  15. #define ST sizeof(T)
  16. #define W while
  17. T B=10000.;
  18. T iB=0.0001;
  19. struct L{L(){l=0;m=new T[MX];MS(m,0,MX*ST);}T*m;I l;};
  20. V P(L&x){I i;printf("%d.",(I)x.m[0]);for(i=1;i<x.l;++i){printf("%.4d",(I)x.m[i]);}}
  21. V A(L&r,L&x,L&y){I c=0,i;r.l=max(x.l,y.l);for(i=r.l-1;i>=0;i--){r.m[i]=x.m[i]+y.m[i]+c;if(r.m[i]<B)c=0;else{c=1;r.m[i]-=B;}}}
  22. V S(L&r,L&x,I q){T b=B;I i=x.l-1;W(x.m[i]==0.)--i;x.l=i+1;r.l=x.l;for(;i>0;i--){r.m[i]=b-x.m[i];b=B-1;}r.m[0]=q-x.m[0]-1;}
  23. V M(L&r,L&x,I q){T c=0.,y;r.l=x.l;for(I i=x.l-1;i>=0;--i){y=x.m[i]*q;y+=c;if(y>=B){c=(I)(y*iB);y-=(c*B);}else c=0;r.m[i]=y;}}
  24. V D(L&r,L&x){T c=0.,y,q;r.l=x.l+1;for(I i=0;i<r.l;++i){y=x.m[i]+c*B;q=I(y*0.125);c=y-q*8;r.m[i]=q;}if(r.m[r.l-1]==0.)--r.l;}
  25. V M(L&r,L&x,L&y){I i,j;T k,t;r.l=x.l+y.l;MS(r.m,0,r.l*ST);for(j=y.l-1;j>=0;--j){if(y.m[j]==0)r.m[(r.l-(j+1))+x.l]=0;
  26. else{k=0;for(i=x.l-1;i>=0;--i){t=x.m[i]*y.m[j]+r.m[i+j]+k;if(t<B){k=0;r.m[i+j]=t;}else{k=I(t*iB);r.m[i+j]=t-k*B;}}
  27. r.m[(r.l-(j+1))+x.l]=k;}}}
  28. struct C{T r,i;C(){}C(T a,T b){r=a;i=b;}inline K C operator+(K C&c)K{R C(r+c.r,i+c.i);}
  29. inline K C operator-(K C&c)K{R C(r-c.r,i-c.i);}inline K C operator*(K C&c)
  30. K{R C(r*c.r-i*c.i,r*c.i+i*c.r);}inline K C operator/(K T&d)K{R C(r/d,i/d);}};
  31. inline K C cj(K C &c){R C(c.r,-c.i);}
  32. #define TV I tl,dx;I td;C pr,rt;
  33. #define IT(l,d)dx=0;tl=(l);td=(d);pr.r=1.0;pr.i=0.0;rt.r=sin(PI/((l)*2.0));rt.r=-2.0*rt.r*rt.r;rt.i=sin(PI/(l))*(d);
  34. #define NT if(((++dx)&15)==0){T ag=(PI*(dx))/tl;pr.r=sin(ag*0.5);pr.r=1.0-2.0*pr.r*pr.r;pr.i=sin(ag)*(td);}else{C t;t=pr;pr=pr*rt;pr=pr+t;}
  35. inline I rn(I r,I n){do{n=n>>1;r=r^n;}W((r&n)==0);R r;}
  36. V ro(C*d,I l){C t;if(l<=2)R;I r=0;for(I x=1;x<l;++x){r=rn(r,l);if(r>x){t=d[x];d[x]=d[r];d[r]=t;}}}
  37. V it(C*d,I n,I i){I s,h,b;TV;s=1;W(s<n){h=s;s*=2;IT(h,i);for(b=0;b<h;++b){I l,r;for(l=b;l<n;l+=s){C a,b;r=l+h;b=d[l];a=d[r];a=a*pr;d[l]=b+a;d[r]=b-a;}NT;}}}
  38. V ft(C*d,I n,I i){I k;TV;if(n<=(16384/sizeof(C))){it(d,n,i);R;}n/=2;IT(n,i);ft(d,n,i);ft(d+n,n,i);for(k=0;k<n;++k){C b,c;b=d[k];c=d[k+n]*pr;d[k]=b+c;d[k+n]=b-c;NT;}}
  39. V rf(T*t,I n,I u){I i,j;C*d=(C*)t;TV;n/=2;if(u>0){ro(d,n);ft(d,n,1);}IT(n,u);NT;for(i=1,j=n-i;i<n/2;i++,j--){C p1,p2,t;t=cj(d[j]);p1=d[i]+t;p2=d[i]-t;p2=p2*pr;
  40. t=C(-u*p2.i,u*p2.r);d[i]=p1-t;d[j]=p1+t;d[j]=cj(d[j]);d[i]=d[i]/2;d[j]=d[j]/2;NT;}{T r,i;r=d[0].r;i=d[0].i;d[0]=C(r+i,r-i);}if(u<0){d[0]=d[0]/2.0;ro(d,n);ft(d,n,-1);}}
  41. T X[MX],Y[MX],*Z=0;
  42. V cn(I l,L&r){T i=1.0/(l/2);T w,p,q=0.;T*c=r.m;I x;for(x=l-1;x>=0;--x){w=(c[x]*i)+q;p=floor(w+0.5);q=I(p*iB);c[x]=T(p-q*B);}x=l;do{x--;}W(c[x]==0.);r.l=x+1;}
  43. V rs(T*x,K T*y,I l){K C*a=(C*)x;K C*b=(C*)y;C*c=(C*)Z;Z[0]=x[0]*y[0];Z[1]=x[1]*y[1];for(I x=1;x<l/2;++x)c[x]=a[x]*b[x];}
  44. V F(L&r,L&x,L&y){Z=r.m;I l=2;W(l<x.l+y.l)l*=2;MC(X,x.m,x.l*ST);MS(&X[x.l],0,(l-x.l)*ST);rf(X,l,1);if(&x==&y)rs(X,X,l);else{MC(Y,y.m,y.l*ST);MS(&Y[y.l],0,(l-y.l)*ST);
  45. rf(Y,l,1);rs(X,Y,l);}rf(Z,l,-1);cn(l,r);}
  46. V U(L&r,L&x,L&y){if(x.l>24&&y.l>24)F(r,x,y);else M(r,x,y);}
  47. V E(L&a,L&b,L&c,L&r,L&d,L&e,L&f,I g){U(c,a,a);M(c,c,2);D(d,a);S(e,c,1);M(f,c,3);S(f,f,7);U(r,d,e);U(d,r,f);A(b,a,d);if(g<b.l){MS(&b.m[g],0,(b.l-g)*ST);b.l=g;}}
  48. I main(){K I G[]={2,5,14,27,78,230,687,2059,6174,18520,55558,166668};L a,b,c,r,d,e,f;a.l=2;a.m[1]=7000;
  49. for(I i=0;i<6;++i){E(a,b,c,r,d,e,f,G[2*i]);E(b,a,c,r,d,e,f,G[2*i+1]);}E(a,b,c,r,d,e,f,500002);M(r,b,2);P(r);R 0;}