Integers.txt 26 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700701702703704705706707708709710711712713714715716717718719720721722723724725726727728729730731732733734735736737738739740741742743744745746747748749750751752753754755756757758759760761762763764765766767768769770771772773774775776777778779780781782783784785786787788789790791792793794795796797798799800801802803804805806807808809810811812813814815816817818819820821822823824825826827828829830831832833834835836837838839840841842843844845846847848
  1. MODULE Integers;
  2. (* THIS IS TEXT COPY OF BlackBox 1.6-rc6 System/Mod/Integers.odc *)
  3. (* DO NOT EDIT *)
  4. IMPORT Files, Math;
  5. CONST
  6. B = 10000; DecPerDig = 4; BinBase = 16 * 1024;
  7. KaratsubaBreak = 41;
  8. TYPE
  9. Index = INTEGER;
  10. Digit = SHORTINT;
  11. DoubleDigit = INTEGER;
  12. IntegerDesc = ARRAY OF Digit; (* to hide internal structure from interface *)
  13. Integer* = POINTER TO IntegerDesc;
  14. Buffer = RECORD
  15. digit: Integer;
  16. beg, len: Index
  17. END;
  18. VAR zero, one, two, buf6: Integer;
  19. PROCEDURE CopyOf (x: Integer; len: Index): Integer;
  20. VAR buf: Integer;
  21. BEGIN
  22. ASSERT(len > 0, 20);
  23. NEW(buf, len);
  24. REPEAT DEC(len); buf[len] := x[len] UNTIL len = 0;
  25. RETURN buf
  26. END CopyOf;
  27. (* Operations on Digits *)
  28. PROCEDURE Add (x, y, sum: Integer; xL, yL: Index; OUT sumL: Index);
  29. VAR i, l: Index; c: Digit;
  30. BEGIN
  31. l := MIN(xL, yL);
  32. i := 0; c := 0;
  33. WHILE i < l DO c := SHORT(c DIV B + x[i] + y[i]); sum[i] := SHORT(c MOD B); INC(i) END;
  34. WHILE i < xL DO c := SHORT(c DIV B + x[i]); sum[i] := SHORT(c MOD B); INC(i) END;
  35. WHILE i < yL DO c := SHORT(c DIV B + y[i]); sum[i] := SHORT(c MOD B); INC(i) END;
  36. IF c >= B THEN sum[i] := SHORT(c DIV B); INC(i) END;
  37. sumL := i
  38. END Add;
  39. PROCEDURE Subtract (x, y, dif: Integer; xL, yL: Index; OUT difL: Index);
  40. VAR i: Index; c, d: Digit;
  41. BEGIN
  42. ASSERT(xL >= yL, 20);
  43. i := 0; difL := 0; c := 0;
  44. WHILE i < yL DO
  45. c := SHORT(c DIV B + x[i] - y[i]); d := SHORT(c MOD B);
  46. IF d # 0 THEN
  47. WHILE difL # i DO dif[difL] := 0; INC(difL) END;
  48. dif[i] := d; INC(difL)
  49. END;
  50. INC(i)
  51. END;
  52. WHILE i < xL DO
  53. c := SHORT(c DIV B + x[i]); d := SHORT(c MOD B);
  54. IF d # 0 THEN
  55. WHILE difL # i DO dif[difL] := 0; INC(difL) END;
  56. dif[i] := d; INC(difL)
  57. END;
  58. INC(i)
  59. END;
  60. ASSERT(c DIV B = 0, 100)
  61. END Subtract;
  62. PROCEDURE OneDigitMult (a, b: Buffer; VAR c: Buffer);
  63. VAR i: Index; carry, factor: DoubleDigit;
  64. BEGIN
  65. ASSERT(a.len = 1, 20);
  66. factor := a.digit[a.beg]; i := 0; carry := 0;
  67. WHILE i # b.len DO
  68. carry := carry DIV B + factor * b.digit[b.beg + i]; c.digit[c.beg + i] := SHORT(carry MOD B);
  69. INC(i)
  70. END;
  71. IF carry >= B THEN c.digit[c.beg + i] := SHORT(carry DIV B); INC(i) END;
  72. c.len := i
  73. END OneDigitMult;
  74. PROCEDURE SimpleMult (a, b: Buffer; VAR c: Buffer);
  75. VAR i, j, k: Index; c0, c1: DoubleDigit;
  76. BEGIN
  77. ASSERT(a.len <= b.len, 20);
  78. c.len := a.len + b.len - 1;
  79. i := 0; c0 := 0; c1 := 0;
  80. REPEAT
  81. IF i < b.len THEN
  82. IF i < a.len THEN j := i; k := 0 ELSE j := a.len - 1; k := i - a.len + 1 END;
  83. REPEAT
  84. c0 := c0 + a.digit[a.beg + j] * b.digit[b.beg + k];
  85. IF c0 > MAX(DoubleDigit) - BinBase * (B - 1) THEN
  86. c1 := c1 + c0 DIV BinBase; c0 := c0 MOD BinBase
  87. END;
  88. DEC(j); INC(k)
  89. UNTIL j < 0
  90. ELSE
  91. j := a.len - 1; k := i - a.len + 1;
  92. REPEAT
  93. c0 := c0 + a.digit[a.beg + j] * b.digit[b.beg + k];
  94. IF c0 > MAX(DoubleDigit) - BinBase * (B - 1) THEN
  95. c1 := c1 + c0 DIV BinBase; c0 := c0 MOD BinBase
  96. END;
  97. DEC(j); INC(k)
  98. UNTIL k = b.len
  99. END;
  100. IF c1 = 0 THEN c.digit[c.beg + i] := SHORT(c0 MOD B); c0 := c0 DIV B
  101. ELSE
  102. c0 := c0 + BinBase * (c1 MOD B);
  103. c.digit[c.beg + i] := SHORT(c0 MOD B); c0 := c0 DIV B; c1 := c1 DIV B
  104. END;
  105. INC(i)
  106. UNTIL i = c.len;
  107. IF c0 # 0 THEN c.digit[c.beg + c.len] := SHORT(c0); INC(c.len) END
  108. END SimpleMult;
  109. PROCEDURE AddBuf (a, b: Buffer; VAR c: Buffer); (* c := a + b *)
  110. VAR i: Index; carry: Digit;
  111. BEGIN
  112. ASSERT(a.len <= b.len, 20);
  113. i := 0; carry := 0;
  114. WHILE i # a.len DO
  115. carry := SHORT(carry DIV B + a.digit[a.beg + i] + b.digit[b.beg + i]);
  116. c.digit[c.beg + i] := SHORT(carry MOD B);
  117. INC(i)
  118. END;
  119. WHILE (i # b.len) & (carry >= B) DO
  120. carry := SHORT(carry DIV B + b.digit[b.beg + i]); c.digit[c.beg + i] := SHORT(carry MOD B);
  121. INC(i)
  122. END;
  123. IF carry >= B THEN c.digit[c.beg + i] := SHORT(carry DIV B); INC(i)
  124. ELSE
  125. WHILE i # b.len DO c.digit[c.beg + i] := b.digit[b.beg + i]; INC(i) END
  126. END;
  127. c.len := i
  128. END AddBuf;
  129. PROCEDURE AddToBuf (VAR a: Buffer; b: Buffer; shift: Index); (* a := a + b * B^shift *)
  130. VAR i, n: Index; carry: Digit;
  131. BEGIN
  132. b.beg := b.beg - shift; b.len := b.len + shift; i := shift; n := MIN(a.len, b.len); carry := 0;
  133. WHILE i # n DO
  134. carry := SHORT(carry DIV B + a.digit[a.beg + i] + b.digit[b.beg + i]);
  135. a.digit[a.beg + i] := SHORT(carry MOD B);
  136. INC(i)
  137. END;
  138. IF i # a.len THEN
  139. WHILE (i # a.len) & (carry >= B) DO
  140. carry := SHORT(carry DIV B + a.digit[a.beg + i]); a.digit[a.beg + i] := SHORT(carry MOD B);
  141. INC(i)
  142. END;
  143. IF carry >= B THEN a.digit[a.beg + i] := SHORT(carry DIV B); INC(i) END
  144. ELSE
  145. WHILE (i # b.len) & (carry >= B) DO
  146. carry := SHORT(carry DIV B + b.digit[b.beg + i]); a.digit[a.beg + i] := SHORT(carry MOD B);
  147. INC(i)
  148. END;
  149. IF carry >= B THEN a.digit[a.beg + i] := SHORT(carry DIV B); INC(i)
  150. ELSE
  151. WHILE i # b.len DO a.digit[a.beg + i] := b.digit[b.beg + i]; INC(i) END
  152. END
  153. END;
  154. a.len := MAX(i, a.len)
  155. END AddToBuf;
  156. PROCEDURE SubtractFromBuf (VAR a: Buffer; b, c: Buffer); (* a := a - b - c *)
  157. VAR i: Index; carry: Digit;
  158. BEGIN
  159. ASSERT(b.len <= c.len, 20);
  160. i := 0; carry := 0;
  161. WHILE i # b.len DO
  162. carry := SHORT(carry DIV B + a.digit[a.beg + i] - b.digit[b.beg + i] - c.digit[c.beg + i]);
  163. a.digit[a.beg + i] := SHORT(carry MOD B);
  164. INC(i)
  165. END;
  166. WHILE i # c.len DO
  167. carry := SHORT(carry DIV B + a.digit[a.beg + i] - c.digit[c.beg + i]);
  168. a.digit[a.beg + i] := SHORT(carry MOD B);
  169. INC(i)
  170. END;
  171. WHILE carry < 0 DO
  172. carry := SHORT(carry DIV B + a.digit[a.beg + i]); a.digit[a.beg + i] := SHORT(carry MOD B);
  173. INC(i)
  174. END;
  175. ASSERT(i <= a.len, 100);
  176. WHILE (a.len # 0) & (a.digit[a.beg + a.len - 1] = 0) DO DEC(a.len) END
  177. END SubtractFromBuf;
  178. PROCEDURE KStep (a, b: Buffer; VAR c: Buffer; stack: Buffer);
  179. VAR n2, i: Index; a0, a1, b0, b1, c0, c1, h: Buffer;
  180. BEGIN
  181. ASSERT(a.len <= b.len, 20);
  182. IF a.len = 0 THEN c.len := 0
  183. ELSIF a.len = 1 THEN OneDigitMult(a, b, c)
  184. ELSIF a.len <= KaratsubaBreak THEN SimpleMult(a, b, c)
  185. ELSE
  186. n2 := b.len DIV 2;
  187. c0.digit := c.digit; c0.beg := c.beg; c1.digit := c.digit; c1.beg := c.beg + 2 * n2;
  188. a0.digit := a.digit; a0.beg := a.beg; a0.len := MIN(a.len, n2);
  189. a1.digit := a.digit; a1.beg := a.beg + n2; a1.len := MAX(0, a.len - n2);
  190. WHILE (a0.len # 0) & (a0.digit[a0.beg + a0.len - 1] = 0) DO DEC(a0.len) END;
  191. b0.digit := b.digit; b0.beg := b.beg; b0.len := MIN(b.len, n2);
  192. b1.digit := b.digit; b1.beg := b.beg + n2; b1.len := MAX(0, b.len - n2);
  193. WHILE (b0.len # 0) & (b0.digit[b0.beg + b0.len - 1] = 0) DO DEC(b0.len) END;
  194. IF (a0.len # 0) OR (b0.len # 0) THEN
  195. IF a0.len <= a1.len THEN AddBuf(a0, a1, c1) ELSE AddBuf(a1, a0, c1) END;
  196. IF b0.len <= b1.len THEN AddBuf(b0, b1, c0) ELSE AddBuf(b1, b0, c0) END;
  197. h.digit := stack.digit; h.beg := stack.beg; stack.beg := stack.beg + c0.len + c1.len;
  198. IF c0.len <= c1.len THEN KStep(c0, c1, h, stack) ELSE KStep(c1, c0, h, stack) END;
  199. IF a0.len <= b0.len THEN KStep(a0, b0, c0, stack) ELSE KStep(b0, a0, c0, stack) END;
  200. KStep(a1, b1, c1, stack);
  201. IF c0.len <= c1.len THEN SubtractFromBuf(h, c0, c1) ELSE SubtractFromBuf(h, c1, c0) END;
  202. IF c1.len # 0 THEN
  203. i := c0.beg + c0.len;
  204. WHILE i < c1.beg DO c.digit[i] := 0; INC(i) END;
  205. c.len := c1.beg + c1.len - c.beg
  206. ELSE
  207. WHILE c0.len < n2 DO c0.digit[c0.beg + c0.len] := 0; INC(c0.len) END;
  208. c.len := c0.len
  209. END;
  210. ASSERT(h.len # 0, 100);
  211. AddToBuf(c, h, n2)
  212. ELSE
  213. KStep(a1, b1, c1, stack); c.len := c1.beg + c1.len - c.beg;
  214. i := c.beg;
  215. WHILE i # c1.beg DO c.digit[i] := 0; INC(i) END
  216. END
  217. END
  218. END KStep;
  219. PROCEDURE Karatsuba (x, y, pro:Integer; xL, yL: Index; OUT proL: Index);
  220. VAR a, b, c, stack: Buffer;
  221. BEGIN
  222. ASSERT(xL <= yL, 20);
  223. a.digit := x; a.beg := 0; a.len := xL; b.digit := y; b.beg := 0; b.len := yL;
  224. c.digit := pro; c.beg := 0;
  225. NEW(stack.digit, 2 * b.len); stack.beg := 0;
  226. KStep(a, b, c, stack);
  227. proL := c.len
  228. END Karatsuba;
  229. PROCEDURE Multiply (x, y, pro: Integer; xL, yL: Index; OUT proL: Index);
  230. VAR i, j, k: Index; c0, c1: DoubleDigit;
  231. BEGIN
  232. ASSERT(xL <= yL, 20);
  233. IF xL > KaratsubaBreak THEN Karatsuba(x, y, pro, xL, yL, proL)
  234. ELSIF xL = 1 THEN
  235. proL := 0; c1 := x[0]; c0 := 0;
  236. WHILE proL < yL DO
  237. c0 := c1 * y[proL] + c0; pro[proL] := SHORT(c0 MOD B);
  238. c0 := c0 DIV B ; INC(proL)
  239. END;
  240. IF c0 # 0 THEN pro[proL] := SHORT(c0); INC(proL) END
  241. ELSE
  242. proL := xL + yL - 1;
  243. i := 0; c0 := 0; c1 := 0;
  244. REPEAT
  245. IF i < yL THEN
  246. IF i < xL THEN j := i; k := 0 ELSE j := xL - 1; k := i - xL + 1 END;
  247. REPEAT
  248. c0 := c0 + x[j] * y[k];
  249. IF c0 > MAX(DoubleDigit) - BinBase * (B - 1) THEN
  250. c1 := c1 + c0 DIV BinBase; c0 := c0 MOD BinBase
  251. END;
  252. DEC(j); INC(k)
  253. UNTIL j < 0
  254. ELSE
  255. j := xL - 1; k := i - xL + 1;
  256. REPEAT
  257. c0 := c0 + x[j] * y[k];
  258. IF c0 > MAX(DoubleDigit) - BinBase * (B - 1) THEN
  259. c1 := c1 + c0 DIV BinBase; c0 := c0 MOD BinBase
  260. END;
  261. DEC(j); INC(k)
  262. UNTIL k = yL
  263. END;
  264. IF c1 = 0 THEN pro[i] := SHORT(c0 MOD B); c0 := c0 DIV B
  265. ELSE c0 := c0 + BinBase * (c1 MOD B); pro[i] := SHORT(c0 MOD B);
  266. c0 := c0 DIV B; c1 := c1 DIV B
  267. END;
  268. INC(i)
  269. UNTIL i = proL;
  270. IF c0 # 0 THEN pro[proL] := SHORT(c0); INC(proL) END
  271. END
  272. END Multiply;
  273. PROCEDURE DecomposeQuoRem (x, y: Integer; xL, yL: Index);
  274. VAR ix, iy, j: Index; d, q, h, yLead, ySecond: DoubleDigit; yBuf: Integer;
  275. BEGIN
  276. ASSERT((yL # 0) & (y[yL - 1] # 0), 20);
  277. IF yL = 1 THEN
  278. j := xL - 1; h := 0; d := y[0];
  279. WHILE j >= 0 DO h := x[j] + h * B; x[j + 1] := SHORT(h DIV d); h := h MOD d; DEC(j) END;
  280. x[0] := SHORT(h)
  281. ELSIF xL >= yL THEN
  282. x[xL] := 0; d := (B DIV 2 - 1) DIV y[yL - 1] + 1; yBuf := CopyOf(y, yL);
  283. IF d # 1 THEN
  284. j := 0; h := 0;
  285. WHILE j < xL DO h := d * x[j] + h DIV B; x[j] := SHORT(h MOD B); INC(j) END;
  286. x[xL] := SHORT(h DIV B);
  287. j := 0; h := 0;
  288. WHILE j < yL DO h := d * yBuf[j] + h DIV B; yBuf[j] := SHORT(h MOD B); INC(j) END;
  289. ASSERT(h DIV B = 0, 100)
  290. END;
  291. yLead := yBuf[yL - 1]; ySecond := yBuf[yL - 2]; j := xL;
  292. WHILE j >= yL DO
  293. IF x[j] # yLead THEN q := (x[j] * B + x[j - 1]) DIV yLead ELSE q := B - 1 END;
  294. WHILE ySecond * q > (x[j] * B + x[j - 1] - yLead * q) * B + x[j - 2] DO
  295. DEC(q)
  296. END;
  297. ix := j - yL; iy := 0; h := 0;
  298. WHILE iy < yL DO
  299. h := x[ix] - q * yBuf[iy] + h DIV B; x[ix] := SHORT(h MOD B); INC(ix); INC(iy)
  300. END;
  301. IF (-x[j]) # (h DIV B) THEN
  302. ix := j - yL; iy := 0; h := 0;
  303. WHILE iy < yL DO
  304. h := h DIV B + x[ix] + yBuf[iy]; x[ix] := SHORT(h MOD B); INC(ix); INC(iy)
  305. END;
  306. x[j] := SHORT(q - 1)
  307. ELSE x[j] := SHORT(q)
  308. END;
  309. DEC(j)
  310. END;
  311. IF d # 1 THEN
  312. j := yL; h := 0;
  313. WHILE j # 0 DO DEC(j); h := h + x[j]; x[j] := SHORT(h DIV d); h := (h MOD d) * B END
  314. END
  315. END
  316. END DecomposeQuoRem;
  317. PROCEDURE GetQuoRem (x, y: Integer; xL, yL: Index; xNeg, yNeg: BOOLEAN;
  318. quo, rem: Integer; OUT quoL, remL: Index; OUT quoNeg, remNeg: BOOLEAN;
  319. doQuo, doRem: BOOLEAN);
  320. VAR i: Index; c: Digit; xBuf: Integer;
  321. BEGIN
  322. ASSERT(xL >= yL, 20);
  323. xBuf := CopyOf(x, xL + 1);
  324. DecomposeQuoRem(xBuf, y, xL, yL);
  325. i := xL;
  326. WHILE (i >= yL) & (xBuf[i] = 0) DO DEC(i) END;
  327. quoL := i - yL + 1;
  328. i := yL - 1;
  329. WHILE (i >= 0) & (xBuf[i] = 0) DO DEC(i) END;
  330. remL := i + 1;
  331. IF doQuo THEN
  332. quoNeg := xNeg # yNeg;
  333. IF quoNeg & (remL # 0) THEN
  334. i := 0; c := 1;
  335. WHILE (i # quoL) & (c # 0) DO
  336. c := SHORT(c + xBuf[i + yL]); quo[i] := SHORT(c MOD B); c := SHORT(c DIV B);
  337. INC(i)
  338. END;
  339. IF c = 0 THEN
  340. WHILE i # quoL DO quo[i] := xBuf[i + yL]; INC(i) END
  341. ELSE quo[i] := c; INC(quoL)
  342. END
  343. ELSE
  344. i := 0;
  345. WHILE i # quoL DO quo[i] := xBuf[i + yL]; INC(i) END
  346. END
  347. END;
  348. IF doRem THEN
  349. remNeg := yNeg & (remL # 0);
  350. IF (xNeg # yNeg) & (remL # 0) THEN Subtract(y, xBuf, rem, yL, remL, remL)
  351. ELSE
  352. i := 0;
  353. WHILE i # remL DO rem[i] := xBuf[i]; INC(i) END
  354. END
  355. END
  356. END GetQuoRem;
  357. PROCEDURE BinPower (x: Integer; exp: INTEGER; y: Integer; xL: Index; OUT yL: Index);
  358. VAR zL: Index; b: INTEGER; z: Integer;
  359. BEGIN
  360. ASSERT(exp > 0, 20); ASSERT(xL # 0, 21);
  361. b := 1;
  362. WHILE 2 * b <= exp DO b := 2 * b END;
  363. y[0] := 1; yL := 1; NEW(z, LEN(y^));
  364. (* y^b * x^exp = const.) & (2 * b > exp) *)
  365. WHILE (exp # 0) OR (b # 1) DO
  366. IF exp >= b THEN
  367. exp := exp - b;
  368. IF xL <= yL THEN Multiply(x, y, z, xL, yL, zL) ELSE Multiply(y, x, z, yL, xL, zL) END
  369. ELSE b := b DIV 2; Multiply(y, y, z, yL, yL, zL)
  370. END;
  371. yL := zL;
  372. REPEAT DEC(zL); y[zL] := z[zL] UNTIL zL = 0
  373. END
  374. END BinPower;
  375. (* Data Format Support *)
  376. PROCEDURE New (nofDigits: Index): Integer;
  377. VAR x: Integer;
  378. BEGIN
  379. NEW(x, nofDigits + 2); RETURN x
  380. END New;
  381. PROCEDURE SetLength (x: Integer; len: Index; negative: BOOLEAN);
  382. VAR low, high: Digit;
  383. BEGIN
  384. ASSERT(len >= 0, 20); ASSERT(~negative OR (len # 0), 21);
  385. IF negative THEN len := -len END;
  386. low := SHORT(len MOD 10000H - 8000H); high := SHORT(len DIV 10000H);
  387. x[LEN(x^) - 1] := low; x[LEN(x^) - 2] := high
  388. END SetLength;
  389. PROCEDURE GetLength (x: Integer; OUT len: Index; OUT negative: BOOLEAN);
  390. VAR low, high: Digit;
  391. BEGIN
  392. low := x[LEN(x^) - 1]; high := x[LEN(x^) - 2];
  393. len := low + 8000H + high * 10000H;
  394. negative := len < 0; len := ABS(len)
  395. END GetLength;
  396. (* Exported Services *)
  397. PROCEDURE Long* (x: LONGINT): Integer;
  398. VAR i: Index; negative: BOOLEAN; int: Integer;
  399. BEGIN
  400. IF x # 0 THEN
  401. negative := x < 0; x := ABS(x);
  402. int := New(5); i := 0;
  403. REPEAT int[i] := SHORT(SHORT(x MOD B)); x := x DIV B; INC(i) UNTIL x = 0;
  404. SetLength(int, i, negative)
  405. ELSE int := zero
  406. END;
  407. RETURN int
  408. END Long;
  409. PROCEDURE Short* (x: Integer): LONGINT;
  410. VAR i: Index; res: LONGINT; negative: BOOLEAN;
  411. BEGIN
  412. res := 0; GetLength(x, i, negative);
  413. WHILE i # 0 DO DEC(i); res := res * B + x[i] END;
  414. IF negative THEN res := -res END;
  415. RETURN res
  416. END Short;
  417. PROCEDURE Entier* (x: REAL): Integer;
  418. VAR mL, yL, i: Index; mx: REAL; ex: INTEGER; neg: BOOLEAN; y, z: Integer;
  419. PROCEDURE Inc(m: Integer; VAR mL: Index);
  420. VAR i: Index;
  421. BEGIN
  422. i := 0;
  423. WHILE m[i] = B - 1 DO m[i] := 0; INC(i) END;
  424. INC(m[i]);
  425. IF i = mL THEN INC(mL); m[mL] := 0 END
  426. END Inc;
  427. PROCEDURE Double (m: Integer; VAR mL: Index);
  428. VAR i: Index; c: Digit;
  429. BEGIN
  430. i := 0; c := 0;
  431. WHILE i < mL DO
  432. c := SHORT(c + m[i] * 2); m[i] := SHORT(c MOD B); c := SHORT(c DIV B);
  433. INC(i)
  434. END;
  435. IF c # 0 THEN INC(mL); m[mL] := 0; m[i] := c END
  436. END Double;
  437. BEGIN
  438. IF (x >= 1) OR (x < 0) THEN
  439. neg := x < 0; x := ABS(x);
  440. mL := 0; buf6[0] := 0; mx := Math.Mantissa(x); ex := Math.Exponent(x);
  441. WHILE (mx # 0) & (ex > 0) DO (* mx * 2^ex + m * 2^ex = const. *)
  442. IF ENTIER(mx) = 1 THEN Inc(buf6, mL); mx := mx - 1
  443. ELSE ASSERT(ENTIER(mx) = 0, 100)
  444. END;
  445. Double(buf6, mL); mx := 2 * mx; DEC(ex)
  446. END;
  447. IF (ENTIER(mx) = 1) & (ex = 0) THEN Inc(buf6, mL); mx := mx - 1 END;
  448. IF ex > 0 THEN
  449. y := New(mL + SHORT(ENTIER(Math.Ln(2) * ex / Math.Ln(B)) + 1));
  450. z := New(SHORT(ENTIER(Math.Ln(2) * ex / Math.Ln(B)) + 1));
  451. BinPower(two, ex, z, 1, yL);
  452. IF mL <= yL THEN Multiply(buf6, z, y, mL, yL, yL) ELSE Multiply(z, buf6, y, yL, mL, yL) END
  453. ELSE
  454. y := New(mL + 1); yL := mL;
  455. i := 0;
  456. WHILE i # mL DO y[i] := buf6[i]; INC(i) END
  457. END;
  458. IF neg & (mx # 0) THEN Inc(y, yL) END;
  459. SetLength(y, yL, neg)
  460. ELSE y := zero
  461. END;
  462. RETURN y
  463. END Entier;
  464. PROCEDURE Float* (x: Integer): REAL;
  465. VAR i: Index; y: REAL; negative: BOOLEAN;
  466. BEGIN
  467. y := 0; GetLength(x, i, negative);
  468. WHILE i # 0 DO DEC(i); y := y * B + x[i] END;
  469. IF negative THEN y := -y END;
  470. RETURN y
  471. END Float;
  472. PROCEDURE Sign* (x: Integer): INTEGER;
  473. VAR len: Index; negative: BOOLEAN;
  474. BEGIN
  475. GetLength(x, len, negative);
  476. IF len = 0 THEN RETURN 0
  477. ELSIF negative THEN RETURN -1
  478. ELSE RETURN 1
  479. END
  480. END Sign;
  481. PROCEDURE Abs* (x: Integer): Integer;
  482. VAR len: Index; negative: BOOLEAN; y: Integer;
  483. BEGIN
  484. GetLength(x, len, negative);
  485. IF negative THEN
  486. y := New(len); SetLength(y, len, FALSE);
  487. REPEAT DEC(len); y[len] := x[len] UNTIL len = 0
  488. ELSE y := x
  489. END;
  490. RETURN y
  491. END Abs;
  492. PROCEDURE Digits10Of* (x: Integer): INTEGER;
  493. VAR i, n: Index; d: Digit; negative: BOOLEAN;
  494. BEGIN
  495. GetLength(x, n, negative);
  496. IF n # 0 THEN
  497. d := x[n - 1]; i := 0;
  498. REPEAT INC(i); d := SHORT(d DIV 10) UNTIL d = 0;
  499. n := DecPerDig * (n - 1) + i
  500. END;
  501. RETURN n
  502. END Digits10Of;
  503. PROCEDURE ThisDigit10* (x: Integer; exp10: INTEGER): CHAR;
  504. VAR i, n: Index; d: Digit; negative: BOOLEAN;
  505. BEGIN
  506. ASSERT(exp10 >= 0, 20);
  507. GetLength(x, n, negative); i := exp10 DIV DecPerDig;
  508. IF n > i THEN
  509. d := x[i]; i := exp10 MOD DecPerDig;
  510. WHILE i # 0 DO d := SHORT(d DIV 10); DEC(i) END;
  511. d := SHORT(d MOD 10)
  512. ELSE d := 0
  513. END;
  514. RETURN CHR(ORD("0") + d)
  515. END ThisDigit10;
  516. PROCEDURE Compare* (x, y: Integer): INTEGER;
  517. VAR xL, yL: Index; res: INTEGER; xNeg, yNeg: BOOLEAN;
  518. BEGIN
  519. GetLength(x, xL, xNeg); GetLength(y, yL, yNeg);
  520. IF xNeg = yNeg THEN
  521. IF (xL = yL) & (xL # 0) THEN
  522. DEC(xL);
  523. WHILE (xL # 0) & (x[xL] = y[xL]) DO DEC(xL) END;
  524. IF x[xL] = y[xL] THEN res := 0 ELSIF (x[xL] < y[xL]) = xNeg THEN res := 1 ELSE res := -1 END
  525. ELSE
  526. IF xL = yL THEN res := 0 ELSIF (xL < yL) = xNeg THEN res := 1 ELSE res := -1 END
  527. END
  528. ELSIF xNeg THEN res := -1
  529. ELSE res := 1
  530. END;
  531. RETURN res
  532. END Compare;
  533. PROCEDURE AddOp (x, y: Integer; subtract: BOOLEAN): Integer;
  534. VAR i, d, xL, yL, intL: Index; xNeg, yNeg: BOOLEAN; int: Integer;
  535. BEGIN
  536. GetLength(x, xL, xNeg); GetLength(y, yL, yNeg);
  537. IF yL = 0 THEN int := x
  538. ELSIF xL = 0 THEN
  539. IF subtract THEN
  540. int := New(yL); SetLength(int, yL, ~yNeg);
  541. REPEAT DEC(yL); int[yL] := y[yL] UNTIL yL = 0
  542. ELSE int := y
  543. END
  544. ELSIF (xNeg = yNeg) # subtract THEN
  545. int := New(MAX(xL, yL) + 1); Add(x, y, int, xL, yL, intL); SetLength(int, intL, xNeg)
  546. ELSE
  547. d := xL - yL;
  548. IF d # 0 THEN i := MAX(xL, yL) - 1
  549. ELSE
  550. i := xL;
  551. REPEAT DEC(i); d := x[i] - y[i] UNTIL (i = 0) OR (d # 0)
  552. END;
  553. IF d > 0 THEN
  554. int := New(i + 1); Subtract(x, y, int, xL, yL, intL); SetLength(int, intL, xNeg)
  555. ELSIF d < 0 THEN
  556. int := New(i + 1); Subtract(y, x, int, yL, xL, intL); SetLength(int, intL, yNeg # subtract)
  557. ELSE int := zero
  558. END
  559. END;
  560. RETURN int
  561. END AddOp;
  562. PROCEDURE Sum* (x, y: Integer): Integer;
  563. BEGIN
  564. RETURN AddOp(x, y, FALSE)
  565. END Sum;
  566. PROCEDURE Difference*(x, y: Integer): Integer;
  567. BEGIN
  568. RETURN AddOp(x, y, TRUE)
  569. END Difference;
  570. PROCEDURE Product* (x, y: Integer): Integer;
  571. VAR xL, yL, intL: Index; neg, xNeg, yNeg: BOOLEAN; int: Integer;
  572. BEGIN
  573. GetLength(x, xL, xNeg); GetLength(y, yL, yNeg); neg := xNeg # yNeg;
  574. IF xL > yL THEN int := x; x := y; y := int; intL := xL; xL := yL; yL := intL; xNeg := yNeg END;
  575. (* x.nofDigits <= y.nofDigits - yNeg no more valid! *)
  576. IF xL = 0 THEN int := zero
  577. ELSIF (xL = 1) & (x[0] = 1) THEN
  578. IF xNeg THEN
  579. int := New(yL); SetLength(int, yL, neg);
  580. REPEAT DEC(yL); int[yL] := y[yL] UNTIL yL = 0
  581. ELSE int := y
  582. END
  583. ELSE
  584. int := New(xL + yL); Multiply(x, y, int, xL, yL, intL); SetLength(int, intL, neg)
  585. END;
  586. RETURN int
  587. END Product;
  588. PROCEDURE Quotient* (x, y: Integer): Integer;
  589. VAR xL, yL, intL, remL: Index; xNeg, yNeg, qNeg, rNeg: BOOLEAN;
  590. int: Integer;
  591. BEGIN
  592. GetLength(x, xL, xNeg); GetLength(y, yL, yNeg);
  593. ASSERT(yL # 0, 20);
  594. IF xL < yL THEN int := zero
  595. ELSIF (yL = 1) & (y[0] = 1) THEN
  596. IF yNeg THEN
  597. int := New(xL); SetLength(int, xL, ~xNeg);
  598. REPEAT DEC(xL); int[xL] := x[xL] UNTIL xL = 0
  599. ELSE int := x
  600. END
  601. ELSE
  602. int := New(xL - yL + 2);
  603. GetQuoRem(x, y, xL, yL, xNeg, yNeg, int, NIL, intL, remL, qNeg, rNeg, TRUE, FALSE);
  604. SetLength(int, intL, qNeg)
  605. END;
  606. RETURN int
  607. END Quotient;
  608. PROCEDURE Remainder* (x, y: Integer): Integer;
  609. VAR xL, yL, intL, quoL: Index; xNeg, yNeg, qNeg, rNeg: BOOLEAN;
  610. int: Integer;
  611. BEGIN
  612. GetLength(x, xL, xNeg); GetLength(y, yL, yNeg);
  613. ASSERT(yL # 0, 20);
  614. IF xL < yL THEN int := x
  615. ELSIF (yL = 1) & (y[0] = 1) THEN int := zero
  616. ELSE
  617. int := New(yL);
  618. GetQuoRem(x, y, xL, yL, xNeg, yNeg, NIL, int, quoL, intL, qNeg, rNeg, FALSE, TRUE);
  619. SetLength(int, intL, rNeg)
  620. END;
  621. RETURN int
  622. END Remainder;
  623. PROCEDURE QuoRem* (x, y: Integer; OUT quo, rem: Integer);
  624. VAR xL, yL, quoL, remL: Index; xNeg, yNeg, qNeg, rNeg: BOOLEAN;
  625. BEGIN
  626. GetLength(x, xL, xNeg); GetLength(y, yL, yNeg);
  627. ASSERT(yL # 0, 20);
  628. IF xL < yL THEN quo := zero; rem := x
  629. ELSIF (yL = 1) & (y[0] = 1) THEN
  630. rem := zero;
  631. IF yNeg THEN
  632. quo := New(xL); SetLength(quo, xL, ~xNeg);
  633. REPEAT DEC(xL); quo[xL] := x[xL] UNTIL xL = 0
  634. ELSE quo := x
  635. END
  636. ELSE
  637. quo := New(xL - yL + 2); rem := New(yL);
  638. GetQuoRem(x, y, xL, yL, xNeg, yNeg, quo, rem, quoL, remL, qNeg, rNeg, TRUE, TRUE);
  639. SetLength(quo, quoL, qNeg); SetLength(rem, remL, rNeg)
  640. END
  641. END QuoRem;
  642. PROCEDURE GCD* (x, y: Integer): Integer;
  643. VAR xL, yL, i: Index; h: Digit; negative: BOOLEAN; xBuf, yBuf, int: Integer;
  644. BEGIN
  645. GetLength(x, xL, negative); GetLength(y, yL, negative);
  646. IF xL = 0 THEN int := y
  647. ELSIF yL = 0 THEN int := x
  648. ELSE
  649. IF xL >= yL THEN xBuf := CopyOf(x, xL + 1); yBuf := CopyOf(y, yL + 1)
  650. ELSE xBuf := CopyOf(y, yL + 1); yBuf := CopyOf(x, xL + 1); i := xL; xL := yL; yL := i
  651. END;
  652. WHILE yL # 0 DO
  653. DecomposeQuoRem(xBuf, yBuf, xL, yL);
  654. xL := yL;
  655. WHILE (xL # 0) & (xBuf[xL - 1] = 0) DO DEC(xL) END;
  656. i := yL;
  657. WHILE i # 0 DO DEC(i); h := xBuf[i]; xBuf[i] := yBuf[i]; yBuf[i] := h END;
  658. i := xL; xL := yL; yL := i
  659. END;
  660. int := New(xL); SetLength(int, xL, FALSE);
  661. WHILE xL # 0 DO DEC(xL); int[xL] := xBuf[xL] END
  662. END;
  663. RETURN int
  664. END GCD;
  665. PROCEDURE Power* (x: Integer; exp: INTEGER): Integer;
  666. VAR xL, intL: Index; negative: BOOLEAN; int: Integer;
  667. BEGIN
  668. ASSERT(exp >= 0, 20);
  669. GetLength(x, xL, negative);
  670. IF xL = 0 THEN int := zero
  671. ELSIF (xL = 1) & (x[0] = 1) THEN
  672. IF negative & ~ODD(exp) THEN
  673. int := New(xL); SetLength(int, xL, FALSE);
  674. REPEAT DEC(xL); int[xL] := x[xL] UNTIL xL = 0
  675. ELSE int := x
  676. END
  677. ELSIF exp = 0 THEN int := one
  678. ELSIF exp = 1 THEN int := x
  679. ELSE
  680. int := New(SHORT((xL - 1) * exp + ENTIER(Math.Ln(x[xL - 1] + 1) * exp / Math.Ln(B)) + 1));
  681. BinPower(x, exp, int, xL, intL); SetLength(int, intL, negative & ODD(exp))
  682. END;
  683. RETURN int
  684. END Power;
  685. (* Read from and Write to String and File *)
  686. PROCEDURE ConvertFromString* (IN s: ARRAY OF CHAR; OUT x: Integer);
  687. VAR i, j, k: INTEGER; dig, b: Digit; ch: CHAR; negative: BOOLEAN; new: Integer;
  688. BEGIN
  689. i := 0; ch := s[0];
  690. WHILE (ch # 0X) & (ch <= " ") DO INC(i); ch := s[i] END;
  691. negative := ch = "-";
  692. IF negative THEN INC(i); ch := s[i] END;
  693. IF ch = "+" THEN INC(i); ch := s[i] END;
  694. WHILE (ch # 0X) & (ch <= " ") DO INC(i); ch := s[i] END;
  695. ASSERT((ch >= "0") & (ch <= "9"), 20);
  696. WHILE ch = "0" DO INC(i); ch := s[i] END;
  697. IF (ch > "0") & (ch <= "9") THEN
  698. j := i;
  699. REPEAT INC(j); ch := s[j] UNTIL (ch < "0") OR (ch > "9");
  700. k := (j - i - 1) DIV DecPerDig + 2;
  701. new := New(k); SetLength(new, k - 1, negative);
  702. k := (j - i) MOD DecPerDig;
  703. IF k # 0 THEN
  704. b := 1; DEC(k);
  705. WHILE k # 0 DO DEC(k); b := SHORT(b * 10) END
  706. ELSE b := B DIV 10
  707. END;
  708. REPEAT
  709. dig := 0;
  710. WHILE b # 0 DO
  711. dig := SHORT(dig + b * (ORD(s[i]) - ORD("0"))); b := SHORT(b DIV 10);
  712. INC(i)
  713. END;
  714. new[(j - i) DIV DecPerDig] := dig; b := B DIV 10
  715. UNTIL i = j;
  716. x := new
  717. ELSE x := zero
  718. END
  719. END ConvertFromString;
  720. PROCEDURE ConvertToString* (x: Integer; OUT s: ARRAY OF CHAR);
  721. VAR j: Index; i: INTEGER; d, b: Digit; negative: BOOLEAN;
  722. BEGIN
  723. GetLength(x, j, negative);
  724. IF negative THEN s[0] := "-"; i := 1 ELSE i := 0 END;
  725. IF j # 0 THEN
  726. DEC(j); d := x[j]; b := B DIV 10;
  727. WHILE d DIV b = 0 DO b := SHORT(b DIV 10) END;
  728. REPEAT
  729. s[i] := CHR(d DIV b + ORD("0")); INC(i); d := SHORT(d MOD b); b := SHORT(b DIV 10)
  730. UNTIL b = 0;
  731. WHILE j # 0 DO
  732. DEC(j); d := x[j]; b := B DIV 10;
  733. REPEAT
  734. s[i] := CHR(d DIV b + ORD("0")); INC(i); d := SHORT(d MOD b); b := SHORT(b DIV 10)
  735. UNTIL b = 0
  736. END
  737. ELSE s[i] := "0"; INC(i)
  738. END;
  739. s[i] := 0X
  740. END ConvertToString;
  741. PROCEDURE Internalize* (r: Files.Reader; OUT x: Integer);
  742. VAR len: Index; n, version: INTEGER; negative: BOOLEAN;
  743. new: Integer; buf: ARRAY 4 OF BYTE;
  744. BEGIN
  745. r.ReadByte(buf[0]); version := buf[0];
  746. ASSERT((version = 0) OR (version >= 128), 20);
  747. IF version = 0 THEN
  748. r.ReadBytes(buf, 0, 4);
  749. len := (((buf[0] MOD 128) * 256 + buf[1] MOD 256) * 256
  750. + buf[2] MOD 256) * 256 + buf[3] MOD 256;
  751. new := New(len); SetLength(new, len, buf[0] < 0);
  752. WHILE len # 0 DO
  753. DEC(len);
  754. r.ReadBytes(buf, 0, 2); new[len] := SHORT((buf[0] MOD 256) * 256 + buf[1] MOD 256)
  755. END;
  756. x := new
  757. ELSE (* version >= 128 *)
  758. r.ReadByte(buf[1]); n := (buf[0] MOD 256) * 256 + buf[1] MOD 256 - 32768;
  759. r.ReadBytes(buf, 0, 2); DEC(n);
  760. len := (buf[0] MOD 256) * 256 + buf[1] MOD 256; negative := len < 0; len := ABS(len);
  761. new := New(len); SetLength(new, len, negative);
  762. WHILE n # len DO DEC(n); r.ReadBytes(buf, 0, 2) END;
  763. WHILE len # 0 DO
  764. DEC(len);
  765. r.ReadBytes(buf, 0, 2); new[len] := SHORT((buf[0] MOD 256) * 256 + buf[1] MOD 256)
  766. END;
  767. x := new
  768. END
  769. END Internalize;
  770. PROCEDURE Externalize* (w: Files.Writer; x: Integer);
  771. VAR len, l: Index; d: Digit; i: INTEGER; negative: BOOLEAN; buf: ARRAY 4 OF BYTE;
  772. PROCEDURE Byte(x: INTEGER): BYTE;
  773. BEGIN
  774. ASSERT((x >= MIN(BYTE)) & (x <= MAX(BYTE) - MIN(BYTE)), 20);
  775. IF x > MAX(BYTE) THEN RETURN SHORT(SHORT(x - 256)) ELSE RETURN SHORT(SHORT(x)) END
  776. END Byte;
  777. BEGIN
  778. GetLength(x, len, negative); l := len; i := 4;
  779. REPEAT DEC(i); buf[i] := Byte(l MOD 256); l := l DIV 256 UNTIL i = 0;
  780. IF negative THEN buf[0] := Byte(128 + buf[0] MOD 256) END;
  781. w.WriteByte(0); w.WriteBytes(buf, 0, 4);
  782. WHILE len # 0 DO
  783. DEC(len);
  784. d := x[len]; buf[0] := Byte(d DIV 256); buf[1] := Byte(d MOD 256); w.WriteBytes(buf, 0, 2)
  785. END
  786. END Externalize;
  787. BEGIN
  788. ASSERT(B <= BinBase, 20);
  789. zero := New(0); SetLength(zero, 0, FALSE);
  790. one := New(1); one[0] := 1; SetLength(one, 1, FALSE);
  791. two := New(1); two[0] := 2; SetLength(two, 1, FALSE);
  792. NEW(buf6, 6)
  793. END Integers.