Назад к лучшим решениям Status: AC, problem ZSQRT2, contest ZEL07. By jedi_knight_ (Ivan Popelyshev), 2007-02-11 14:18:27.
  1. {$IFDEF ONLINE_JUDGE}
  2. {$ASMMODE Intel}
  3. {$ELSE}
  4. {$APPTYPE CONSOLE}
  5. {$ENDIF}
  6. const
  7. base10=4;
  8. base=10000;
  9. ndigs2=20;
  10. ndigs=1 shl ndigs2;
  11. mask=1000;
  12. iter=2000000;
  13. type
  14. int=longint;
  15. real=double;
  16. tlong=array [0..ndigs+1] of int;
  17. tcomp=record re,im:real;end;
  18. tfur=array[0..ndigs-1] of tcomp;
  19. var
  20. a,b,c:tlong;
  21. aa,bb,cc: tfur;
  22. rev:array[0..ndigs-1]of int;
  23. len,len2:int;
  24. cl:int;
  25. procedure add(var a, b: tlong;var c: tlong);
  26. var
  27. i,t:int;
  28. begin
  29. t:=0;
  30. for i:=0 to 2*cl do
  31. begin
  32. t:=t+a[i]+b[i];
  33. c[i]:=t mod base;
  34. t:=t div base;
  35. end
  36. end;
  37.  
  38. procedure sub(var a,b:tlong;var c:tlong);
  39. var
  40. i,t:int;
  41. begin
  42. t:=0;
  43. for i:=0 to 2*cl do
  44. begin
  45. t:=t+a[i]-b[i];
  46. if (t<0) then
  47. begin
  48. c[i]:=t+base;
  49. t:=-1;
  50. end else
  51. begin
  52. c[i]:=t;
  53. t:=0;
  54. end
  55. end
  56. end;
  57. procedure calc;
  58. var i,j,k:int;
  59. begin
  60. j:=0;
  61. for i:=0 to len-2 do
  62. begin
  63. rev[i]:=j;
  64. k:=len;
  65. repeat
  66. k:=k shr 1;
  67. j:=j xor k;
  68. until (j and k)>0;
  69. end;
  70. rev[len-1]:=len-1;
  71. end;
  72. procedure copy1(var source:tlong; var dest: tfur);
  73. var
  74. i: int;
  75. begin
  76. for i:=0 to len-1 do
  77. begin
  78. dest[i].re:=source[rev[i]];
  79. dest[i].im:=0;
  80. end
  81. end;
  82. type
  83. tct=array of tcomp;
  84. var
  85. tb:array[0..21] of tct;
  86. procedure init;
  87. var
  88. i,j:int;
  89. sina,cosa:real;
  90. begin
  91. for i:=0 to 20 do
  92. begin
  93. setlength(tb[i],1 shl i);
  94. cosa:=cos(pi/(1 shl i));
  95. sina:=sin(pi/(1 shl i));
  96. for j:=0 to (1 shl i)-1 do
  97. begin
  98. if (j mod mask)=0 then
  99. begin
  100. tb[i][j].re:=cos((pi*j)/(1 shl i));
  101. tb[i][j].im:=sin((pi*j)/(1 shl i));
  102. continue;
  103. end;
  104. tb[i][j].re:=tb[i][j-1].re*cosa-tb[i][j-1].im*sina;
  105. tb[i][j].im:=tb[i][j-1].re*sina+tb[i][j-1].im*cosa;
  106. end
  107. end
  108. end;
  109. procedure fur(var res:tfur);
  110. var
  111. i,j,j1,k,l:int;
  112. c2:tcomp;
  113. cosa, sina: real;
  114. v:tct;
  115. begin
  116. for i:=1 to len2 do
  117. begin
  118. j:=1 shl (i-1);
  119. j1:=1 shl i;
  120. k:=0;
  121. v:=tb[i-1];
  122. while (k<len) do
  123. begin
  124. for l:=0 to j-1 do
  125. begin
  126. cosa:=v[l].re;
  127. sina:=v[l].im;
  128. c2.re:=res[k+l+j].re*cosa-res[k+l+j].im*sina;
  129. c2.im:=res[k+l+j].im*cosa+res[k+l+j].re*sina;
  130. res[k+l+j].re:=res[k+l].re-c2.re;
  131. res[k+l+j].im:=res[k+l].im-c2.im;
  132. res[k+l].re:=res[k+l].re+c2.re;
  133. res[k+l].im:=res[k+l].im+c2.im;
  134. end;
  135. inc(k,j1);
  136. end
  137. end
  138. end;
  139.  
  140. procedure mult1(var aa, bb, res: tfur);
  141. var
  142. i:int;
  143. c: tcomp;
  144. begin
  145. for i:=0 to len-1 do
  146. begin
  147. c.re:=(aa[i].re*bb[i].re-aa[i].im*bb[i].im);
  148. c.im:=(aa[i].im*bb[i].re+aa[i].re*bb[i].im);
  149. res[rev[i]].re:=c.re;
  150. res[rev[i]].im:=-c.im;
  151. end
  152. end;
  153. procedure round(var f: tfur; var res: tlong; s: int);
  154. var
  155. i:int;
  156. cf:int64;
  157. label L1;
  158. begin
  159. cf:=0;
  160. for i:=0 to len-1 do
  161. begin
  162. f[i].re:=f[i].re/len;
  163. cf:=cf+trunc((f[i].re)+0.4);
  164. asm
  165. mov eax,dword ptr cf+4
  166. xor edx,edx
  167. mov ecx,base
  168. div ecx
  169. mov dword ptr cf+4,eax
  170. mov eax,dword ptr cf
  171. div ecx
  172. mov dword ptr cf,eax
  173. mov eax,i
  174. sub eax,s
  175. jl L1
  176. mov ecx,res
  177. mov [ecx+4*eax],edx
  178. L1:
  179. end
  180. end
  181. end;
  182. procedure shift(var aa:tlong;sh:int);
  183. var
  184. i:int;
  185. begin
  186. if sh>0 then
  187. for i:=cl*2 downto 0 do
  188. begin
  189. a[i+sh]:=aa[i];aa[i]:=0
  190. end
  191. else
  192. for i:=0 to cl*2 do
  193. begin
  194. a[i]:=aa[i-sh];aa[i-sh]:=0
  195. end
  196. end;
  197. procedure solve;
  198. var
  199. i, k: int;
  200. s,s1:shortstring;
  201. begin
  202. fillchar(a,sizeof(a),0);
  203. fillchar(b,sizeof(b),0);
  204. fillchar(c,sizeof(c),0);
  205. a[0]:=7071;
  206. cl:=1;
  207. len:=2;len2:=1;
  208. while cl*base10<iter do
  209. begin
  210. if (cl=1024) then
  211. begin
  212. shift(a,-24);
  213. cl:=1000;
  214. end;
  215. while (cl*2>len) or (cl*3>len) and (cl<1000) do
  216. begin
  217. Inc(len2);
  218. len:=len*2;
  219. end;
  220. calc;
  221. copy1(a,aa);
  222. fur(aa);
  223. mult1(aa, aa, cc);
  224. fur(cc);
  225. round(cc,b,0);
  226. c[cl*2-1]:=base div 2;
  227. sub(c, b, b);
  228. c[cl*2-1]:=0;
  229. copy1(b,bb);
  230. fur(bb);
  231. mult1(aa,bb,cc);
  232. fur(cc);
  233. round(cc,b,cl);
  234. shift(a,cl);
  235. add(a,b,a);
  236. cl:=cl*2;
  237. end;
  238. add(a,a,a);
  239. i:=ndigs;
  240. while (a[i]=0) do dec(i);
  241. k:=i;
  242. s:='';s1:='';
  243. write('1.');
  244. while (i>=0) do
  245. begin
  246. str(a[i], s);
  247. if (i<k) then
  248. begin
  249. while (length(s)<base10) do
  250. s:='0'+s;
  251. s1:=s1+s;
  252. if length(s1)>100 then
  253. begin
  254. write(s1);
  255. s1:='';
  256. end
  257. end else
  258. begin
  259. delete(s, 1, 1);
  260. write(s);
  261. end;
  262. dec(i);
  263. end;
  264. writeln(s1);
  265. end;
  266.  
  267. begin
  268. {$ifndef ONLINE_JUDGE}
  269. assignfile(output, 'output.txt');
  270. rewrite(output);
  271. {$endif}
  272. init;
  273. solve;
  274. end.