![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to b1nuI.c CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [CSRG BSD Unix] / 43BSD / contrib / B / src / bint|
Return to b1nuI.c CVS log | Up to [CSRG BSD Unix] / 43BSD / contrib / B / src / bint |
1.1 root 1: /* Copyright (c) Stichting Mathematisch Centrum, Amsterdam, 1985. */
2:
3: /*
4: $Header: b1nuI.c,v 1.4 85/08/22 16:51:13 timo Exp $
5: */
6:
7: /* Multi-precision integer arithmetic */
8:
9: #include "b.h"
10: #include "b1obj.h"
11: #include "b1num.h"
12: #include "b0con.h"
13: #include "b3err.h"
14:
15: /*
16: * Number representation:
17: * ======================
18: *
19: * (Think of BASE = 10 for ordinary decimal notation.)
20: * A number is a sequence of N "digits" b1, b2, ..., bN
21: * where each bi is in {0..BASE-1}, except for negative numbers,
22: * where bN = -1.
23: * The number represented by b1, ..., bN is
24: * b1*BASE**(N-1) + b2*BASE**(N-2) + ... + bN .
25: * The base BASE is chosen so that multiplication of two positive
26: * integers up to BASE-1 can be multiplied exactly using double
27: * precision floating point arithmetic.
28: * Also it must be possible to add two long integers between
29: * -BASE and +BASE (exclusive), giving a result between -2BASE and
30: * +2BASE.
31: * BASE must be even (so we can easily decide whether the whole
32: * number is even), and positive (to avoid all kinds of other trouble).
33: * Presently, it is restricted to a power of 10 by the I/O-conversion
34: * routines (file "b1nuC.c").
35: *
36: * Canonical representation:
37: * bN is never zero (for the number zero itself, N is zero).
38: * If bN is -1, b[N-1] is never BASE-1 .
39: * All operands are assumed te be in canonical representation.
40: * Routine "int_canon" brings a number in canonical representation.
41: *
42: * Mapping to C objects:
43: * A "digit" is an integer of type "digit", probably an "int".
44: * A number is represented as a "B-integer", i.e. something
45: * of type "integer" (which is actually a pointer to some struct).
46: * The number of digits N is extracted through the macro Length(v).
47: * The i-th digit is extracted through the macro Digit(v,N-i).
48: * (So in C, we count in a backwards direction from 0 ... n-1 !)
49: * A number is created through a call to grab_num(N), which sets
50: * N zero digits (thus not in canonical form!).
51: */
52:
53:
54: /*
55: * Bring an integer into canonical form.
56: * Make a SmallInt if at all possible.
57: * NB: Work done by int_canon is duplicated by mk_integer for optimization;
58: * if the strategy here changes, look at mk_integer, too!
59: */
60:
61: Visible integer int_canon(v) integer v; {
62: register int i;
63:
64: if (IsSmallInt(v)) return v;
65:
66: for (i = Length(v) - 1; i >= 0 && Digit(v,i) == 0; --i)
67: ;
68:
69: if (i < 0) {
70: release((value) v);
71: return int_0;
72: }
73:
74: if (i == 0) {
75: digit dig = Digit(v,0);
76: release((value) v);
77: return (integer) MkSmallInt(dig);
78: }
79:
80: if (i > 0 && Digit(v,i) == -1) {
81: while (i > 0 && Digit(v, i-1) == BASE-1) --i;
82: if (i == 0) {
83: release((value) v);
84: return (integer) MkSmallInt(-1);
85: }
86: if (i == 1) {
87: digit dig = Digit(v,0) - BASE;
88: release((value) v);
89: return (integer) MkSmallInt(dig);
90: }
91: Digit(v,i) = -1;
92: }
93:
94: if (i+1 < Length(v)) return (integer) regrab_num((value) v, i+1);
95:
96: return v;
97: }
98:
99:
100: /* General add/subtract subroutine */
101:
102: typedef double twodigit; /* Might be long on 16 bit machines */
103: /* Should be in b0con.h */
104:
105: Hidden twodigit fmodulo(x, y) twodigit x, y; {
106: return x - y * (twodigit) floor((double)x / (double)y);
107: }
108:
109: Visible Procedure dig_gadd(to, nto, from, nfrom, ffactor)
110: digit *to, *from; intlet nto, nfrom; digit ffactor; {
111: twodigit carry= 0;
112: twodigit factor= ffactor;
113: digit save;
114:
115: nto -= nfrom;
116: if (nto < 0)
117: syserr(MESS(1000, "dig_gadd: nto < nfrom"));
118: for (; nfrom > 0; ++to, ++from, --nfrom) {
119: carry += *to + *from * factor;
120: *to= save= fmodulo(carry, (twodigit)BASE);
121: carry= (carry-save) / BASE;
122: }
123: for (; nto > 0; ++to, --nto) {
124: if (carry == 0)
125: return;
126: carry += *to;
127: *to= save= fmodulo(carry, (twodigit)BASE);
128: carry= (carry-save) / BASE;
129: }
130: if (carry != 0)
131: to[-1] += carry*BASE; /* Assume it's -1 */
132: }
133:
134:
135: /* Sum or difference of two integers */
136: /* Should have its own version of dig-gadd without double precision */
137:
138: Visible integer int_gadd(v, w, factor) integer v, w; intlet factor; {
139: struct integer vv, ww;
140: integer s;
141: int len, lenv, i;
142:
143: FreezeSmallInt(v, vv);
144: FreezeSmallInt(w, ww);
145: lenv= len= Length(v);
146: if (Length(w) > len)
147: len= Length(w);
148: ++len;
149: s= (integer) grab_num(len);
150: for (i= 0; i < lenv; ++i)
151: Digit(s, i)= Digit(v, i);
152: for (; i < len; ++i)
153: Digit(s, i)= 0;
154: dig_gadd(&Digit(s, 0), len, &Digit(w, 0), Length(w), (digit)factor);
155: return int_canon(s);
156: }
157:
158:
159: /* Product of two integers */
160:
161: Visible integer int_prod(v, w) integer v, w; {
162: int i;
163: integer a;
164: struct integer vv, ww;
165:
166: if (v == int_0 || w == int_0) return int_0;
167: if (v == int_1) return (integer) Copy(w);
168: if (w == int_1) return (integer) Copy(v);
169:
170: FreezeSmallInt(v, vv);
171: FreezeSmallInt(w, ww);
172:
173: a = (integer) grab_num(Length(v) + Length(w));
174:
175: for (i= Length(a)-1; i >= 0; --i)
176: Digit(a, i)= 0;
177: for (i = 0; i < Length(v) && !interrupted; ++i)
178: dig_gadd(&Digit(a, i), Length(w)+1, &Digit(w, 0), Length(w),
179: Digit(v, i));
180:
181: return int_canon(a);
182: }
183:
184:
185: /* Compare two integers */
186:
187: Visible relation int_comp(v, w) integer v, w; {
188: int sv, sw;
189: register int i;
190: struct integer vv, ww;
191:
192: /* 1. Compare pointers and equal SmallInts */
193: if (v == w) return 0;
194:
195: /* 1a. Handle SmallInts */
196: if (IsSmallInt(v) && IsSmallInt(w))
197: return SmallIntVal(v) - SmallIntVal(w);
198: FreezeSmallInt(v, vv);
199: FreezeSmallInt(w, ww);
200:
201: /* 2. Extract signs */
202: sv = Length(v)==0 ? 0 : Digit(v,Length(v)-1)<0 ? -1 : 1;
203: sw = Length(w)==0 ? 0 : Digit(w,Length(w)-1)<0 ? -1 : 1;
204:
205: /* 3. Compare signs */
206: if (sv != sw) return (sv>sw) - (sv<sw);
207:
208: /* 4. Compare sizes */
209: if (Length(v) != Length(w))
210: return sv * ( (Length(v)>Length(w)) - (Length(v)<Length(w)) );
211:
212: /* 5. Compare individual digits */
213: for (i = Length(v)-1; i >= 0 && Digit(v,i) == Digit(w,i); --i)
214: ;
215:
216: /* 6. All digits equal? */
217: if (i < 0) return 0; /* Yes */
218:
219: /* 7. Compare leftmost different digits */
220: if (Digit(v,i) < Digit(w,i)) return -1;
221:
222: return 1;
223: }
224:
225:
226: /* Construct an integer out of a floating point number */
227:
228: #define GRAN 8 /* Granularity used when requesting more storage */
229: /* MOVE TO MEM! */
230: Visible integer mk_int(x) double x; {
231: register integer a;
232: integer b;
233: register int i, j;
234: int negate;
235:
236: if (MinSmallInt <= x && x <= MaxSmallInt)
237: return (integer) MkSmallInt((int)x);
238:
239: a = (integer) grab_num(1);
240: negate = x < 0 ? 1 : 0;
241: if (negate) x = -x;
242:
243: for (i = 0; x != 0; ++i) {
244: double z = floor(x/BASE);
245: digit save = Modulo((digit)(x-z*BASE), BASE);
246: if (i >= Length(a)) {
247: a = (integer) regrab_num((value) a, Length(a)+GRAN);
248: for (j = Length(a)-1; j > i; --j)
249: Digit(a,j) = 0; /* clear higher digits */
250: }
251: Digit(a,i) = save;
252: x = floor((x-save)/BASE);
253: }
254:
255: if (negate) {
256: b = int_neg(a);
257: release((value) a);
258: return b;
259: }
260:
261: return int_canon(a);
262: }
263:
264: /* Construct an integer out of a C int. Like mk_int, but optimized. */
265:
266: Visible value mk_integer(x) int x; {
267: if (MinSmallInt <= x && x <= MaxSmallInt) return MkSmallInt(x);
268: return (value) mk_int((double)x);
269: }
270:
271:
272: /* Efficiently compute 10**n as a B integer, where n is a C int >= 0 */
273:
274: Visible integer int_tento(n) int n; {
275: integer i;
276: digit msd = 1;
277: if (n < 0) syserr(MESS(1001, "int_tento(-n)"));
278: if (n < tenlogBASE) {
279: while (n != 0) msd *= 10, --n;
280: return (integer) MkSmallInt(msd);
281: }
282: i = (integer) grab_num(1 + (int)(n/tenlogBASE));
283: n %= tenlogBASE;
284: while (n != 0) msd *= 10, --n;
285: Digit(i, Length(i)-1) = msd;
286: return i;
287: }
288:
289: #ifdef NOT_USED
290: /* Approximate ceiling(10 log abs(u/v)), as C int.
291: It only works for v > 0, u, v both integers.
292: The result may be one too large or too small */
293:
294: Visible int scale(u, v) integer u, v; {
295: int s;
296: double z;
297: struct integer uu, vv;
298:
299: if (Msd(v) <= 0) syserr(MESS(1002, "scale(u,v<=0)"));
300: if (u == int_0) return 0; /* `Don't care' case */
301: FreezeSmallInt(u, uu);
302: FreezeSmallInt(v, vv);
303: s = (Length(u) - Length(v)) * tenlogBASE;
304: if (Digit(u, Length(u)-1) >= 0) z = Digit(u, Length(u)-1);
305: else {
306: s -= tenlogBASE;
307: if (Length(u) == 1) z = 1;
308: else z = BASE - Digit(u, Length(u)-2);
309: }
310: z /= Digit(v, Length(v)-1);
311: while (z >= 10) z /= 10, ++s;
312: while (z < 1) z *= 10, --s;
313: return s;
314: }
315: #endif NOT_USED
This archive runs on limited infrastructure. Preserving old code on modern bandwidth. Automated agents are requested to crawl responsibly.