123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700701702703704705706707708709710711712713714715716717718719720721722723724725726727728729730731732733734735736737738739740741742743744745746747748749750751752753754755756757758759760761762763764765766767768769770771772773774775776777778779780781782783784785786787788789790791792793794795796797798799800801802803804805806807808809810811812813814815816817818819820821822823824825826827828829830831832833834835836837838839840841842843844845846847848 |
- MODULE Integers;
- (* THIS IS TEXT COPY OF BlackBox 1.6-rc6 System/Mod/Integers.odc *)
- (* DO NOT EDIT *)
- IMPORT Files, Math;
- CONST
- B = 10000; DecPerDig = 4; BinBase = 16 * 1024;
- KaratsubaBreak = 41;
- TYPE
- Index = INTEGER;
- Digit = SHORTINT;
- DoubleDigit = INTEGER;
- IntegerDesc = ARRAY OF Digit; (* to hide internal structure from interface *)
- Integer* = POINTER TO IntegerDesc;
- Buffer = RECORD
- digit: Integer;
- beg, len: Index
- END;
- VAR zero, one, two, buf6: Integer;
- PROCEDURE CopyOf (x: Integer; len: Index): Integer;
- VAR buf: Integer;
- BEGIN
- ASSERT(len > 0, 20);
- NEW(buf, len);
- REPEAT DEC(len); buf[len] := x[len] UNTIL len = 0;
- RETURN buf
- END CopyOf;
- (* Operations on Digits *)
- PROCEDURE Add (x, y, sum: Integer; xL, yL: Index; OUT sumL: Index);
- VAR i, l: Index; c: Digit;
- BEGIN
- l := MIN(xL, yL);
- i := 0; c := 0;
- WHILE i < l DO c := SHORT(c DIV B + x[i] + y[i]); sum[i] := SHORT(c MOD B); INC(i) END;
- WHILE i < xL DO c := SHORT(c DIV B + x[i]); sum[i] := SHORT(c MOD B); INC(i) END;
- WHILE i < yL DO c := SHORT(c DIV B + y[i]); sum[i] := SHORT(c MOD B); INC(i) END;
- IF c >= B THEN sum[i] := SHORT(c DIV B); INC(i) END;
- sumL := i
- END Add;
- PROCEDURE Subtract (x, y, dif: Integer; xL, yL: Index; OUT difL: Index);
- VAR i: Index; c, d: Digit;
- BEGIN
- ASSERT(xL >= yL, 20);
- i := 0; difL := 0; c := 0;
- WHILE i < yL DO
- c := SHORT(c DIV B + x[i] - y[i]); d := SHORT(c MOD B);
- IF d # 0 THEN
- WHILE difL # i DO dif[difL] := 0; INC(difL) END;
- dif[i] := d; INC(difL)
- END;
- INC(i)
- END;
- WHILE i < xL DO
- c := SHORT(c DIV B + x[i]); d := SHORT(c MOD B);
- IF d # 0 THEN
- WHILE difL # i DO dif[difL] := 0; INC(difL) END;
- dif[i] := d; INC(difL)
- END;
- INC(i)
- END;
- ASSERT(c DIV B = 0, 100)
- END Subtract;
- PROCEDURE OneDigitMult (a, b: Buffer; VAR c: Buffer);
- VAR i: Index; carry, factor: DoubleDigit;
- BEGIN
- ASSERT(a.len = 1, 20);
- factor := a.digit[a.beg]; i := 0; carry := 0;
- WHILE i # b.len DO
- carry := carry DIV B + factor * b.digit[b.beg + i]; c.digit[c.beg + i] := SHORT(carry MOD B);
- INC(i)
- END;
- IF carry >= B THEN c.digit[c.beg + i] := SHORT(carry DIV B); INC(i) END;
- c.len := i
- END OneDigitMult;
- PROCEDURE SimpleMult (a, b: Buffer; VAR c: Buffer);
- VAR i, j, k: Index; c0, c1: DoubleDigit;
- BEGIN
- ASSERT(a.len <= b.len, 20);
- c.len := a.len + b.len - 1;
- i := 0; c0 := 0; c1 := 0;
- REPEAT
- IF i < b.len THEN
- IF i < a.len THEN j := i; k := 0 ELSE j := a.len - 1; k := i - a.len + 1 END;
- REPEAT
- c0 := c0 + a.digit[a.beg + j] * b.digit[b.beg + k];
- IF c0 > MAX(DoubleDigit) - BinBase * (B - 1) THEN
- c1 := c1 + c0 DIV BinBase; c0 := c0 MOD BinBase
- END;
- DEC(j); INC(k)
- UNTIL j < 0
- ELSE
- j := a.len - 1; k := i - a.len + 1;
- REPEAT
- c0 := c0 + a.digit[a.beg + j] * b.digit[b.beg + k];
- IF c0 > MAX(DoubleDigit) - BinBase * (B - 1) THEN
- c1 := c1 + c0 DIV BinBase; c0 := c0 MOD BinBase
- END;
- DEC(j); INC(k)
- UNTIL k = b.len
- END;
- IF c1 = 0 THEN c.digit[c.beg + i] := SHORT(c0 MOD B); c0 := c0 DIV B
- ELSE
- c0 := c0 + BinBase * (c1 MOD B);
- c.digit[c.beg + i] := SHORT(c0 MOD B); c0 := c0 DIV B; c1 := c1 DIV B
- END;
- INC(i)
- UNTIL i = c.len;
- IF c0 # 0 THEN c.digit[c.beg + c.len] := SHORT(c0); INC(c.len) END
- END SimpleMult;
- PROCEDURE AddBuf (a, b: Buffer; VAR c: Buffer); (* c := a + b *)
- VAR i: Index; carry: Digit;
- BEGIN
- ASSERT(a.len <= b.len, 20);
- i := 0; carry := 0;
- WHILE i # a.len DO
- carry := SHORT(carry DIV B + a.digit[a.beg + i] + b.digit[b.beg + i]);
- c.digit[c.beg + i] := SHORT(carry MOD B);
- INC(i)
- END;
- WHILE (i # b.len) & (carry >= B) DO
- carry := SHORT(carry DIV B + b.digit[b.beg + i]); c.digit[c.beg + i] := SHORT(carry MOD B);
- INC(i)
- END;
- IF carry >= B THEN c.digit[c.beg + i] := SHORT(carry DIV B); INC(i)
- ELSE
- WHILE i # b.len DO c.digit[c.beg + i] := b.digit[b.beg + i]; INC(i) END
- END;
- c.len := i
- END AddBuf;
- PROCEDURE AddToBuf (VAR a: Buffer; b: Buffer; shift: Index); (* a := a + b * B^shift *)
- VAR i, n: Index; carry: Digit;
- BEGIN
- b.beg := b.beg - shift; b.len := b.len + shift; i := shift; n := MIN(a.len, b.len); carry := 0;
- WHILE i # n DO
- carry := SHORT(carry DIV B + a.digit[a.beg + i] + b.digit[b.beg + i]);
- a.digit[a.beg + i] := SHORT(carry MOD B);
- INC(i)
- END;
- IF i # a.len THEN
- WHILE (i # a.len) & (carry >= B) DO
- carry := SHORT(carry DIV B + a.digit[a.beg + i]); a.digit[a.beg + i] := SHORT(carry MOD B);
- INC(i)
- END;
- IF carry >= B THEN a.digit[a.beg + i] := SHORT(carry DIV B); INC(i) END
- ELSE
- WHILE (i # b.len) & (carry >= B) DO
- carry := SHORT(carry DIV B + b.digit[b.beg + i]); a.digit[a.beg + i] := SHORT(carry MOD B);
- INC(i)
- END;
- IF carry >= B THEN a.digit[a.beg + i] := SHORT(carry DIV B); INC(i)
- ELSE
- WHILE i # b.len DO a.digit[a.beg + i] := b.digit[b.beg + i]; INC(i) END
- END
- END;
- a.len := MAX(i, a.len)
- END AddToBuf;
- PROCEDURE SubtractFromBuf (VAR a: Buffer; b, c: Buffer); (* a := a - b - c *)
- VAR i: Index; carry: Digit;
- BEGIN
- ASSERT(b.len <= c.len, 20);
- i := 0; carry := 0;
- WHILE i # b.len DO
- carry := SHORT(carry DIV B + a.digit[a.beg + i] - b.digit[b.beg + i] - c.digit[c.beg + i]);
- a.digit[a.beg + i] := SHORT(carry MOD B);
- INC(i)
- END;
- WHILE i # c.len DO
- carry := SHORT(carry DIV B + a.digit[a.beg + i] - c.digit[c.beg + i]);
- a.digit[a.beg + i] := SHORT(carry MOD B);
- INC(i)
- END;
- WHILE carry < 0 DO
- carry := SHORT(carry DIV B + a.digit[a.beg + i]); a.digit[a.beg + i] := SHORT(carry MOD B);
- INC(i)
- END;
- ASSERT(i <= a.len, 100);
- WHILE (a.len # 0) & (a.digit[a.beg + a.len - 1] = 0) DO DEC(a.len) END
- END SubtractFromBuf;
- PROCEDURE KStep (a, b: Buffer; VAR c: Buffer; stack: Buffer);
- VAR n2, i: Index; a0, a1, b0, b1, c0, c1, h: Buffer;
- BEGIN
- ASSERT(a.len <= b.len, 20);
- IF a.len = 0 THEN c.len := 0
- ELSIF a.len = 1 THEN OneDigitMult(a, b, c)
- ELSIF a.len <= KaratsubaBreak THEN SimpleMult(a, b, c)
- ELSE
- n2 := b.len DIV 2;
- c0.digit := c.digit; c0.beg := c.beg; c1.digit := c.digit; c1.beg := c.beg + 2 * n2;
- a0.digit := a.digit; a0.beg := a.beg; a0.len := MIN(a.len, n2);
- a1.digit := a.digit; a1.beg := a.beg + n2; a1.len := MAX(0, a.len - n2);
- WHILE (a0.len # 0) & (a0.digit[a0.beg + a0.len - 1] = 0) DO DEC(a0.len) END;
- b0.digit := b.digit; b0.beg := b.beg; b0.len := MIN(b.len, n2);
- b1.digit := b.digit; b1.beg := b.beg + n2; b1.len := MAX(0, b.len - n2);
- WHILE (b0.len # 0) & (b0.digit[b0.beg + b0.len - 1] = 0) DO DEC(b0.len) END;
- IF (a0.len # 0) OR (b0.len # 0) THEN
- IF a0.len <= a1.len THEN AddBuf(a0, a1, c1) ELSE AddBuf(a1, a0, c1) END;
- IF b0.len <= b1.len THEN AddBuf(b0, b1, c0) ELSE AddBuf(b1, b0, c0) END;
- h.digit := stack.digit; h.beg := stack.beg; stack.beg := stack.beg + c0.len + c1.len;
- IF c0.len <= c1.len THEN KStep(c0, c1, h, stack) ELSE KStep(c1, c0, h, stack) END;
- IF a0.len <= b0.len THEN KStep(a0, b0, c0, stack) ELSE KStep(b0, a0, c0, stack) END;
- KStep(a1, b1, c1, stack);
- IF c0.len <= c1.len THEN SubtractFromBuf(h, c0, c1) ELSE SubtractFromBuf(h, c1, c0) END;
- IF c1.len # 0 THEN
- i := c0.beg + c0.len;
- WHILE i < c1.beg DO c.digit[i] := 0; INC(i) END;
- c.len := c1.beg + c1.len - c.beg
- ELSE
- WHILE c0.len < n2 DO c0.digit[c0.beg + c0.len] := 0; INC(c0.len) END;
- c.len := c0.len
- END;
- ASSERT(h.len # 0, 100);
- AddToBuf(c, h, n2)
- ELSE
- KStep(a1, b1, c1, stack); c.len := c1.beg + c1.len - c.beg;
- i := c.beg;
- WHILE i # c1.beg DO c.digit[i] := 0; INC(i) END
- END
- END
- END KStep;
- PROCEDURE Karatsuba (x, y, pro:Integer; xL, yL: Index; OUT proL: Index);
- VAR a, b, c, stack: Buffer;
- BEGIN
- ASSERT(xL <= yL, 20);
- a.digit := x; a.beg := 0; a.len := xL; b.digit := y; b.beg := 0; b.len := yL;
- c.digit := pro; c.beg := 0;
- NEW(stack.digit, 2 * b.len); stack.beg := 0;
- KStep(a, b, c, stack);
- proL := c.len
- END Karatsuba;
-
- PROCEDURE Multiply (x, y, pro: Integer; xL, yL: Index; OUT proL: Index);
- VAR i, j, k: Index; c0, c1: DoubleDigit;
- BEGIN
- ASSERT(xL <= yL, 20);
- IF xL > KaratsubaBreak THEN Karatsuba(x, y, pro, xL, yL, proL)
- ELSIF xL = 1 THEN
- proL := 0; c1 := x[0]; c0 := 0;
- WHILE proL < yL DO
- c0 := c1 * y[proL] + c0; pro[proL] := SHORT(c0 MOD B);
- c0 := c0 DIV B ; INC(proL)
- END;
- IF c0 # 0 THEN pro[proL] := SHORT(c0); INC(proL) END
- ELSE
- proL := xL + yL - 1;
- i := 0; c0 := 0; c1 := 0;
- REPEAT
- IF i < yL THEN
- IF i < xL THEN j := i; k := 0 ELSE j := xL - 1; k := i - xL + 1 END;
- REPEAT
- c0 := c0 + x[j] * y[k];
- IF c0 > MAX(DoubleDigit) - BinBase * (B - 1) THEN
- c1 := c1 + c0 DIV BinBase; c0 := c0 MOD BinBase
- END;
- DEC(j); INC(k)
- UNTIL j < 0
- ELSE
- j := xL - 1; k := i - xL + 1;
- REPEAT
- c0 := c0 + x[j] * y[k];
- IF c0 > MAX(DoubleDigit) - BinBase * (B - 1) THEN
- c1 := c1 + c0 DIV BinBase; c0 := c0 MOD BinBase
- END;
- DEC(j); INC(k)
- UNTIL k = yL
- END;
- IF c1 = 0 THEN pro[i] := SHORT(c0 MOD B); c0 := c0 DIV B
- ELSE c0 := c0 + BinBase * (c1 MOD B); pro[i] := SHORT(c0 MOD B);
- c0 := c0 DIV B; c1 := c1 DIV B
- END;
- INC(i)
- UNTIL i = proL;
- IF c0 # 0 THEN pro[proL] := SHORT(c0); INC(proL) END
- END
- END Multiply;
- PROCEDURE DecomposeQuoRem (x, y: Integer; xL, yL: Index);
- VAR ix, iy, j: Index; d, q, h, yLead, ySecond: DoubleDigit; yBuf: Integer;
- BEGIN
- ASSERT((yL # 0) & (y[yL - 1] # 0), 20);
- IF yL = 1 THEN
- j := xL - 1; h := 0; d := y[0];
- WHILE j >= 0 DO h := x[j] + h * B; x[j + 1] := SHORT(h DIV d); h := h MOD d; DEC(j) END;
- x[0] := SHORT(h)
- ELSIF xL >= yL THEN
- x[xL] := 0; d := (B DIV 2 - 1) DIV y[yL - 1] + 1; yBuf := CopyOf(y, yL);
- IF d # 1 THEN
- j := 0; h := 0;
- WHILE j < xL DO h := d * x[j] + h DIV B; x[j] := SHORT(h MOD B); INC(j) END;
- x[xL] := SHORT(h DIV B);
- j := 0; h := 0;
- WHILE j < yL DO h := d * yBuf[j] + h DIV B; yBuf[j] := SHORT(h MOD B); INC(j) END;
- ASSERT(h DIV B = 0, 100)
- END;
- yLead := yBuf[yL - 1]; ySecond := yBuf[yL - 2]; j := xL;
- WHILE j >= yL DO
- IF x[j] # yLead THEN q := (x[j] * B + x[j - 1]) DIV yLead ELSE q := B - 1 END;
- WHILE ySecond * q > (x[j] * B + x[j - 1] - yLead * q) * B + x[j - 2] DO
- DEC(q)
- END;
- ix := j - yL; iy := 0; h := 0;
- WHILE iy < yL DO
- h := x[ix] - q * yBuf[iy] + h DIV B; x[ix] := SHORT(h MOD B); INC(ix); INC(iy)
- END;
- IF (-x[j]) # (h DIV B) THEN
- ix := j - yL; iy := 0; h := 0;
- WHILE iy < yL DO
- h := h DIV B + x[ix] + yBuf[iy]; x[ix] := SHORT(h MOD B); INC(ix); INC(iy)
- END;
- x[j] := SHORT(q - 1)
- ELSE x[j] := SHORT(q)
- END;
- DEC(j)
- END;
- IF d # 1 THEN
- j := yL; h := 0;
- WHILE j # 0 DO DEC(j); h := h + x[j]; x[j] := SHORT(h DIV d); h := (h MOD d) * B END
- END
- END
- END DecomposeQuoRem;
- PROCEDURE GetQuoRem (x, y: Integer; xL, yL: Index; xNeg, yNeg: BOOLEAN;
- quo, rem: Integer; OUT quoL, remL: Index; OUT quoNeg, remNeg: BOOLEAN;
- doQuo, doRem: BOOLEAN);
- VAR i: Index; c: Digit; xBuf: Integer;
- BEGIN
- ASSERT(xL >= yL, 20);
- xBuf := CopyOf(x, xL + 1);
- DecomposeQuoRem(xBuf, y, xL, yL);
- i := xL;
- WHILE (i >= yL) & (xBuf[i] = 0) DO DEC(i) END;
- quoL := i - yL + 1;
- i := yL - 1;
- WHILE (i >= 0) & (xBuf[i] = 0) DO DEC(i) END;
- remL := i + 1;
- IF doQuo THEN
- quoNeg := xNeg # yNeg;
- IF quoNeg & (remL # 0) THEN
- i := 0; c := 1;
- WHILE (i # quoL) & (c # 0) DO
- c := SHORT(c + xBuf[i + yL]); quo[i] := SHORT(c MOD B); c := SHORT(c DIV B);
- INC(i)
- END;
- IF c = 0 THEN
- WHILE i # quoL DO quo[i] := xBuf[i + yL]; INC(i) END
- ELSE quo[i] := c; INC(quoL)
- END
- ELSE
- i := 0;
- WHILE i # quoL DO quo[i] := xBuf[i + yL]; INC(i) END
- END
- END;
- IF doRem THEN
- remNeg := yNeg & (remL # 0);
- IF (xNeg # yNeg) & (remL # 0) THEN Subtract(y, xBuf, rem, yL, remL, remL)
- ELSE
- i := 0;
- WHILE i # remL DO rem[i] := xBuf[i]; INC(i) END
- END
- END
- END GetQuoRem;
- PROCEDURE BinPower (x: Integer; exp: INTEGER; y: Integer; xL: Index; OUT yL: Index);
- VAR zL: Index; b: INTEGER; z: Integer;
- BEGIN
- ASSERT(exp > 0, 20); ASSERT(xL # 0, 21);
- b := 1;
- WHILE 2 * b <= exp DO b := 2 * b END;
- y[0] := 1; yL := 1; NEW(z, LEN(y^));
- (* y^b * x^exp = const.) & (2 * b > exp) *)
- WHILE (exp # 0) OR (b # 1) DO
- IF exp >= b THEN
- exp := exp - b;
- IF xL <= yL THEN Multiply(x, y, z, xL, yL, zL) ELSE Multiply(y, x, z, yL, xL, zL) END
- ELSE b := b DIV 2; Multiply(y, y, z, yL, yL, zL)
- END;
- yL := zL;
- REPEAT DEC(zL); y[zL] := z[zL] UNTIL zL = 0
- END
- END BinPower;
- (* Data Format Support *)
- PROCEDURE New (nofDigits: Index): Integer;
- VAR x: Integer;
- BEGIN
- NEW(x, nofDigits + 2); RETURN x
- END New;
- PROCEDURE SetLength (x: Integer; len: Index; negative: BOOLEAN);
- VAR low, high: Digit;
- BEGIN
- ASSERT(len >= 0, 20); ASSERT(~negative OR (len # 0), 21);
- IF negative THEN len := -len END;
- low := SHORT(len MOD 10000H - 8000H); high := SHORT(len DIV 10000H);
- x[LEN(x^) - 1] := low; x[LEN(x^) - 2] := high
- END SetLength;
- PROCEDURE GetLength (x: Integer; OUT len: Index; OUT negative: BOOLEAN);
- VAR low, high: Digit;
- BEGIN
- low := x[LEN(x^) - 1]; high := x[LEN(x^) - 2];
- len := low + 8000H + high * 10000H;
- negative := len < 0; len := ABS(len)
- END GetLength;
- (* Exported Services *)
- PROCEDURE Long* (x: LONGINT): Integer;
- VAR i: Index; negative: BOOLEAN; int: Integer;
- BEGIN
- IF x # 0 THEN
- negative := x < 0; x := ABS(x);
- int := New(5); i := 0;
- REPEAT int[i] := SHORT(SHORT(x MOD B)); x := x DIV B; INC(i) UNTIL x = 0;
- SetLength(int, i, negative)
- ELSE int := zero
- END;
- RETURN int
- END Long;
- PROCEDURE Short* (x: Integer): LONGINT;
- VAR i: Index; res: LONGINT; negative: BOOLEAN;
- BEGIN
- res := 0; GetLength(x, i, negative);
- WHILE i # 0 DO DEC(i); res := res * B + x[i] END;
- IF negative THEN res := -res END;
- RETURN res
- END Short;
- PROCEDURE Entier* (x: REAL): Integer;
- VAR mL, yL, i: Index; mx: REAL; ex: INTEGER; neg: BOOLEAN; y, z: Integer;
- PROCEDURE Inc(m: Integer; VAR mL: Index);
- VAR i: Index;
- BEGIN
- i := 0;
- WHILE m[i] = B - 1 DO m[i] := 0; INC(i) END;
- INC(m[i]);
- IF i = mL THEN INC(mL); m[mL] := 0 END
- END Inc;
- PROCEDURE Double (m: Integer; VAR mL: Index);
- VAR i: Index; c: Digit;
- BEGIN
- i := 0; c := 0;
- WHILE i < mL DO
- c := SHORT(c + m[i] * 2); m[i] := SHORT(c MOD B); c := SHORT(c DIV B);
- INC(i)
- END;
- IF c # 0 THEN INC(mL); m[mL] := 0; m[i] := c END
- END Double;
- BEGIN
- IF (x >= 1) OR (x < 0) THEN
- neg := x < 0; x := ABS(x);
- mL := 0; buf6[0] := 0; mx := Math.Mantissa(x); ex := Math.Exponent(x);
- WHILE (mx # 0) & (ex > 0) DO (* mx * 2^ex + m * 2^ex = const. *)
- IF ENTIER(mx) = 1 THEN Inc(buf6, mL); mx := mx - 1
- ELSE ASSERT(ENTIER(mx) = 0, 100)
- END;
- Double(buf6, mL); mx := 2 * mx; DEC(ex)
- END;
- IF (ENTIER(mx) = 1) & (ex = 0) THEN Inc(buf6, mL); mx := mx - 1 END;
- IF ex > 0 THEN
- y := New(mL + SHORT(ENTIER(Math.Ln(2) * ex / Math.Ln(B)) + 1));
- z := New(SHORT(ENTIER(Math.Ln(2) * ex / Math.Ln(B)) + 1));
- BinPower(two, ex, z, 1, yL);
- IF mL <= yL THEN Multiply(buf6, z, y, mL, yL, yL) ELSE Multiply(z, buf6, y, yL, mL, yL) END
- ELSE
- y := New(mL + 1); yL := mL;
- i := 0;
- WHILE i # mL DO y[i] := buf6[i]; INC(i) END
- END;
- IF neg & (mx # 0) THEN Inc(y, yL) END;
- SetLength(y, yL, neg)
- ELSE y := zero
- END;
- RETURN y
- END Entier;
- PROCEDURE Float* (x: Integer): REAL;
- VAR i: Index; y: REAL; negative: BOOLEAN;
- BEGIN
- y := 0; GetLength(x, i, negative);
- WHILE i # 0 DO DEC(i); y := y * B + x[i] END;
- IF negative THEN y := -y END;
- RETURN y
- END Float;
- PROCEDURE Sign* (x: Integer): INTEGER;
- VAR len: Index; negative: BOOLEAN;
- BEGIN
- GetLength(x, len, negative);
- IF len = 0 THEN RETURN 0
- ELSIF negative THEN RETURN -1
- ELSE RETURN 1
- END
- END Sign;
- PROCEDURE Abs* (x: Integer): Integer;
- VAR len: Index; negative: BOOLEAN; y: Integer;
- BEGIN
- GetLength(x, len, negative);
- IF negative THEN
- y := New(len); SetLength(y, len, FALSE);
- REPEAT DEC(len); y[len] := x[len] UNTIL len = 0
- ELSE y := x
- END;
- RETURN y
- END Abs;
- PROCEDURE Digits10Of* (x: Integer): INTEGER;
- VAR i, n: Index; d: Digit; negative: BOOLEAN;
- BEGIN
- GetLength(x, n, negative);
- IF n # 0 THEN
- d := x[n - 1]; i := 0;
- REPEAT INC(i); d := SHORT(d DIV 10) UNTIL d = 0;
- n := DecPerDig * (n - 1) + i
- END;
- RETURN n
- END Digits10Of;
- PROCEDURE ThisDigit10* (x: Integer; exp10: INTEGER): CHAR;
- VAR i, n: Index; d: Digit; negative: BOOLEAN;
- BEGIN
- ASSERT(exp10 >= 0, 20);
- GetLength(x, n, negative); i := exp10 DIV DecPerDig;
- IF n > i THEN
- d := x[i]; i := exp10 MOD DecPerDig;
- WHILE i # 0 DO d := SHORT(d DIV 10); DEC(i) END;
- d := SHORT(d MOD 10)
- ELSE d := 0
- END;
- RETURN CHR(ORD("0") + d)
- END ThisDigit10;
- PROCEDURE Compare* (x, y: Integer): INTEGER;
- VAR xL, yL: Index; res: INTEGER; xNeg, yNeg: BOOLEAN;
- BEGIN
- GetLength(x, xL, xNeg); GetLength(y, yL, yNeg);
- IF xNeg = yNeg THEN
- IF (xL = yL) & (xL # 0) THEN
- DEC(xL);
- WHILE (xL # 0) & (x[xL] = y[xL]) DO DEC(xL) END;
- IF x[xL] = y[xL] THEN res := 0 ELSIF (x[xL] < y[xL]) = xNeg THEN res := 1 ELSE res := -1 END
- ELSE
- IF xL = yL THEN res := 0 ELSIF (xL < yL) = xNeg THEN res := 1 ELSE res := -1 END
- END
- ELSIF xNeg THEN res := -1
- ELSE res := 1
- END;
- RETURN res
- END Compare;
- PROCEDURE AddOp (x, y: Integer; subtract: BOOLEAN): Integer;
- VAR i, d, xL, yL, intL: Index; xNeg, yNeg: BOOLEAN; int: Integer;
- BEGIN
- GetLength(x, xL, xNeg); GetLength(y, yL, yNeg);
- IF yL = 0 THEN int := x
- ELSIF xL = 0 THEN
- IF subtract THEN
- int := New(yL); SetLength(int, yL, ~yNeg);
- REPEAT DEC(yL); int[yL] := y[yL] UNTIL yL = 0
- ELSE int := y
- END
- ELSIF (xNeg = yNeg) # subtract THEN
- int := New(MAX(xL, yL) + 1); Add(x, y, int, xL, yL, intL); SetLength(int, intL, xNeg)
- ELSE
- d := xL - yL;
- IF d # 0 THEN i := MAX(xL, yL) - 1
- ELSE
- i := xL;
- REPEAT DEC(i); d := x[i] - y[i] UNTIL (i = 0) OR (d # 0)
- END;
- IF d > 0 THEN
- int := New(i + 1); Subtract(x, y, int, xL, yL, intL); SetLength(int, intL, xNeg)
- ELSIF d < 0 THEN
- int := New(i + 1); Subtract(y, x, int, yL, xL, intL); SetLength(int, intL, yNeg # subtract)
- ELSE int := zero
- END
- END;
- RETURN int
- END AddOp;
- PROCEDURE Sum* (x, y: Integer): Integer;
- BEGIN
- RETURN AddOp(x, y, FALSE)
- END Sum;
- PROCEDURE Difference*(x, y: Integer): Integer;
- BEGIN
- RETURN AddOp(x, y, TRUE)
- END Difference;
- PROCEDURE Product* (x, y: Integer): Integer;
- VAR xL, yL, intL: Index; neg, xNeg, yNeg: BOOLEAN; int: Integer;
- BEGIN
- GetLength(x, xL, xNeg); GetLength(y, yL, yNeg); neg := xNeg # yNeg;
- IF xL > yL THEN int := x; x := y; y := int; intL := xL; xL := yL; yL := intL; xNeg := yNeg END;
- (* x.nofDigits <= y.nofDigits - yNeg no more valid! *)
- IF xL = 0 THEN int := zero
- ELSIF (xL = 1) & (x[0] = 1) THEN
- IF xNeg THEN
- int := New(yL); SetLength(int, yL, neg);
- REPEAT DEC(yL); int[yL] := y[yL] UNTIL yL = 0
- ELSE int := y
- END
- ELSE
- int := New(xL + yL); Multiply(x, y, int, xL, yL, intL); SetLength(int, intL, neg)
- END;
- RETURN int
- END Product;
- PROCEDURE Quotient* (x, y: Integer): Integer;
- VAR xL, yL, intL, remL: Index; xNeg, yNeg, qNeg, rNeg: BOOLEAN;
- int: Integer;
- BEGIN
- GetLength(x, xL, xNeg); GetLength(y, yL, yNeg);
- ASSERT(yL # 0, 20);
- IF xL < yL THEN int := zero
- ELSIF (yL = 1) & (y[0] = 1) THEN
- IF yNeg THEN
- int := New(xL); SetLength(int, xL, ~xNeg);
- REPEAT DEC(xL); int[xL] := x[xL] UNTIL xL = 0
- ELSE int := x
- END
- ELSE
- int := New(xL - yL + 2);
- GetQuoRem(x, y, xL, yL, xNeg, yNeg, int, NIL, intL, remL, qNeg, rNeg, TRUE, FALSE);
- SetLength(int, intL, qNeg)
- END;
- RETURN int
- END Quotient;
- PROCEDURE Remainder* (x, y: Integer): Integer;
- VAR xL, yL, intL, quoL: Index; xNeg, yNeg, qNeg, rNeg: BOOLEAN;
- int: Integer;
- BEGIN
- GetLength(x, xL, xNeg); GetLength(y, yL, yNeg);
- ASSERT(yL # 0, 20);
- IF xL < yL THEN int := x
- ELSIF (yL = 1) & (y[0] = 1) THEN int := zero
- ELSE
- int := New(yL);
- GetQuoRem(x, y, xL, yL, xNeg, yNeg, NIL, int, quoL, intL, qNeg, rNeg, FALSE, TRUE);
- SetLength(int, intL, rNeg)
- END;
- RETURN int
- END Remainder;
- PROCEDURE QuoRem* (x, y: Integer; OUT quo, rem: Integer);
- VAR xL, yL, quoL, remL: Index; xNeg, yNeg, qNeg, rNeg: BOOLEAN;
- BEGIN
- GetLength(x, xL, xNeg); GetLength(y, yL, yNeg);
- ASSERT(yL # 0, 20);
- IF xL < yL THEN quo := zero; rem := x
- ELSIF (yL = 1) & (y[0] = 1) THEN
- rem := zero;
- IF yNeg THEN
- quo := New(xL); SetLength(quo, xL, ~xNeg);
- REPEAT DEC(xL); quo[xL] := x[xL] UNTIL xL = 0
- ELSE quo := x
- END
- ELSE
- quo := New(xL - yL + 2); rem := New(yL);
- GetQuoRem(x, y, xL, yL, xNeg, yNeg, quo, rem, quoL, remL, qNeg, rNeg, TRUE, TRUE);
- SetLength(quo, quoL, qNeg); SetLength(rem, remL, rNeg)
- END
- END QuoRem;
- PROCEDURE GCD* (x, y: Integer): Integer;
- VAR xL, yL, i: Index; h: Digit; negative: BOOLEAN; xBuf, yBuf, int: Integer;
- BEGIN
- GetLength(x, xL, negative); GetLength(y, yL, negative);
- IF xL = 0 THEN int := y
- ELSIF yL = 0 THEN int := x
- ELSE
- IF xL >= yL THEN xBuf := CopyOf(x, xL + 1); yBuf := CopyOf(y, yL + 1)
- ELSE xBuf := CopyOf(y, yL + 1); yBuf := CopyOf(x, xL + 1); i := xL; xL := yL; yL := i
- END;
- WHILE yL # 0 DO
- DecomposeQuoRem(xBuf, yBuf, xL, yL);
- xL := yL;
- WHILE (xL # 0) & (xBuf[xL - 1] = 0) DO DEC(xL) END;
- i := yL;
- WHILE i # 0 DO DEC(i); h := xBuf[i]; xBuf[i] := yBuf[i]; yBuf[i] := h END;
- i := xL; xL := yL; yL := i
- END;
- int := New(xL); SetLength(int, xL, FALSE);
- WHILE xL # 0 DO DEC(xL); int[xL] := xBuf[xL] END
- END;
- RETURN int
- END GCD;
- PROCEDURE Power* (x: Integer; exp: INTEGER): Integer;
- VAR xL, intL: Index; negative: BOOLEAN; int: Integer;
- BEGIN
- ASSERT(exp >= 0, 20);
- GetLength(x, xL, negative);
- IF xL = 0 THEN int := zero
- ELSIF (xL = 1) & (x[0] = 1) THEN
- IF negative & ~ODD(exp) THEN
- int := New(xL); SetLength(int, xL, FALSE);
- REPEAT DEC(xL); int[xL] := x[xL] UNTIL xL = 0
- ELSE int := x
- END
- ELSIF exp = 0 THEN int := one
- ELSIF exp = 1 THEN int := x
- ELSE
- int := New(SHORT((xL - 1) * exp + ENTIER(Math.Ln(x[xL - 1] + 1) * exp / Math.Ln(B)) + 1));
- BinPower(x, exp, int, xL, intL); SetLength(int, intL, negative & ODD(exp))
- END;
- RETURN int
- END Power;
- (* Read from and Write to String and File *)
- PROCEDURE ConvertFromString* (IN s: ARRAY OF CHAR; OUT x: Integer);
- VAR i, j, k: INTEGER; dig, b: Digit; ch: CHAR; negative: BOOLEAN; new: Integer;
- BEGIN
- i := 0; ch := s[0];
- WHILE (ch # 0X) & (ch <= " ") DO INC(i); ch := s[i] END;
- negative := ch = "-";
- IF negative THEN INC(i); ch := s[i] END;
- IF ch = "+" THEN INC(i); ch := s[i] END;
- WHILE (ch # 0X) & (ch <= " ") DO INC(i); ch := s[i] END;
- ASSERT((ch >= "0") & (ch <= "9"), 20);
- WHILE ch = "0" DO INC(i); ch := s[i] END;
- IF (ch > "0") & (ch <= "9") THEN
- j := i;
- REPEAT INC(j); ch := s[j] UNTIL (ch < "0") OR (ch > "9");
- k := (j - i - 1) DIV DecPerDig + 2;
- new := New(k); SetLength(new, k - 1, negative);
- k := (j - i) MOD DecPerDig;
- IF k # 0 THEN
- b := 1; DEC(k);
- WHILE k # 0 DO DEC(k); b := SHORT(b * 10) END
- ELSE b := B DIV 10
- END;
- REPEAT
- dig := 0;
- WHILE b # 0 DO
- dig := SHORT(dig + b * (ORD(s[i]) - ORD("0"))); b := SHORT(b DIV 10);
- INC(i)
- END;
- new[(j - i) DIV DecPerDig] := dig; b := B DIV 10
- UNTIL i = j;
- x := new
- ELSE x := zero
- END
- END ConvertFromString;
- PROCEDURE ConvertToString* (x: Integer; OUT s: ARRAY OF CHAR);
- VAR j: Index; i: INTEGER; d, b: Digit; negative: BOOLEAN;
- BEGIN
- GetLength(x, j, negative);
- IF negative THEN s[0] := "-"; i := 1 ELSE i := 0 END;
- IF j # 0 THEN
- DEC(j); d := x[j]; b := B DIV 10;
- WHILE d DIV b = 0 DO b := SHORT(b DIV 10) END;
- REPEAT
- s[i] := CHR(d DIV b + ORD("0")); INC(i); d := SHORT(d MOD b); b := SHORT(b DIV 10)
- UNTIL b = 0;
- WHILE j # 0 DO
- DEC(j); d := x[j]; b := B DIV 10;
- REPEAT
- s[i] := CHR(d DIV b + ORD("0")); INC(i); d := SHORT(d MOD b); b := SHORT(b DIV 10)
- UNTIL b = 0
- END
- ELSE s[i] := "0"; INC(i)
- END;
- s[i] := 0X
- END ConvertToString;
- PROCEDURE Internalize* (r: Files.Reader; OUT x: Integer);
- VAR len: Index; n, version: INTEGER; negative: BOOLEAN;
- new: Integer; buf: ARRAY 4 OF BYTE;
- BEGIN
- r.ReadByte(buf[0]); version := buf[0];
- ASSERT((version = 0) OR (version >= 128), 20);
- IF version = 0 THEN
- r.ReadBytes(buf, 0, 4);
- len := (((buf[0] MOD 128) * 256 + buf[1] MOD 256) * 256
- + buf[2] MOD 256) * 256 + buf[3] MOD 256;
- new := New(len); SetLength(new, len, buf[0] < 0);
- WHILE len # 0 DO
- DEC(len);
- r.ReadBytes(buf, 0, 2); new[len] := SHORT((buf[0] MOD 256) * 256 + buf[1] MOD 256)
- END;
- x := new
- ELSE (* version >= 128 *)
- r.ReadByte(buf[1]); n := (buf[0] MOD 256) * 256 + buf[1] MOD 256 - 32768;
- r.ReadBytes(buf, 0, 2); DEC(n);
- len := (buf[0] MOD 256) * 256 + buf[1] MOD 256; negative := len < 0; len := ABS(len);
- new := New(len); SetLength(new, len, negative);
- WHILE n # len DO DEC(n); r.ReadBytes(buf, 0, 2) END;
- WHILE len # 0 DO
- DEC(len);
- r.ReadBytes(buf, 0, 2); new[len] := SHORT((buf[0] MOD 256) * 256 + buf[1] MOD 256)
- END;
- x := new
- END
- END Internalize;
- PROCEDURE Externalize* (w: Files.Writer; x: Integer);
- VAR len, l: Index; d: Digit; i: INTEGER; negative: BOOLEAN; buf: ARRAY 4 OF BYTE;
- PROCEDURE Byte(x: INTEGER): BYTE;
- BEGIN
- ASSERT((x >= MIN(BYTE)) & (x <= MAX(BYTE) - MIN(BYTE)), 20);
- IF x > MAX(BYTE) THEN RETURN SHORT(SHORT(x - 256)) ELSE RETURN SHORT(SHORT(x)) END
- END Byte;
- BEGIN
- GetLength(x, len, negative); l := len; i := 4;
- REPEAT DEC(i); buf[i] := Byte(l MOD 256); l := l DIV 256 UNTIL i = 0;
- IF negative THEN buf[0] := Byte(128 + buf[0] MOD 256) END;
- w.WriteByte(0); w.WriteBytes(buf, 0, 4);
- WHILE len # 0 DO
- DEC(len);
- d := x[len]; buf[0] := Byte(d DIV 256); buf[1] := Byte(d MOD 256); w.WriteBytes(buf, 0, 2)
- END
- END Externalize;
- BEGIN
- ASSERT(B <= BinBase, 20);
- zero := New(0); SetLength(zero, 0, FALSE);
- one := New(1); one[0] := 1; SetLength(one, 1, FALSE);
- two := New(1); two[0] := 2; SetLength(two, 1, FALSE);
- NEW(buf6, 6)
- END Integers.
|