blob: b696563095e3486d15272446751785d27e76a7a0 [file] [log] [blame]
Robert Griesemerdb3bf9c2009-08-14 11:53:27 -07001// Copyright 2009 The Go Authors. All rights reserved.
2// Use of this source code is governed by a BSD-style
3// license that can be found in the LICENSE file.
4
5// This file contains operations on unsigned multi-precision integers.
6// These are the building blocks for the operations on signed integers
7// and rationals.
8
9package big
10
11// An unsigned integer x of the form
12//
13// x = x[n-1]*_B^(n-1) + x[n-2]*_B^(n-2) + ... + x[1]*_B + x[0]
14//
15// with 0 <= x[i] < _B and 0 <= i < n is stored in a slice of length n,
16// with the digits x[i] as the slice elements.
17//
18// A number is normalized if the slice contains no leading 0 digits.
19// During arithmetic operations, denormalized values may occur but are
20// always normalized before returning the final result. The normalized
21// representation of 0 is the empty or nil slice (length = 0).
22
Robert Griesemere5874222009-08-15 11:43:54 -070023// TODO(gri) - convert these routines into methods for type 'nat'
24// - decide if type 'nat' should be exported
25
Robert Griesemerdb3bf9c2009-08-14 11:53:27 -070026func normN(z []Word) []Word {
27 i := len(z);
28 for i > 0 && z[i-1] == 0 {
29 i--;
30 }
31 z = z[0 : i];
32 return z;
33}
34
35
36func makeN(z []Word, m int) []Word {
37 if len(z) > m {
38 z = z[0 : m]; // has at least one extra word for a carry, if any
39 return z; // reuse z
40 }
41 c := 4; // minimum capacity
42 if m > c {
43 c = m;
44 }
45 return make([]Word, m, c+1); // +1: extra word for a carry, if any
46}
47
48
49func newN(z []Word, x uint64) []Word {
50 if x == 0 {
Robert Griesemere5874222009-08-15 11:43:54 -070051 return makeN(z, 0);
Robert Griesemerdb3bf9c2009-08-14 11:53:27 -070052 }
53
54 // single-digit values
55 if x == uint64(Word(x)) {
56 z = makeN(z, 1);
57 z[0] = Word(x);
58 return z;
59 }
60
61 // compute number of words n required to represent x
62 n := 0;
63 for t := x; t > 0; t >>= _W {
64 n++;
65 }
66
67 // split x into n words
68 z = makeN(z, n);
69 for i := 0; i < n; i++ {
70 z[i] = Word(x & _M);
71 x >>= _W;
72 }
73
74 return z;
75}
76
77
78func setN(z, x []Word) []Word {
79 z = makeN(z, len(x));
80 for i, d := range x {
81 z[i] = d;
82 }
83 return z;
84}
85
86
87func addNN(z, x, y []Word) []Word {
88 m := len(x);
89 n := len(y);
90
91 switch {
92 case m < n:
93 return addNN(z, y, x);
94 case m == 0:
95 // n == 0 because m >= n; result is 0
96 return makeN(z, 0);
97 case n == 0:
98 // result is x
99 return setN(z, x);
100 }
Robert Griesemere5874222009-08-15 11:43:54 -0700101 // m > 0
Robert Griesemerdb3bf9c2009-08-14 11:53:27 -0700102
103 z = makeN(z, m);
104 c := addVV(&z[0], &x[0], &y[0], n);
105 if m > n {
106 c = addVW(&z[n], &x[n], c, m-n);
107 }
108 if c > 0 {
109 z = z[0 : m+1];
110 z[m] = c;
111 }
112
113 return z;
114}
115
116
117func subNN(z, x, y []Word) []Word {
118 m := len(x);
119 n := len(y);
120
121 switch {
122 case m < n:
123 panic("underflow");
124 case m == 0:
125 // n == 0 because m >= n; result is 0
126 return makeN(z, 0);
127 case n == 0:
128 // result is x
129 return setN(z, x);
130 }
Robert Griesemere5874222009-08-15 11:43:54 -0700131 // m > 0
Robert Griesemerdb3bf9c2009-08-14 11:53:27 -0700132
133 z = makeN(z, m);
134 c := subVV(&z[0], &x[0], &y[0], n);
135 if m > n {
136 c = subVW(&z[n], &x[n], c, m-n);
137 }
138 if c != 0 {
139 panic("underflow");
140 }
Robert Griesemerdb3bf9c2009-08-14 11:53:27 -0700141 z = normN(z);
Robert Griesemere5874222009-08-15 11:43:54 -0700142
Robert Griesemerdb3bf9c2009-08-14 11:53:27 -0700143 return z;
144}
145
146
147func cmpNN(x, y []Word) int {
148 m := len(x);
149 n := len(y);
150 if m != n || m == 0 {
151 return m-n;
152 }
153
154 i := m-1;
155 for i > 0 && x[i] == y[i] {
156 i--;
157 }
158
159 z := 0;
160 switch {
161 case x[i] < y[i]: z = -1;
162 case x[i] > y[i]: z = 1;
163 }
164 return z;
165}
166
167
Robert Griesemere5874222009-08-15 11:43:54 -0700168func mulAddNWW(z, x []Word, y, r Word) []Word {
Robert Griesemerdb3bf9c2009-08-14 11:53:27 -0700169 m := len(x);
Robert Griesemere5874222009-08-15 11:43:54 -0700170 if m == 0 || y == 0 {
171 return newN(z, uint64(r)); // result is r
Robert Griesemerdb3bf9c2009-08-14 11:53:27 -0700172 }
173 // m > 0
Robert Griesemere5874222009-08-15 11:43:54 -0700174
175 z = makeN(z, m);
176 c := mulAddVWW(&z[0], &x[0], y, r, m);
Robert Griesemerdb3bf9c2009-08-14 11:53:27 -0700177 if c > 0 {
178 z = z[0 : m+1];
179 z[m] = c;
180 }
Robert Griesemere5874222009-08-15 11:43:54 -0700181
Robert Griesemerdb3bf9c2009-08-14 11:53:27 -0700182 return z;
183}
184
185
186func mulNN(z, x, y []Word) []Word {
Robert Griesemere5874222009-08-15 11:43:54 -0700187 m := len(x);
188 n := len(y);
189
190 switch {
191 case m < n:
192 return mulNN(z, x, y);
193 case m == 0 || n == 0:
194 return makeN(z, 0);
195 }
196 // m > 0 && n > 0 && m >= n
197
Robert Griesemerdb3bf9c2009-08-14 11:53:27 -0700198 panic("mulNN unimplemented");
Robert Griesemere5874222009-08-15 11:43:54 -0700199
Robert Griesemerdb3bf9c2009-08-14 11:53:27 -0700200 return z
201}
202
203
204// q = (x-r)/y, with 0 <= r < y
205func divNW(z, x []Word, y Word) (q []Word, r Word) {
206 m := len(x);
207 switch {
208 case y == 0:
209 panic("division by zero");
210 case y == 1:
211 q = setN(z, x); // result is x
212 return;
213 case m == 0:
214 q = setN(z, nil); // result is 0
215 return;
216 }
217 // m > 0
218 z = makeN(z, m);
219 r = divWVW(&z[0], 0, &x[0], y, m);
220 q = normN(z);
221 return;
222}
223
224
225// log2 computes the binary logarithm of x.
226// The result is the integer n for which 2^n <= x < 2^(n+1).
227// If x == 0, the result is < 0.
228func log2(x Word) int {
229 n := 0;
230 for ; x > 0; x >>= 1 {
231 n++;
232 }
233 return n-1;
234}
235
236
237// log2N computes the binary logarithm of x.
238// The result is the integer n for which 2^n <= x < 2^(n+1).
239// If x == 0, the result is < 0.
240func log2N(x []Word) int {
241 m := len(x);
242 if m > 0 {
243 return (m-1)*int(_W) + log2(x[m-1]);
244 }
245 return -1;
246}
247
248
249func hexValue(ch byte) int {
250 var d byte;
251 switch {
252 case '0' <= ch && ch <= '9': d = ch - '0';
253 case 'a' <= ch && ch <= 'f': d = ch - 'a' + 10;
254 case 'A' <= ch && ch <= 'F': d = ch - 'A' + 10;
255 default: return -1;
256 }
257 return int(d);
258}
259
260
261// scanN returns the natural number corresponding to the
262// longest possible prefix of s representing a natural number in a
263// given conversion base, the actual conversion base used, and the
264// prefix length. The syntax of natural numbers follows the syntax
265// of unsigned integer literals in Go.
266//
267// If the base argument is 0, the string prefix determines the actual
268// conversion base. A prefix of ``0x'' or ``0X'' selects base 16; the
269// ``0'' prefix selects base 8. Otherwise the selected base is 10.
270//
271func scanN(z []Word, s string, base int) ([]Word, int, int) {
272 // determine base if necessary
273 i, n := 0, len(s);
274 if base == 0 {
275 base = 10;
276 if n > 0 && s[0] == '0' {
277 if n > 1 && (s[1] == 'x' || s[1] == 'X') {
278 base, i = 16, 2;
279 } else {
280 base, i = 8, 1;
281 }
282 }
283 }
284 if base < 2 || 16 < base {
285 panic("illegal base");
286 }
287
288 // convert string
289 z = makeN(z, len(z));
290 for ; i < n; i++ {
291 d := hexValue(s[i]);
292 if 0 <= d && d < base {
Robert Griesemere5874222009-08-15 11:43:54 -0700293 z = mulAddNWW(z, z, Word(base), Word(d));
Robert Griesemerdb3bf9c2009-08-14 11:53:27 -0700294 } else {
295 break;
296 }
297 }
298
299 return z, base, i;
300}
301
302
303// string converts x to a string for a given base, with 2 <= base <= 16.
304// TODO(gri) in the style of the other routines, perhaps this should take
305// a []byte buffer and return it
306func stringN(x []Word, base int) string {
307 if base < 2 || 16 < base {
308 panic("illegal base");
309 }
310
311 if len(x) == 0 {
312 return "0";
313 }
314
315 // allocate buffer for conversion
316 i := (log2N(x) + 1) / log2(Word(base)) + 1; // +1: round up
317 s := make([]byte, i);
318
319 // don't destroy x
320 q := setN(nil, x);
321
322 // convert
323 for len(q) > 0 {
324 i--;
325 var r Word;
326 q, r = divNW(q, q, 10);
327 s[i] = "0123456789abcdef"[r];
328 };
329
330 return string(s[i : len(s)]);
331}