Contact
CoCalc Logo Icon
StoreFeaturesDocsShareSupport News AboutSign UpSign In
| Download

GAP 4.8.9 installation with standard packages -- copy to your CoCalc project to get it

Path: gap4r8 / src / integer.c
Views: 415070
1
/****************************************************************************
2
**
3
*W integer.c GAP source Martin Schönert
4
** & Alice Niemeyer
5
** & Werner Nickel
6
**
7
**
8
*Y Copyright (C) 1996, Lehrstuhl D für Mathematik, RWTH Aachen, Germany
9
*Y (C) 1998 School Math and Comp. Sci., University of St Andrews, Scotland
10
*Y Copyright (C) 2002 The GAP Group
11
**
12
** This file implements the functions handling arbitrary size integers.
13
**
14
** There are three integer types in GAP: 'T_INT', 'T_INTPOS' and 'T_INTNEG'.
15
** Each integer has a unique representation, e.g., an integer that can be
16
** represented as 'T_INT' is never represented as 'T_INTPOS' or 'T_INTNEG'.
17
**
18
** 'T_INT' is the type of those integers small enough to fit into 29 bits.
19
** Therefore the value range of this small integers is $-2^{28}...2^{28}-1$.
20
** This range contains about 99\% of all integers that usually occur in GAP.
21
** (I just made up this number, obviously it depends on the application :-)
22
** Only these small integers can be used as index expression into sequences.
23
**
24
** Small integers are represented by an immediate integer handle, containing
25
** the value instead of pointing to it, which has the following form:
26
**
27
** +-------+-------+-------+-------+- - - -+-------+-------+-------+
28
** | guard | sign | bit | bit | | bit | tag | tag |
29
** | bit | bit | 27 | 26 | | 0 | 0 | 1 |
30
** +-------+-------+-------+-------+- - - -+-------+-------+-------+
31
**
32
** Immediate integers handles carry the tag 'T_INT', i.e. the last bit is 1.
33
** This distinguishes immediate integers from other handles which point to
34
** structures aligned on 4 byte boundaries and therefore have last bit zero.
35
** (The second bit is reserved as tag to allow extensions of this scheme.)
36
** Using immediates as pointers and dereferencing them gives address errors.
37
**
38
** To aid overflow check the most significant two bits must always be equal,
39
** that is to say that the sign bit of immediate integers has a guard bit.
40
**
41
** The macros 'INTOBJ_INT' and 'INT_INTOBJ' should be used to convert between
42
** a small integer value and its representation as immediate integer handle.
43
**
44
** 'T_INTPOS' and 'T_INTNEG' are the types of positive (respectively, negative)
45
** integer values that can not be represented by immediate integers.
46
**
47
** This large integers values are represented in signed base 65536 notation.
48
** That means that the bag of a large integer has the following form:
49
**
50
** +-------+-------+-------+-------+- - - -+-------+-------+-------+
51
** | digit | digit | digit | digit | | digit | digit | digit |
52
** | 0 | 1 | 2 | 3 | | <n>-2 | <n>-1 | <n> |
53
** +-------+-------+-------+-------+- - - -+-------+-------+-------+
54
**
55
** The value of this is: $d0 + d1 65536 + d2 65536^2 + ... + d_n 65536^n$,
56
** respectively the negative of this if the type of this object is 'T_INTNEG'.
57
**
58
** Each digit is of course stored as a 16 bit wide unsigned short.
59
** Note that base 65536 allows us to multiply 2 digits and add a carry digit
60
** without overflow in 32 bit long arithmetic, available on most processors.
61
**
62
** The number of digits in every large integer is a multiple of four.
63
** Therefore the leading three digits of some values will actually be zero.
64
** Note that the uniqueness of representation implies that not four or more
65
** leading digits may be zero, since |d0|d1|d2|d3| and |d0|d1|d2|d3|0|0|0|0|
66
** have the same value only one, the first, can be a legal representation.
67
**
68
** Because of this it is possible to do a little bit of loop unrolling.
69
** Thus instead of looping <n> times, handling one digit in each iteration,
70
** we can loop <n>/4 times, handling four digits during each iteration.
71
** This reduces the overhead of the loop by a factor of approximately four.
72
**
73
** Using base 65536 representation has advantages over using other bases.
74
** Integers in base 65536 representation can be packed dense and therefore
75
** use roughly 20\% less space than integers in base 10000 representation.
76
** 'SumInt' is 20\% and 'ProdInt' is 40\% faster for 65536 than for 10000,
77
** as their runtime is linear respectively quadratic in the number of digits.
78
** Dividing by 65536 and computing the remainder mod 65536 can be done fast
79
** by shifting 16 bit to the right and by taking the lower 16 bits.
80
** Larger bases are difficult because the product of two digits will not fit
81
** into 32 bit, which is the word size of most modern micro processors.
82
** Base 10000 would have the advantage that printing is very much easier,
83
** but 'PrInt' keeps a terminal at 9600 baud busy for almost all integers.
84
*/
85
86
#include "system.h" /* Ints, UInts */
87
88
#ifndef USE_GMP /* use this file, otherwise ignore what follows */
89
90
#include "gasman.h" /* garbage collector */
91
#include "objects.h" /* objects */
92
#include "scanner.h" /* scanner */
93
94
#include "gvars.h" /* global variables */
95
96
#include "calls.h" /* generic call mechanism */
97
#include "opers.h" /* generic operations */
98
99
#include "ariths.h" /* basic arithmetic */
100
101
#include "bool.h" /* booleans */
102
103
#include "integer.h" /* integers */
104
105
#include "gap.h" /* error handling, initialisation */
106
107
#include "records.h" /* generic records */
108
#include "precord.h" /* plain records */
109
110
#include "lists.h" /* generic lists */
111
#include "string.h" /* strings */
112
113
#include "saveload.h" /* saving and loading */
114
115
#include "intfuncs.h"
116
117
#include "code.h" /* coder */
118
#include "thread.h" /* threads */
119
#include "tls.h" /* thread-local storage */
120
121
#include <stdio.h>
122
123
/* for fallbacks to library */
124
Obj String;
125
126
127
/****************************************************************************
128
**
129
*F TypeInt(<int>) . . . . . . . . . . . . . . . . . . . . . type of integer
130
**
131
** 'TypeInt' returns the type of the integer <int>.
132
**
133
** 'TypeInt' is the function in 'TypeObjFuncs' for integers.
134
*/
135
Obj TYPE_INT_SMALL_ZERO;
136
Obj TYPE_INT_SMALL_POS;
137
Obj TYPE_INT_SMALL_NEG;
138
Obj TYPE_INT_LARGE_POS;
139
Obj TYPE_INT_LARGE_NEG;
140
141
Obj TypeIntSmall (
142
Obj val )
143
{
144
if ( 0 == INT_INTOBJ(val) ) {
145
return TYPE_INT_SMALL_ZERO;
146
}
147
else if ( 0 < INT_INTOBJ(val) ) {
148
return TYPE_INT_SMALL_POS;
149
}
150
else /* if ( 0 > INT_INTOBJ(val) ) */ {
151
return TYPE_INT_SMALL_NEG;
152
}
153
}
154
155
Obj TypeIntLargePos (
156
Obj val )
157
{
158
return TYPE_INT_LARGE_POS;
159
}
160
161
Obj TypeIntLargeNeg (
162
Obj val )
163
{
164
return TYPE_INT_LARGE_NEG;
165
}
166
167
/**************************************************************************
168
** The following two functions convert a C Int or UInt respectively into
169
** a GAP integer, either an immediate, small integer if possible or
170
** otherwise a new GAP bag with TNUM T_INTPOS or T_INTNEG.
171
**
172
*F ObjInt_Int(Int i)
173
*F ObjInt_UInt(UInt i)
174
**
175
****************************************************************************/
176
177
Obj ObjInt_Int(Int i)
178
{
179
Obj n;
180
Int bound = 1L << NR_SMALL_INT_BITS;
181
if (i >= bound) {
182
/* We have to make a big integer */
183
n = NewBag(T_INTPOS,4*sizeof(TypDigit));
184
ADDR_INT(n)[0] = (TypDigit) (i & ((Int) INTBASE - 1L));
185
ADDR_INT(n)[1] = (TypDigit) (i >> NR_DIGIT_BITS);
186
ADDR_INT(n)[2] = 0;
187
ADDR_INT(n)[3] = 0;
188
return n;
189
} else if (-i > bound) {
190
n = NewBag(T_INTNEG,4*sizeof(TypDigit));
191
ADDR_INT(n)[0] = (TypDigit) ((-i) & ((Int) INTBASE - 1L));
192
ADDR_INT(n)[1] = (TypDigit) ((-i) >> NR_DIGIT_BITS);
193
ADDR_INT(n)[2] = 0;
194
ADDR_INT(n)[3] = 0;
195
return n;
196
} else {
197
return INTOBJ_INT(i);
198
}
199
}
200
201
Obj ObjInt_UInt(UInt i)
202
{
203
Obj n;
204
UInt bound = 1UL << NR_SMALL_INT_BITS;
205
if (i >= bound) {
206
/* We have to make a big integer */
207
n = NewBag(T_INTPOS,4*sizeof(TypDigit));
208
ADDR_INT(n)[0] = (TypDigit) (i & ((UInt) INTBASE - 1L));
209
ADDR_INT(n)[1] = (TypDigit) (i >> NR_DIGIT_BITS);
210
ADDR_INT(n)[2] = 0;
211
ADDR_INT(n)[3] = 0;
212
return n;
213
} else {
214
return INTOBJ_INT(i);
215
}
216
}
217
218
219
220
/****************************************************************************
221
**
222
*F PrintInt( <int> ) . . . . . . . . . . . . . . . print an integer constant
223
**
224
** 'PrintInt' prints the integer <int> in the usual decimal notation.
225
** 'PrintInt' handles objects of type 'T_INT', 'T_INTPOS' and 'T_INTNEG'.
226
**
227
** Large integers are first converted into base 10000 and then printed.
228
** The time for a conversion depends quadratically on the number of digits.
229
** For 2000 decimal digit integers, a screenfull, it is reasonable fast.
230
**
231
** The number of digits needed in PrIntD[] is the ceiling of the logarithm
232
** with respect to base PRINT_BASE of
233
**
234
** ( (1<<NR_DIGIT_BITS) )^1000 - 1.
235
**
236
** The latter is the largest number that can be represented with 1000 digits
237
** of type TypDigit.
238
**
239
** If NR_DIGIT_BITS is 16, we get 1205.
240
** If NR_DIGIT_BITS is 32, we get 1071.
241
**
242
** The subsidiary function IntToPrintBase converts an integer into base
243
** PRINT_BASE, leaving the result in base PrIntD. It returns the index of the
244
** most significant digits. It assumes that the argument is a large
245
** integer small enough to fit.
246
*/
247
248
TypDigit PrIntC [1000]; /* copy of integer to be printed */
249
250
#ifdef SYS_IS_64_BIT
251
252
#define PRINT_BASE 1000000000L /* 10^9 */
253
#define PRINT_FORMAT "%09d" /* print 9 decimals at a time */
254
#define CHARS_PER_PRINT_BASE 9
255
TypDigit PrIntD [1071]; /* integer converted to base 10^9 */
256
#define NR_HEX_DIGITS 8
257
258
#else
259
260
#define PRINT_BASE 10000
261
#define PRINT_FORMAT "%04d" /* print 4 decimals at a time */
262
#define CHARS_PER_PRINT_BASE 4
263
TypDigit PrIntD [1205]; /* integer converted to base 10000 */
264
#define NR_HEX_DIGITS 4
265
266
#endif
267
/****************************************************************************
268
**
269
*F FuncHexStringInt( <self>, <int> ) . . . . . . . . hex string for integer
270
*F FuncIntHexString( <self>, <string> ) . . . . . . integer from hex string
271
**
272
** The function `FuncHexStringInt' constructs from an integer the
273
** corresponding string in hexadecimal notation. It has a leading '-'
274
** for negative numbers and the digits 10..15 are written as A..F.
275
**
276
** The function `FuncIntHexString' does the converse, but here the
277
** letters a..f are also allowed in <string> instead of A..F.
278
**
279
*/
280
Obj FuncHexStringInt( Obj self, Obj integer )
281
{
282
Int len, i, j, n;
283
UInt nf;
284
TypDigit d, f;
285
UInt1 *p, a;
286
Obj res;
287
288
/* immediate integers */
289
if (IS_INTOBJ(integer)) {
290
n = INT_INTOBJ(integer);
291
/* 0 is special */
292
if (n == 0) {
293
res = NEW_STRING(1);
294
CHARS_STRING(res)[0] = '0';
295
return res;
296
}
297
298
/* else we create a string big enough for any immediate integer */
299
res = NEW_STRING(2 * NR_HEX_DIGITS + 1);
300
p = CHARS_STRING(res);
301
/* handle sign */
302
if (n<0) {
303
p[0] = '-';
304
n = -n;
305
p++;
306
}
307
else
308
SET_LEN_STRING(res, GET_LEN_STRING(res)-1);
309
/* collect digits, skipping leading zeros */
310
j = 0;
311
nf = ((UInt)15) << (4*(2*NR_HEX_DIGITS-1));
312
for (i = 2*NR_HEX_DIGITS; i; i-- ) {
313
a = ((UInt)n & nf) >> (4*(i-1));
314
if (j==0 && a==0) SET_LEN_STRING(res, GET_LEN_STRING(res)-1);
315
else if (a<10) p[j++] = a + '0';
316
else p[j++] = a - 10 + 'A';
317
nf = nf >> 4;
318
}
319
/* final null character */
320
p[j] = 0;
321
return res;
322
}
323
else if (TNUM_OBJ(integer) == T_INTNEG || TNUM_OBJ(integer) == T_INTPOS) {
324
/* nr of digits */
325
len = SIZE_INT(integer);
326
for (; ADDR_INT(integer)[len-1] == 0; len--);
327
328
/* result string and sign */
329
if (TNUM_OBJ(integer) == T_INTNEG) {
330
res = NEW_STRING(len * NR_HEX_DIGITS + 1);
331
p = CHARS_STRING(res);
332
p[0] = '-';
333
p++;
334
}
335
else {
336
res = NEW_STRING(len * NR_HEX_DIGITS);
337
p = CHARS_STRING(res);
338
}
339
/* collect digits */
340
j = 0;
341
for (; len; len--) {
342
d = ADDR_INT(integer)[len-1];
343
f = 15L << (4*(NR_HEX_DIGITS-1));
344
for (i = NR_HEX_DIGITS; i; i-- ) {
345
a = (d & f) >> (4*(i-1));
346
if (j==0 && a==0) SET_LEN_STRING(res, GET_LEN_STRING(res)-1);
347
else if (a<10) p[j++] = a + '0';
348
else p[j++] = a - 10 + 'A';
349
f = f >> 4;
350
}
351
}
352
/* final null character */
353
p[j] = 0;
354
return res;
355
}
356
else
357
ErrorReturnObj("HexStringInt: argument must be integer, (not a %s)",
358
(Int)TNAM_OBJ(integer), 0L,
359
"");
360
return (Obj) 0L; /* please picky cc */
361
}
362
363
Obj FuncIntHexString( Obj self, Obj str )
364
{
365
Obj res;
366
Int i, j, s, ii, len, sign, nd;
367
UInt n;
368
UInt1 *p, a;
369
TypDigit d;
370
371
if (! IsStringConv(str))
372
ErrorReturnObj("IntHexString: argument must be string (not a %s)",
373
(Int)TNAM_OBJ(str), 0L,
374
"");
375
376
/* number of hex digits and sign */
377
len = GET_LEN_STRING(str);
378
if (len == 0) {
379
res = INTOBJ_INT(0);
380
return res;
381
}
382
if (*(CHARS_STRING(str)) == '-') {
383
sign = -1;
384
i = 1;
385
}
386
else {
387
sign = 1;
388
i = 0;
389
}
390
391
/* skip leading zeros */
392
while ((CHARS_STRING(str))[i] == '0' && i < len)
393
i++;
394
395
/* small int case */
396
if ((len-i)*4 <= NR_SMALL_INT_BITS) {
397
n = 0;
398
p = CHARS_STRING(str);
399
for (; i<len; i++) {
400
a = p[i];
401
if (a>96)
402
a -= 87;
403
else if (a>64)
404
a -= 55;
405
else
406
a -= 48;
407
if (a > 15)
408
ErrorReturnObj("IntHexString: non-valid character in hex-string",
409
0L, 0L, "");
410
n = (n << 4) + a;
411
}
412
res = INTOBJ_INT(sign * n);
413
return res;
414
}
415
else {
416
/* number of Digits */
417
nd = (len-i)/NR_HEX_DIGITS;
418
if (nd * NR_HEX_DIGITS < (len-i)) nd++;
419
nd += ((3*nd) % 4);
420
if (sign == 1)
421
res = NewBag( T_INTPOS, nd*sizeof(TypDigit) );
422
else
423
res = NewBag( T_INTNEG, nd*sizeof(TypDigit) );
424
/* collect digits, easiest to start from the end */
425
p = CHARS_STRING(str);
426
for (j=0; j < nd; j++) {
427
d = 0;
428
for (s=0, ii=len-j*NR_HEX_DIGITS-1;
429
ii>=i && ii>len-(j+1)*NR_HEX_DIGITS-1;
430
s+=4, ii--) {
431
a = p[ii];
432
if (a>96)
433
a -= 87;
434
else if (a>64)
435
a -= 55;
436
else
437
a -= 48;
438
if (a > 15)
439
ErrorReturnObj("IntHexString: non-valid character in hex-string",
440
0L, 0L, "");
441
442
d += (a<<s);
443
}
444
ADDR_INT(res)[j] = d;
445
}
446
return res;
447
}
448
}
449
450
451
452
Int IntToPrintBase ( Obj op )
453
{
454
UInt i, k; /* loop counter */
455
TypDigit * p; /* loop pointer */
456
UInt c; /* carry in division step */
457
458
i = 0;
459
for ( k = 0; k < SIZE_INT(op); k++ )
460
PrIntC[k] = ADDR_INT(op)[k];
461
while ( k > 0 && PrIntC[k-1] == 0 ) k--;
462
while ( k > 0 ) {
463
for ( c = 0, p = PrIntC+k-1; p >= PrIntC; p-- ) {
464
c = (c<<NR_DIGIT_BITS) + *p;
465
*p = (TypDigit)(c / PRINT_BASE);
466
c = c - PRINT_BASE * *p;
467
}
468
PrIntD[i++] = (TypDigit)c;
469
while ( k > 0 && PrIntC[k-1] == 0 ) k--;
470
}
471
return i-1;
472
473
}
474
475
void PrintInt (
476
Obj op )
477
{
478
Int i; /* loop counter */
479
Obj str; /* fallback to lib for large ints */
480
481
/* print a small integer */
482
if ( IS_INTOBJ(op) ) {
483
Pr( "%>%d%<", INT_INTOBJ(op), 0L );
484
}
485
486
/* print a large integer */
487
else if ( SIZE_INT(op) < 1000 ) {
488
489
/* start printing, %> means insert '\' before a linebreak */
490
Pr("%>",0L,0L);
491
492
if ( TNUM_OBJ(op) == T_INTNEG )
493
Pr("-",0L,0L);
494
495
/* convert the integer into base PRINT_BASE */
496
i = IntToPrintBase(op);
497
498
/* print the base PRINT_BASE digits */
499
Pr( "%d", (Int)PrIntD[i], 0L );
500
while ( i > 0 )
501
Pr( PRINT_FORMAT, (Int)PrIntD[--i], 0L );
502
Pr("%<",0L,0L);
503
504
}
505
506
else {
507
str = CALL_1ARGS( String, op );
508
Pr("%>%s%<",(Int)(CHARS_STRING(str)), 0);
509
/* for a long time Print of large ints did not follow the general idea
510
* that Print should produce something that can be read back into GAP:
511
Pr("<<an integer too large to be printed>>",0L,0L); */
512
}
513
}
514
515
/****************************************************************************
516
**
517
** Implementation of Log2Int for C integers.
518
*/
519
520
Int CLog2Int(Int a)
521
{
522
Int res, mask;
523
if (a < 0) a = -a;
524
if (a < 1) return -1;
525
if (a < 65536) {
526
for(mask = 2, res = 0; ;mask *= 2, res += 1) {
527
if(a < mask) return res;
528
}
529
}
530
for(mask = 65536, res = 15; ;mask *= 2, res += 1) {
531
if(a < mask) return res;
532
}
533
}
534
535
/****************************************************************************
536
**
537
*F FuncLog2Int( <self>, <int> ) . . . . . . . . . nr of bits of integer - 1
538
**
539
** Given to GAP-Level as "Log2Int".
540
*/
541
Obj FuncLog2Int( Obj self, Obj integer)
542
{
543
Int d;
544
Int a, len;
545
TypDigit dmask;
546
547
/* case of small ints */
548
if (IS_INTOBJ(integer)) {
549
return INTOBJ_INT(CLog2Int(INT_INTOBJ(integer)));
550
}
551
552
/* case of long ints */
553
if (TNUM_OBJ(integer) == T_INTNEG || TNUM_OBJ(integer) == T_INTPOS) {
554
for (len = SIZE_INT(integer); ADDR_INT(integer)[len-1] == 0; len--);
555
/* Instead of computing
556
res = len * NR_DIGIT_BITS - d;
557
we keep len and d separate, because on 32 bit systems res may
558
not fit into an Int (and not into an immediate integer). */
559
d = 1;
560
a = (TypDigit)(ADDR_INT(integer)[len-1]);
561
for(dmask = (TypDigit)1 << (NR_DIGIT_BITS - 1);
562
(dmask & a) == 0 && dmask != (TypDigit)0;
563
dmask = dmask >> 1, d++);
564
return DiffInt(ProdInt(INTOBJ_INT(len), INTOBJ_INT(NR_DIGIT_BITS)),
565
INTOBJ_INT(d));
566
}
567
else {
568
ErrorReturnObj("Log2Int: argument must be integer, (not a %s)",
569
(Int)TNAM_OBJ(integer), 0L,
570
"");
571
return (Obj) 0L; /* please picky cc */
572
}
573
}
574
575
/****************************************************************************
576
**
577
*F FuncSTRING_INT( <self>, <int> ) . . . . . convert an integer to a string
578
**
579
** `FuncSTRING_INT' returns an immutable string representing the integer
580
** <int>
581
**
582
*/
583
584
Obj FuncSTRING_INT( Obj self, Obj integer )
585
{
586
Int x;
587
Obj str;
588
Int len;
589
Int i;
590
Char c;
591
Int j,top, chunk, neg;
592
593
/* handle a small integer */
594
if ( IS_INTOBJ(integer) ) {
595
x = INT_INTOBJ(integer);
596
str = NEW_STRING( (NR_SMALL_INT_BITS+5)/3 );
597
RetypeBag(str, T_STRING+IMMUTABLE);
598
len = 0;
599
/* Case of zero */
600
if (x == 0)
601
{
602
CHARS_STRING(str)[0] = '0';
603
CHARS_STRING(str)[1] = '\0';
604
ResizeBag(str, SIZEBAG_STRINGLEN(1));
605
SET_LEN_STRING(str, 1);
606
607
return str;
608
}
609
/* Negative numbers */
610
if (x < 0)
611
{
612
CHARS_STRING(str)[len++] = '-';
613
x = -x;
614
neg = 1;
615
}
616
else
617
neg = 0;
618
619
/* Now the main case */
620
while (x != 0)
621
{
622
CHARS_STRING(str)[len++] = '0'+ x % 10;
623
x /= 10;
624
}
625
CHARS_STRING(str)[len] = '\0';
626
627
/* finally, reverse the digits in place */
628
for (i = neg; i < (neg+len)/2; i++)
629
{
630
c = CHARS_STRING(str)[neg+len-1-i];
631
CHARS_STRING(str)[neg+len-1-i] = CHARS_STRING(str)[i];
632
CHARS_STRING(str)[i] = c;
633
}
634
635
ResizeBag(str, SIZEBAG_STRINGLEN(len));
636
SET_LEN_STRING(str, len);
637
return str;
638
}
639
640
/* handle a large integer */
641
else if ( SIZE_INT(integer) < 1000 ) {
642
643
/* convert the integer into base PRINT_BASE */
644
len = IntToPrintBase(integer);
645
str = NEW_STRING(CHARS_PER_PRINT_BASE*(len+1)+2);
646
RetypeBag(str, T_STRING+IMMUTABLE);
647
648
/* sort out the length of the top group */
649
j = 1;
650
top = (Int)PrIntD[len];
651
while ( top >= j)
652
{
653
j *= 10;
654
}
655
656
/* Start filling in the string */
657
i = 0;
658
if ( TNUM_OBJ(integer) == T_INTNEG ) {
659
CHARS_STRING(str)[i++] = '-';
660
}
661
662
while (j > 1)
663
{
664
j /= 10;
665
CHARS_STRING(str)[i++] = '0' + (top / j) % 10;
666
}
667
668
/* Now the rest of the base PRINT_BASE digits are easy */
669
while( len > 0)
670
{
671
chunk = (Int)PrIntD[--len];
672
j = PRINT_BASE/10;
673
while (j > 0)
674
{
675
CHARS_STRING(str)[i++] = '0' + (chunk / j) % 10;
676
j /= 10;
677
}
678
}
679
680
CHARS_STRING(str)[i] = '\0';
681
ResizeBag(str, SIZEBAG_STRINGLEN(i));
682
SET_LEN_STRING(str, i);
683
return str;
684
}
685
else {
686
687
/* Very large integer, fall back on the GAP function */
688
return CALL_1ARGS( String, integer);
689
}
690
}
691
692
693
694
/****************************************************************************
695
**
696
*F EqInt( <intL>, <intR> ) . . . . . . . . . test if two integers are equal
697
**
698
** 'EqInt' returns 1 if the two integer arguments <intL> and <intR> are
699
** equal and 0 otherwise.
700
*/
701
Int EqInt (
702
Obj opL,
703
Obj opR )
704
{
705
Int k; /* loop counter */
706
TypDigit * l; /* pointer into the left operand */
707
TypDigit * r; /* pointer into the right operand */
708
709
/* compare two small integers */
710
if ( ARE_INTOBJS( opL, opR ) ) {
711
if ( INT_INTOBJ(opL) == INT_INTOBJ(opR) ) return 1L;
712
else return 0L;
713
}
714
715
/* compare a small and a large integer */
716
else if ( IS_INTOBJ(opL) ) {
717
return 0L;
718
}
719
else if ( IS_INTOBJ(opR) ) {
720
return 0L;
721
}
722
723
/* compare two large integers */
724
else {
725
726
/* compare the sign and size */
727
if ( TNUM_OBJ(opL) != TNUM_OBJ(opR)
728
|| SIZE_INT(opL) != SIZE_INT(opR) )
729
return 0L;
730
731
/* set up the pointers */
732
l = ADDR_INT(opL);
733
r = ADDR_INT(opR);
734
735
/* run through the digits, four at a time */
736
for ( k = SIZE_INT(opL)/4-1; k >= 0; k-- ) {
737
if ( *l++ != *r++ ) return 0L;
738
if ( *l++ != *r++ ) return 0L;
739
if ( *l++ != *r++ ) return 0L;
740
if ( *l++ != *r++ ) return 0L;
741
}
742
743
/* no differences found, so they must be equal */
744
return 1L;
745
746
}
747
}
748
749
750
/****************************************************************************
751
**
752
*F LtInt( <intL>, <intR> ) . . . . . test if an integer is less than another
753
**
754
** 'LtInt' returns 1 if the integer <intL> is strictly less than the integer
755
** <intR> and 0 otherwise.
756
*/
757
Int LtInt (
758
Obj opL,
759
Obj opR )
760
{
761
Int k; /* loop counter */
762
TypDigit * l; /* pointer into the left operand */
763
TypDigit * r; /* pointer into the right operand */
764
765
/* compare two small integers */
766
if ( ARE_INTOBJS( opL, opR ) ) {
767
if ( INT_INTOBJ(opL) < INT_INTOBJ(opR) ) return 1L;
768
else return 0L;
769
}
770
771
/* compare a small and a large integer */
772
else if ( IS_INTOBJ(opL) ) {
773
if ( TNUM_OBJ(opR) == T_INTPOS ) return 1L;
774
else return 0L;
775
}
776
else if ( IS_INTOBJ(opR) ) {
777
if ( TNUM_OBJ(opL) == T_INTPOS ) return 0L;
778
else return 1L;
779
}
780
781
/* compare two large integers */
782
else {
783
784
/* compare the sign and size */
785
if ( TNUM_OBJ(opL) == T_INTNEG
786
&& TNUM_OBJ(opR) == T_INTPOS )
787
return 1L;
788
else if ( TNUM_OBJ(opL) == T_INTPOS
789
&& TNUM_OBJ(opR) == T_INTNEG )
790
return 0L;
791
else if ( (TNUM_OBJ(opL) == T_INTPOS
792
&& SIZE_INT(opL) < SIZE_INT(opR))
793
|| (TNUM_OBJ(opL) == T_INTNEG
794
&& SIZE_INT(opL) > SIZE_INT(opR)) )
795
return 1L;
796
else if ( (TNUM_OBJ(opL) == T_INTPOS
797
&& SIZE_INT(opL) > SIZE_INT(opR))
798
|| (TNUM_OBJ(opL) == T_INTNEG
799
&& SIZE_INT(opL) < SIZE_INT(opR)) )
800
return 0L;
801
802
/* set up the pointers */
803
l = ADDR_INT(opL);
804
r = ADDR_INT(opR);
805
806
/* run through the digits, from the end downwards */
807
for ( k = SIZE_INT(opL)-1; k >= 0; k-- ) {
808
if ( l[k] != r[k] ) {
809
if ( (TNUM_OBJ(opL) == T_INTPOS
810
&& l[k] < r[k])
811
|| (TNUM_OBJ(opL) == T_INTNEG
812
&& l[k] > r[k]) )
813
return 1L;
814
else
815
return 0L;
816
}
817
}
818
819
/* no differences found, so they must be equal */
820
return 0L;
821
822
}
823
}
824
825
826
/****************************************************************************
827
**
828
*F SumInt( <intL>, <intR> ) . . . . . . . . . . . . . . sum of two integers
829
**
830
** 'SumInt' returns the sum of the two integer arguments <intL> and <intR>.
831
** 'SumInt' handles operands of type 'T_INT', 'T_INTPOS' and 'T_INTNEG'.
832
**
833
** It can also be used in the cases that both operands are small integers
834
** and the result is a small integer too, i.e., that no overflow occurs.
835
** This case is usually already handled in 'EvalSum' for a better efficiency.
836
**
837
** Is called from the 'EvalSum' binop so both operands are already evaluated.
838
**
839
** 'SumInt' is a little bit difficult since there are 16 different cases to
840
** handle, each operand can be positive or negative, small or large integer.
841
** If the operands have opposite sign 'SumInt' calls 'DiffInt', this helps
842
** reduce the total amount of code by a factor of two.
843
*/
844
Obj SumInt (
845
Obj opL,
846
Obj opR )
847
{
848
Int i; /* loop variable */
849
Int k; /* loop variable */
850
Int cs; /* sum of two smalls */
851
UInt c; /* sum of two digits */
852
TypDigit * l; /* pointer into the left operand */
853
TypDigit * r; /* pointer into the right operand */
854
TypDigit * s; /* pointer into the sum */
855
UInt * l2; /* pointer to get 2 digits at once */
856
UInt * s2; /* pointer to put 2 digits at once */
857
Obj sum; /* handle of the result bag */
858
859
/* adding two small integers */
860
if ( ARE_INTOBJS( opL, opR ) ) {
861
862
/* add two small integers with a small sum */
863
/* add and compare top two bits to check that no overflow occured */
864
if ( SUM_INTOBJS( sum, opL, opR ) ) {
865
return sum;
866
}
867
868
/* add two small integers with a large sum */
869
cs = INT_INTOBJ(opL) + INT_INTOBJ(opR);
870
if ( 0 < cs ) {
871
sum = NewBag( T_INTPOS, 4*sizeof(TypDigit) );
872
ADDR_INT(sum)[0] = (TypDigit)cs;
873
ADDR_INT(sum)[1] = (TypDigit)(cs >> NR_DIGIT_BITS);
874
}
875
else {
876
sum = NewBag( T_INTNEG, 4*sizeof(TypDigit) );
877
ADDR_INT(sum)[0] = (TypDigit)(-cs);
878
ADDR_INT(sum)[1] = (TypDigit)((-cs) >> NR_DIGIT_BITS);
879
}
880
881
}
882
883
/* adding one large integer and one small integer */
884
else if ( IS_INTOBJ(opL) || IS_INTOBJ(opR) ) {
885
886
/* make the right operand the small one */
887
if ( IS_INTOBJ(opL) ) {
888
sum = opL; opL = opR; opR = sum;
889
}
890
891
/* if the integers have different sign, let 'DiffInt' do the work */
892
if ( (TNUM_OBJ(opL) == T_INTNEG && 0 <= INT_INTOBJ(opR))
893
|| (TNUM_OBJ(opL) == T_INTPOS && INT_INTOBJ(opR) < 0) ) {
894
if ( TNUM_OBJ(opL) == T_INTPOS ) RetypeBag( opL, T_INTNEG );
895
else RetypeBag( opL, T_INTPOS );
896
sum = DiffInt( opR, opL );
897
if ( TNUM_OBJ(opL) == T_INTPOS ) RetypeBag( opL, T_INTNEG );
898
else RetypeBag( opL, T_INTPOS );
899
return sum;
900
}
901
902
/* allocate the result bag and set up the pointers */
903
if ( TNUM_OBJ(opL) == T_INTPOS ) {
904
i = INT_INTOBJ(opR);
905
sum = NewBag( T_INTPOS, (SIZE_INT(opL)+4)*sizeof(TypDigit) );
906
}
907
else {
908
i = -INT_INTOBJ(opR);
909
sum = NewBag( T_INTNEG, (SIZE_INT(opL)+4)*sizeof(TypDigit) );
910
}
911
l = ADDR_INT(opL);
912
s = ADDR_INT(sum);
913
914
/* add the first four digits,the right operand has only two digits */
915
c = (UInt)*l++ + (TypDigit)i; *s++ = (TypDigit)c;
916
c = (UInt)*l++ + (i>>NR_DIGIT_BITS) + (c>>NR_DIGIT_BITS); *s++ = (TypDigit)c;
917
c = (UInt)*l++ + (c>>NR_DIGIT_BITS); *s++ = (TypDigit)c;
918
c = (UInt)*l++ + (c>>NR_DIGIT_BITS); *s++ = (TypDigit)c;
919
920
/* propagate the carry, this loop is almost never executed */
921
for ( k = SIZE_INT(opL)/4-1; k != 0 && (c>>NR_DIGIT_BITS) != 0; k-- ) {
922
c = (UInt)*l++ + (c>>NR_DIGIT_BITS); *s++ = (TypDigit)c;
923
c = (UInt)*l++ + (c>>NR_DIGIT_BITS); *s++ = (TypDigit)c;
924
c = (UInt)*l++ + (c>>NR_DIGIT_BITS); *s++ = (TypDigit)c;
925
c = (UInt)*l++ + (c>>NR_DIGIT_BITS); *s++ = (TypDigit)c;
926
}
927
928
/* just copy the remaining digits, do it two digits at once */
929
for ( l2 = (UInt*)l, s2 = (UInt*)s; k != 0; k-- ) {
930
*s2++ = *l2++;
931
*s2++ = *l2++;
932
}
933
934
/* if there is a carry, enter it, otherwise shrink the sum */
935
if ( (c>>NR_DIGIT_BITS) != 0 )
936
*s++ = (TypDigit)(c>>NR_DIGIT_BITS);
937
else
938
ResizeBag( sum, (SIZE_INT(sum)-4)*sizeof(TypDigit) );
939
940
}
941
942
/* add two large integers */
943
else {
944
945
/* if the integers have different sign, let 'DiffInt' do the work */
946
if ( (TNUM_OBJ(opL) == T_INTPOS && TNUM_OBJ(opR) == T_INTNEG)
947
|| (TNUM_OBJ(opL) == T_INTNEG && TNUM_OBJ(opR) == T_INTPOS) ) {
948
if ( TNUM_OBJ(opL) == T_INTPOS ) RetypeBag( opL, T_INTNEG );
949
else RetypeBag( opL, T_INTPOS );
950
sum = DiffInt( opR, opL );
951
if ( TNUM_OBJ(opL) == T_INTPOS ) RetypeBag( opL, T_INTNEG );
952
else RetypeBag( opL, T_INTPOS );
953
return sum;
954
}
955
956
/* make the right operand the smaller one */
957
if ( SIZE_INT(opL) < SIZE_INT(opR) ) {
958
sum = opL; opL = opR; opR = sum;
959
}
960
961
/* allocate the result bag and set up the pointers */
962
if ( TNUM_OBJ(opL) == T_INTPOS ) {
963
sum = NewBag( T_INTPOS, (SIZE_INT(opL)+4)*sizeof(TypDigit) );
964
}
965
else {
966
sum = NewBag( T_INTNEG, (SIZE_INT(opL)+4)*sizeof(TypDigit) );
967
}
968
l = ADDR_INT(opL);
969
r = ADDR_INT(opR);
970
s = ADDR_INT(sum);
971
972
/* add the digits, convert to UInt to get maximum precision */
973
c = 0;
974
for ( k = SIZE_INT(opR)/4; k != 0; k-- ) {
975
c = (UInt)*l++ + (UInt)*r++ + (c>>NR_DIGIT_BITS); *s++ = (TypDigit)c;
976
c = (UInt)*l++ + (UInt)*r++ + (c>>NR_DIGIT_BITS); *s++ = (TypDigit)c;
977
c = (UInt)*l++ + (UInt)*r++ + (c>>NR_DIGIT_BITS); *s++ = (TypDigit)c;
978
c = (UInt)*l++ + (UInt)*r++ + (c>>NR_DIGIT_BITS); *s++ = (TypDigit)c;
979
}
980
981
/* propagate the carry, this loop is almost never executed */
982
for ( k=(SIZE_INT(opL)-SIZE_INT(opR))/4;
983
k!=0 && (c>>NR_DIGIT_BITS)!=0; k-- ) {
984
c = (UInt)*l++ + (c>>NR_DIGIT_BITS); *s++ = (TypDigit)c;
985
c = (UInt)*l++ + (c>>NR_DIGIT_BITS); *s++ = (TypDigit)c;
986
c = (UInt)*l++ + (c>>NR_DIGIT_BITS); *s++ = (TypDigit)c;
987
c = (UInt)*l++ + (c>>NR_DIGIT_BITS); *s++ = (TypDigit)c;
988
}
989
990
/* just copy the remaining digits, do it two digits at once */
991
for ( l2 = (UInt*)l, s2 = (UInt*)s; k != 0; k-- ) {
992
*s2++ = *l2++;
993
*s2++ = *l2++;
994
}
995
996
/* if there is a carry, enter it, otherwise shrink the sum */
997
if ( (c>>NR_DIGIT_BITS) != 0 )
998
*s++ = (TypDigit)(c>>NR_DIGIT_BITS);
999
else
1000
ResizeBag( sum, (SIZE_INT(sum)-4)*sizeof(TypDigit) );
1001
1002
}
1003
1004
/* return the sum */
1005
return sum;
1006
}
1007
1008
1009
/****************************************************************************
1010
**
1011
*F ZeroInt(<int>) . . . . . . . . . . . . . . . . . . . . zero of integers
1012
*/
1013
Obj ZeroInt (
1014
Obj op )
1015
{
1016
return INTOBJ_INT(0);
1017
}
1018
1019
1020
/****************************************************************************
1021
**
1022
*F AInvInt(<int>) . . . . . . . . . . . . . additive inverse of an integer
1023
*/
1024
Obj AInvInt (
1025
Obj op )
1026
{
1027
Obj inv;
1028
UInt i;
1029
1030
/* handle small integer */
1031
if ( IS_INTOBJ( op ) ) {
1032
1033
/* special case (ugh) */
1034
if ( op == INTOBJ_INT( -(1L<<NR_SMALL_INT_BITS) ) ) {
1035
inv = NewBag( T_INTPOS, 4*sizeof(TypDigit) );
1036
ADDR_INT(inv)[0] = 0;
1037
ADDR_INT(inv)[1] = (TypDigit)(1L<<(NR_SMALL_INT_BITS-NR_DIGIT_BITS));
1038
}
1039
1040
/* general case */
1041
else {
1042
inv = INTOBJ_INT( - INT_INTOBJ( op ) );
1043
}
1044
1045
}
1046
1047
/* invert a large integer */
1048
else {
1049
1050
/* special case (ugh) */
1051
if ( TNUM_OBJ(op) == T_INTPOS && SIZE_INT(op) == 4
1052
&& ADDR_INT(op)[3] == 0
1053
&& ADDR_INT(op)[2] == 0
1054
&& ADDR_INT(op)[1] == (1L<<(NR_SMALL_INT_BITS-NR_DIGIT_BITS))
1055
&& ADDR_INT(op)[0] == 0 ) {
1056
inv = INTOBJ_INT( -(1L<<NR_SMALL_INT_BITS) );
1057
}
1058
1059
/* general case */
1060
else {
1061
if ( TNUM_OBJ(op) == T_INTPOS ) {
1062
inv = NewBag( T_INTNEG, SIZE_OBJ(op) );
1063
}
1064
else {
1065
inv = NewBag( T_INTPOS, SIZE_OBJ(op) );
1066
}
1067
for ( i = 0; i < SIZE_INT(op); i++ ) {
1068
ADDR_INT(inv)[i] = ADDR_INT(op)[i];
1069
}
1070
}
1071
1072
}
1073
1074
/* return the inverse */
1075
return inv;
1076
}
1077
1078
1079
/****************************************************************************
1080
**
1081
*F DiffInt( <intL>, <intR> ) . . . . . . . . . . difference of two integers
1082
**
1083
** 'DiffInt' returns the difference of the two integer arguments <intL> and
1084
** <intR>. 'DiffInt' handles operands of type 'T_INT', 'T_INTPOS' and
1085
** 'T_INTNEG'.
1086
**
1087
** It can also be used in the cases that both operands are small integers
1088
** and the result is a small integer too, i.e., that no overflow occurs.
1089
** This case is usually already handled in 'EvalDiff' for a better efficiency.
1090
**
1091
** Is called from the 'EvalDiff' binop so both operands are already evaluated.
1092
**
1093
** 'DiffInt' is a little bit difficult since there are 16 different cases to
1094
** handle, each operand can be positive or negative, small or large integer.
1095
** If the operands have opposite sign 'DiffInt' calls 'SumInt', this helps
1096
** reduce the total amount of code by a factor of two.
1097
*/
1098
Obj DiffInt (
1099
Obj opL,
1100
Obj opR )
1101
{
1102
Int i; /* loop variable */
1103
Int k; /* loop variable */
1104
Int c; /* difference of two digits */
1105
TypDigit * l; /* pointer into the left operand */
1106
TypDigit * r; /* pointer into the right operand */
1107
TypDigit * d; /* pointer into the difference */
1108
UInt * l2; /* pointer to get 2 digits at once */
1109
UInt * d2; /* pointer to put 2 digits at once */
1110
Obj dif; /* handle of the result bag */
1111
1112
/* subtracting two small integers */
1113
if ( ARE_INTOBJS( opL, opR ) ) {
1114
1115
/* subtract two small integers with a small difference */
1116
/* sub and compare top two bits to check that no overflow occured */
1117
if ( DIFF_INTOBJS( dif, opL, opR ) ) {
1118
return dif;
1119
}
1120
1121
/* subtract two small integers with a large difference */
1122
c = INT_INTOBJ(opL) - INT_INTOBJ(opR);
1123
if ( 0 < c ) {
1124
dif = NewBag( T_INTPOS, 4*sizeof(TypDigit) );
1125
ADDR_INT(dif)[0] = (TypDigit)c;
1126
ADDR_INT(dif)[1] = (TypDigit)(c >> NR_DIGIT_BITS);
1127
}
1128
else {
1129
dif = NewBag( T_INTNEG, 4*sizeof(TypDigit) );
1130
ADDR_INT(dif)[0] = (TypDigit)(-c);
1131
ADDR_INT(dif)[1] = (TypDigit)((-c) >> NR_DIGIT_BITS);
1132
}
1133
1134
}
1135
1136
/* subtracting one small integer and one large integer */
1137
else if ( IS_INTOBJ( opL ) || IS_INTOBJ( opR ) ) {
1138
1139
/* make the right operand the small one */
1140
if ( IS_INTOBJ( opL ) ) {
1141
dif = opL; opL = opR; opR = dif;
1142
c = -1;
1143
}
1144
else {
1145
c = 1;
1146
}
1147
1148
/* if the integers have different sign, let 'SumInt' do the work */
1149
if ( (TNUM_OBJ(opL) == T_INTNEG && 0 <= INT_INTOBJ(opR))
1150
|| (TNUM_OBJ(opL) == T_INTPOS && INT_INTOBJ(opR) < 0) ) {
1151
if ( TNUM_OBJ(opL) == T_INTPOS ) RetypeBag( opL, T_INTNEG );
1152
else RetypeBag( opL, T_INTPOS );
1153
dif = SumInt( opL, opR );
1154
if ( TNUM_OBJ(opL) == T_INTPOS ) RetypeBag( opL, T_INTNEG );
1155
else RetypeBag( opL, T_INTPOS );
1156
if ( c == 1 ) {
1157
if ( TNUM_OBJ(dif) == T_INTPOS ) RetypeBag( dif, T_INTNEG );
1158
else RetypeBag( dif, T_INTPOS );
1159
}
1160
return dif;
1161
}
1162
1163
/* allocate the result bag and set up the pointers */
1164
if ( TNUM_OBJ(opL) == T_INTPOS ) {
1165
i = INT_INTOBJ(opR);
1166
if ( c == 1 ) dif = NewBag( T_INTPOS, SIZE_OBJ(opL) );
1167
else dif = NewBag( T_INTNEG, SIZE_OBJ(opL) );
1168
}
1169
else {
1170
i = - INT_INTOBJ(opR);
1171
if ( c == 1 ) dif = NewBag( T_INTNEG, SIZE_OBJ(opL) );
1172
else dif = NewBag( T_INTPOS, SIZE_OBJ(opL) );
1173
}
1174
l = ADDR_INT(opL);
1175
d = ADDR_INT(dif);
1176
1177
/* sub the first four digit, note the left operand has only two */
1178
/*N (c>>16<) need not work, replace by (c<0?-1:0) */
1179
c = (Int)*l++ - (TypDigit)i; *d++ = (TypDigit)c;
1180
c = (Int)*l++ - (TypDigit)(i/(1L<<NR_DIGIT_BITS)) + (c<0?-1:0); *d++ = (TypDigit)c;
1181
c = (Int)*l++ + (c<0?-1:0); *d++ = (TypDigit)c;
1182
c = (Int)*l++ + (c<0?-1:0); *d++ = (TypDigit)c;
1183
1184
/* propagate the carry, this loop is almost never executed */
1185
for ( k = SIZE_INT(opL)/4-1; k != 0 && c < 0; k-- ) {
1186
c = (Int)*l++ + (c < 0 ? -1 : 0); *d++ = (TypDigit)c;
1187
c = (Int)*l++ + (c < 0 ? -1 : 0); *d++ = (TypDigit)c;
1188
c = (Int)*l++ + (c < 0 ? -1 : 0); *d++ = (TypDigit)c;
1189
c = (Int)*l++ + (c < 0 ? -1 : 0); *d++ = (TypDigit)c;
1190
}
1191
1192
/* just copy the remaining digits, do it two digits at once */
1193
for ( l2 = (UInt*)l, d2 = (UInt*)d; k != 0; k-- ) {
1194
*d2++ = *l2++;
1195
*d2++ = *l2++;
1196
}
1197
1198
/* no underflow since we subtracted a small int from a large one */
1199
/* but there may be leading zeroes in the result, get rid of them */
1200
/* occurs almost never, so it doesn't matter that it is expensive */
1201
if ( ((UInt*)d == d2
1202
&& d[-4] == 0 && d[-3] == 0 && d[-2] == 0 && d[-1] == 0)
1203
|| (SIZE_INT(dif) == 4 && d[-2] == 0 && d[-1] == 0) ) {
1204
1205
/* find the number of significant digits */
1206
d = ADDR_INT(dif);
1207
for ( k = SIZE_INT(dif); k != 0; k-- ) {
1208
if ( d[k-1] != 0 )
1209
break;
1210
}
1211
1212
/* reduce to small integer if possible, otherwise shrink bag */
1213
if ( k <= 2 && TNUM_OBJ(dif) == T_INTPOS
1214
&& (UInt)(INTBASE*d[1]+d[0])<(1L<<NR_SMALL_INT_BITS) )
1215
dif = INTOBJ_INT( INTBASE*d[1]+d[0] );
1216
else if ( k <= 2 && TNUM_OBJ(dif) == T_INTNEG
1217
&& (UInt)(INTBASE*d[1]+d[0])<=(1L<<NR_SMALL_INT_BITS) )
1218
dif = INTOBJ_INT( -(Int)(INTBASE*d[1]+d[0]) );
1219
else
1220
ResizeBag( dif, (((k + 3) / 4) * 4) * sizeof(TypDigit) );
1221
}
1222
1223
}
1224
1225
/* subtracting two large integers */
1226
else {
1227
1228
/* if the integers have different sign, let 'SumInt' do the work */
1229
if ( (TNUM_OBJ(opL) == T_INTPOS && TNUM_OBJ(opR) == T_INTNEG)
1230
|| (TNUM_OBJ(opL) == T_INTNEG && TNUM_OBJ(opR) == T_INTPOS) ) {
1231
if ( TNUM_OBJ(opR) == T_INTPOS ) RetypeBag( opR, T_INTNEG );
1232
else RetypeBag( opR, T_INTPOS );
1233
dif = SumInt( opL, opR );
1234
if ( TNUM_OBJ(opR) == T_INTPOS ) RetypeBag( opR, T_INTNEG );
1235
else RetypeBag( opR, T_INTPOS );
1236
return dif;
1237
}
1238
1239
/* make the right operand the smaller one */
1240
if ( SIZE_INT(opL) < SIZE_INT(opR)
1241
|| (TNUM_OBJ(opL) == T_INTPOS && LtInt(opL,opR) )
1242
|| (TNUM_OBJ(opL) == T_INTNEG && LtInt(opR,opL) ) ) {
1243
dif = opL; opL = opR; opR = dif; c = -1;
1244
}
1245
else {
1246
c = 1;
1247
}
1248
1249
/* allocate the result bag and set up the pointers */
1250
if ( (TNUM_OBJ(opL) == T_INTPOS && c == 1)
1251
|| (TNUM_OBJ(opL) == T_INTNEG && c == -1) )
1252
dif = NewBag( T_INTPOS, SIZE_OBJ(opL) );
1253
else
1254
dif = NewBag( T_INTNEG, SIZE_OBJ(opL) );
1255
l = ADDR_INT(opL);
1256
r = ADDR_INT(opR);
1257
d = ADDR_INT(dif);
1258
1259
/* subtract the digits */
1260
c = 0;
1261
for ( k = SIZE_INT(opR)/4; k != 0; k-- ) {
1262
c = (Int)*l++ - (Int)*r++ + (c < 0 ? -1 : 0); *d++ = (TypDigit)c;
1263
c = (Int)*l++ - (Int)*r++ + (c < 0 ? -1 : 0); *d++ = (TypDigit)c;
1264
c = (Int)*l++ - (Int)*r++ + (c < 0 ? -1 : 0); *d++ = (TypDigit)c;
1265
c = (Int)*l++ - (Int)*r++ + (c < 0 ? -1 : 0); *d++ = (TypDigit)c;
1266
}
1267
1268
/* propagate the carry, this loop is almost never executed */
1269
for ( k=(SIZE_INT(opL)-SIZE_INT(opR))/4;
1270
k!=0 && c < 0; k-- ) {
1271
c = (Int)*l++ + (c < 0 ? -1 : 0); *d++ = (TypDigit)c;
1272
c = (Int)*l++ + (c < 0 ? -1 : 0); *d++ = (TypDigit)c;
1273
c = (Int)*l++ + (c < 0 ? -1 : 0); *d++ = (TypDigit)c;
1274
c = (Int)*l++ + (c < 0 ? -1 : 0); *d++ = (TypDigit)c;
1275
}
1276
1277
/* just copy the remaining digits, do it two digits at once */
1278
for ( d2 = (UInt*)d, l2 = (UInt*)l; k != 0; k-- ) {
1279
*d2++ = *l2++;
1280
*d2++ = *l2++;
1281
}
1282
1283
/* no underflow since we subtracted a small int from a large one */
1284
/* but there may be leading zeroes in the result, get rid of them */
1285
/* occurs almost never, so it doesn't matter that it is expensive */
1286
if ( ((UInt*)d == d2
1287
&& d[-4] == 0 && d[-3] == 0 && d[-2] == 0 && d[-1] == 0)
1288
|| (SIZE_INT(dif) == 4 && d[-2] == 0 && d[-1] == 0) ) {
1289
1290
/* find the number of significant digits */
1291
d = ADDR_INT(dif);
1292
for ( k = SIZE_INT(dif); k != 0; k-- ) {
1293
if ( d[k-1] != 0 )
1294
break;
1295
}
1296
1297
/* reduce to small integer if possible, otherwise shrink bag */
1298
if ( k <= 2 && TNUM_OBJ(dif) == T_INTPOS
1299
&& (UInt)(INTBASE*d[1]+d[0]) < (1L<<NR_SMALL_INT_BITS) )
1300
dif = INTOBJ_INT( INTBASE*d[1]+d[0] );
1301
else if ( k <= 2 && TNUM_OBJ(dif) == T_INTNEG
1302
&& (UInt)(INTBASE*d[1]+d[0])<=(1L<<NR_SMALL_INT_BITS))
1303
dif = INTOBJ_INT( -(Int)(INTBASE*d[1]+d[0]) );
1304
else
1305
ResizeBag( dif, (((k + 3) / 4) * 4) * sizeof(TypDigit) );
1306
1307
}
1308
1309
}
1310
1311
/* return the difference */
1312
return dif;
1313
}
1314
1315
1316
/****************************************************************************
1317
**
1318
*F ProdInt( <intL>, <intR> ) . . . . . . . . . . . . product of two integers
1319
**
1320
** 'ProdInt' returns the product of the two integer arguments <intL> and
1321
** <intR>. 'ProdInt' handles operands of type 'T_INT', 'T_INTPOS' and
1322
** 'T_INTNEG'.
1323
**
1324
** It can also be used in the cases that both operands are small integers
1325
** and the result is a small integer too, i.e., that no overflow occurs.
1326
** This case is usually already handled in 'EvalProd' for a better efficiency.
1327
**
1328
** Is called from the 'EvalProd' binop so both operands are already evaluated.
1329
**
1330
** The only difficulty about this function is the fact that is has to handle
1331
** 3 different situations, depending on how many arguments are small ints.
1332
*/
1333
Obj ProdInt (
1334
Obj opL,
1335
Obj opR )
1336
{
1337
Int i; /* loop count, value for small int */
1338
Int k; /* loop count, value for small int */
1339
UInt c; /* product of two digits */
1340
TypDigit l; /* one digit of the left operand */
1341
TypDigit * r; /* pointer into the right operand */
1342
TypDigit * p; /* pointer into the product */
1343
Obj prd; /* handle of the result bag */
1344
1345
/* multiplying two small integers */
1346
if ( ARE_INTOBJS( opL, opR ) ) {
1347
1348
/* multiply two small integers with a small product */
1349
/* multiply and divide back to check that no overflow occured */
1350
if ( PROD_INTOBJS( prd, opL, opR ) ) {
1351
return prd;
1352
}
1353
1354
/* get the integer values */
1355
i = INT_INTOBJ(opL);
1356
k = INT_INTOBJ(opR);
1357
1358
/* allocate the product bag */
1359
if ( (0 < i && 0 < k) || (i < 0 && k < 0) )
1360
prd = NewBag( T_INTPOS, 4*sizeof(TypDigit) );
1361
else
1362
prd = NewBag( T_INTNEG, 4*sizeof(TypDigit) );
1363
p = ADDR_INT(prd);
1364
1365
/* make both operands positive */
1366
if ( i < 0 ) i = -i;
1367
if ( k < 0 ) k = -k;
1368
1369
/* multiply digitwise */
1370
c = (UInt)(TypDigit)i * (TypDigit)k; p[0] = (TypDigit)c;
1371
c = (UInt)(TypDigit)i * (((UInt)k)>>NR_DIGIT_BITS)
1372
+ (c>>NR_DIGIT_BITS); p[1] = (TypDigit)c;
1373
p[2] = c>>NR_DIGIT_BITS;
1374
1375
c = (UInt)(TypDigit)(((UInt)i)>>NR_DIGIT_BITS) * (TypDigit)k
1376
+ p[1]; p[1] = (TypDigit)c;
1377
c = (UInt)(TypDigit)(((UInt)i)>>NR_DIGIT_BITS) * (TypDigit)(((UInt)k)>>NR_DIGIT_BITS)
1378
+ p[2] + (c>>NR_DIGIT_BITS); p[2] = (TypDigit)c;
1379
p[3] = (TypDigit)(c>>NR_DIGIT_BITS);
1380
1381
}
1382
1383
/* multiply a small and a large integer */
1384
else if ( IS_INTOBJ(opL) || IS_INTOBJ(opR) ) {
1385
1386
/* make the left operand the small one */
1387
if ( IS_INTOBJ(opR) ) {
1388
i = INT_INTOBJ(opR); opR = opL;
1389
}
1390
else {
1391
i = INT_INTOBJ(opL);
1392
}
1393
1394
/* handle trivial cases first */
1395
if ( i == 0 )
1396
return INTOBJ_INT(0);
1397
if ( i == 1 )
1398
return opR;
1399
1400
/* the large integer 1<<28 times -1 is the small integer -(1<<28) */
1401
if ( i == -1
1402
&& TNUM_OBJ(opR) == T_INTPOS && SIZE_INT(opR) == 4
1403
&& ADDR_INT(opR)[3] == 0
1404
&& ADDR_INT(opR)[2] == 0
1405
&& ADDR_INT(opR)[1] == (1L<<(NR_SMALL_INT_BITS-NR_DIGIT_BITS))
1406
&& ADDR_INT(opR)[0] == 0 )
1407
return INTOBJ_INT( -(Int)(1L<<NR_SMALL_INT_BITS) );
1408
1409
/* multiplication by -1 is easy, just switch the sign and copy */
1410
if ( i == -1 ) {
1411
if ( TNUM_OBJ(opR) == T_INTPOS )
1412
prd = NewBag( T_INTNEG, SIZE_OBJ(opR) );
1413
else
1414
prd = NewBag( T_INTPOS, SIZE_OBJ(opR) );
1415
r = ADDR_INT(opR);
1416
p = ADDR_INT(prd);
1417
for ( k = SIZE_INT(opR)/4; k != 0; k-- ) {
1418
/*N should be: *p2++=*r2++; *p2++=*r2++; */
1419
*p++ = *r++; *p++ = *r++; *p++ = *r++; *p++ = *r++;
1420
}
1421
return prd;
1422
}
1423
1424
/* allocate a bag for the result */
1425
if ( (0 < i && TNUM_OBJ(opR) == T_INTPOS)
1426
|| (i < 0 && TNUM_OBJ(opR) == T_INTNEG) )
1427
prd = NewBag( T_INTPOS, (SIZE_INT(opR)+4)*sizeof(TypDigit) );
1428
else
1429
prd = NewBag( T_INTNEG, (SIZE_INT(opR)+4)*sizeof(TypDigit) );
1430
if ( i < 0 ) i = -i;
1431
1432
/* multiply with the lower digit of the left operand */
1433
l = (TypDigit)i;
1434
if ( l != 0 ) {
1435
1436
r = ADDR_INT(opR);
1437
p = ADDR_INT(prd);
1438
c = 0;
1439
1440
/* multiply the right with this digit and store in the product */
1441
for ( k = SIZE_INT(opR)/4; k != 0; k-- ) {
1442
c = (UInt)l * (UInt)*r++ + (c>>NR_DIGIT_BITS); *p++ = (TypDigit)c;
1443
c = (UInt)l * (UInt)*r++ + (c>>NR_DIGIT_BITS); *p++ = (TypDigit)c;
1444
c = (UInt)l * (UInt)*r++ + (c>>NR_DIGIT_BITS); *p++ = (TypDigit)c;
1445
c = (UInt)l * (UInt)*r++ + (c>>NR_DIGIT_BITS); *p++ = (TypDigit)c;
1446
}
1447
*p = (TypDigit)(c>>NR_DIGIT_BITS);
1448
}
1449
1450
/* multiply with the larger digit of the left operand */
1451
l = ((UInt)i) >> NR_DIGIT_BITS;
1452
if ( l != 0 ) {
1453
1454
r = ADDR_INT(opR);
1455
p = ADDR_INT(prd) + 1;
1456
c = 0;
1457
1458
/* multiply the right with this digit and add into the product */
1459
for ( k = SIZE_INT(opR)/4; k != 0; k-- ) {
1460
c = (UInt)l * (UInt)*r++ + (UInt)*p + (c>>NR_DIGIT_BITS); *p++ = (TypDigit)c;
1461
c = (UInt)l * (UInt)*r++ + (UInt)*p + (c>>NR_DIGIT_BITS); *p++ = (TypDigit)c;
1462
c = (UInt)l * (UInt)*r++ + (UInt)*p + (c>>NR_DIGIT_BITS); *p++ = (TypDigit)c;
1463
c = (UInt)l * (UInt)*r++ + (UInt)*p + (c>>NR_DIGIT_BITS); *p++ = (TypDigit)c;
1464
}
1465
*p = (TypDigit)(c>>NR_DIGIT_BITS);
1466
}
1467
1468
/* remove the leading zeroes, note that there can't be more than 6 */
1469
p = ADDR_INT(prd) + SIZE_INT(prd);
1470
if ( p[-4] == 0 && p[-3] == 0 && p[-2] == 0 && p[-1] == 0 ) {
1471
ResizeBag( prd, (SIZE_INT(prd)-4)*sizeof(TypDigit) );
1472
}
1473
1474
}
1475
1476
/* multiply two large integers */
1477
else {
1478
1479
/* make the left operand the smaller one, for performance */
1480
if ( SIZE_INT(opL) > SIZE_INT(opR) ) {
1481
prd = opR; opR = opL; opL = prd;
1482
}
1483
1484
/* allocate a bag for the result */
1485
if ( TNUM_OBJ(opL) == TNUM_OBJ(opR) )
1486
prd = NewBag( T_INTPOS, SIZE_OBJ(opL)+SIZE_OBJ(opR) );
1487
else
1488
prd = NewBag( T_INTNEG, SIZE_OBJ(opL)+SIZE_OBJ(opR) );
1489
1490
/* run through the digits of the left operand */
1491
for ( i = 0; i < (Int)SIZE_INT(opL); i++ ) {
1492
1493
/* set up pointer for one loop iteration */
1494
l = ADDR_INT(opL)[i];
1495
if ( l == 0 ) continue;
1496
r = ADDR_INT(opR);
1497
p = ADDR_INT(prd) + i;
1498
c = 0;
1499
1500
/* multiply the right with this digit and add into the product */
1501
for ( k = SIZE_INT(opR)/4; k != 0; k-- ) {
1502
c = (UInt)l * (UInt)*r++ + (UInt)*p + (c>>NR_DIGIT_BITS); *p++ = (TypDigit)c;
1503
c = (UInt)l * (UInt)*r++ + (UInt)*p + (c>>NR_DIGIT_BITS); *p++ = (TypDigit)c;
1504
c = (UInt)l * (UInt)*r++ + (UInt)*p + (c>>NR_DIGIT_BITS); *p++ = (TypDigit)c;
1505
c = (UInt)l * (UInt)*r++ + (UInt)*p + (c>>NR_DIGIT_BITS); *p++ = (TypDigit)c;
1506
}
1507
*p = (TypDigit)(c>>NR_DIGIT_BITS);
1508
}
1509
1510
/* remove the leading zeroes, note that there can't be more than 7 */
1511
p = ADDR_INT(prd) + SIZE_INT(prd);
1512
if ( p[-4] == 0 && p[-3] == 0 && p[-2] == 0 && p[-1] == 0 ) {
1513
ResizeBag( prd, (SIZE_INT(prd)-4)*sizeof(TypDigit) );
1514
}
1515
1516
}
1517
1518
/* return the product */
1519
return prd;
1520
}
1521
1522
1523
/****************************************************************************
1524
**
1525
*F ProdIntObj(<n>,<op>) . . . . . . . . product of an integer and an object
1526
*/
1527
Obj ProdIntObj (
1528
Obj n,
1529
Obj op )
1530
{
1531
Obj res = 0; /* result */
1532
UInt i, k, l; /* loop variables */
1533
1534
/* if the integer is zero, return the neutral element of the operand */
1535
if ( TNUM_OBJ(n) == T_INT && INT_INTOBJ(n) == 0 ) {
1536
res = ZERO( op );
1537
}
1538
1539
/* if the integer is one, return the object if immutable
1540
if mutable, add the object to its ZeroSameMutability to
1541
ensure correct mutability propagation */
1542
else if ( TNUM_OBJ(n) == T_INT && INT_INTOBJ(n) == 1 ) {
1543
if (IS_MUTABLE_OBJ(op))
1544
res = SUM(ZERO(op),op);
1545
else
1546
res = op;
1547
}
1548
1549
/* if the integer is minus one, return the inverse of the operand */
1550
else if ( TNUM_OBJ(n) == T_INT && INT_INTOBJ(n) == -1 ) {
1551
res = AINV( op );
1552
}
1553
1554
/* if the integer is negative, invert the operand and the integer */
1555
else if ( TNUM_OBJ(n) == T_INT && INT_INTOBJ(n) < -1 ) {
1556
res = AINV( op );
1557
if ( res == Fail ) {
1558
return ErrorReturnObj(
1559
"Operations: <obj> must have an additive inverse",
1560
0L, 0L,
1561
"you can supply an inverse <inv> for <obj> via 'return <inv>;'" );
1562
}
1563
res = PROD( AINV( n ), res );
1564
}
1565
1566
/* if the integer is negative, invert the operand and the integer */
1567
else if ( TNUM_OBJ(n) == T_INTNEG ) {
1568
res = AINV( op );
1569
if ( res == Fail ) {
1570
return ErrorReturnObj(
1571
"Operations: <obj> must have an additive inverse",
1572
0L, 0L,
1573
"you can supply an inverse <inv> for <obj> via 'return <inv>;'" );
1574
}
1575
res = PROD( AINV( n ), res );
1576
}
1577
1578
/* if the integer is small, compute the product by repeated doubling */
1579
/* the loop invariant is <result> = <k>*<res> + <l>*<op>, <l> < <k> */
1580
/* <res> = 0 means that <res> is the neutral element */
1581
else if ( TNUM_OBJ(n) == T_INT && INT_INTOBJ(n) > 1 ) {
1582
res = 0;
1583
k = 1L << (NR_SMALL_INT_BITS+1);
1584
l = INT_INTOBJ(n);
1585
while ( 1 < k ) {
1586
res = (res == 0 ? res : SUM( res, res ));
1587
k = k / 2;
1588
if ( k <= l ) {
1589
res = (res == 0 ? op : SUM( res, op ));
1590
l = l - k;
1591
}
1592
}
1593
}
1594
1595
/* if the integer is large, compute the product by repeated doubling */
1596
else if ( TNUM_OBJ(n) == T_INTPOS ) {
1597
res = 0;
1598
for ( i = SIZE_OBJ(n)/sizeof(TypDigit); 0 < i; i-- ) {
1599
k = 1L << (8*sizeof(TypDigit));
1600
l = ((TypDigit*) ADDR_OBJ(n))[i-1];
1601
while ( 1 < k ) {
1602
res = (res == 0 ? res : SUM( res, res ));
1603
k = k / 2;
1604
if ( k <= l ) {
1605
res = (res == 0 ? op : SUM( res, op ));
1606
l = l - k;
1607
}
1608
}
1609
}
1610
}
1611
1612
/* return the result */
1613
return res;
1614
}
1615
1616
Obj ProdIntObjFunc;
1617
1618
Obj FuncPROD_INT_OBJ (
1619
Obj self,
1620
Obj opL,
1621
Obj opR )
1622
{
1623
return ProdIntObj( opL, opR );
1624
}
1625
1626
1627
/****************************************************************************
1628
**
1629
*F OneInt(<int>) . . . . . . . . . . . . . . . . . . . . . one of an integer
1630
*/
1631
Obj OneInt (
1632
Obj op )
1633
{
1634
return INTOBJ_INT( 1L );
1635
}
1636
1637
1638
/****************************************************************************
1639
**
1640
*F PowInt( <intL>, <intR> ) . . . . . . . . . . . . . . power of an integer
1641
**
1642
** 'PowInt' returns the <intR>-th (an integer) power of the integer <intL>.
1643
** 'PowInt' handles operands of type 'T_INT', 'T_INTPOS' and 'T_INTNEG'.
1644
**
1645
** It can also be used in the cases that both operands are small integers
1646
** and the result is a small integer too, i.e., that no overflow occurs.
1647
** This case is usually already handled in 'EvalPow' for a better efficiency.
1648
**
1649
** Is called from the 'EvalPow' binop so both operands are already evaluated.
1650
*/
1651
Obj PowInt (
1652
Obj opL,
1653
Obj opR )
1654
{
1655
Int i;
1656
Obj pow;
1657
1658
/* power with a large exponent */
1659
ReadGuard(opL);
1660
ReadGuard(opR);
1661
if ( ! IS_INTOBJ(opR) ) {
1662
if ( opL == INTOBJ_INT(0) )
1663
pow = INTOBJ_INT(0);
1664
else if ( opL == INTOBJ_INT(1) )
1665
pow = INTOBJ_INT(1);
1666
else if ( opL == INTOBJ_INT(-1) && ADDR_INT(opR)[0] % 2 == 0 )
1667
pow = INTOBJ_INT(1);
1668
else if ( opL == INTOBJ_INT(-1) && ADDR_INT(opR)[0] % 2 != 0 )
1669
pow = INTOBJ_INT(-1);
1670
else {
1671
opR = ErrorReturnObj(
1672
"Integer operands: <exponent> is too large",
1673
0L, 0L,
1674
"you can replace the integer <exponent> via 'return <exponent>;'" );
1675
return POW( opL, opR );
1676
}
1677
}
1678
1679
/* power with a negative exponent */
1680
else if ( INT_INTOBJ(opR) < 0 ) {
1681
if ( opL == INTOBJ_INT(0) ) {
1682
opL = ErrorReturnObj(
1683
"Integer operands: <base> must not be zero",
1684
0L, 0L,
1685
"you can replace the integer <base> via 'return <base>;'" );
1686
return POW( opL, opR );
1687
}
1688
else if ( opL == INTOBJ_INT(1) )
1689
pow = INTOBJ_INT(1);
1690
else if ( opL == INTOBJ_INT(-1) && INT_INTOBJ(opR) % 2 == 0 )
1691
pow = INTOBJ_INT(1);
1692
else if ( opL == INTOBJ_INT(-1) && INT_INTOBJ(opR) % 2 != 0 )
1693
pow = INTOBJ_INT(-1);
1694
else
1695
pow = QUO( INTOBJ_INT(1),
1696
PowInt( opL, INTOBJ_INT( -INT_INTOBJ(opR)) ) );
1697
}
1698
1699
/* power with a small positive exponent, do it by a repeated squaring */
1700
else {
1701
pow = INTOBJ_INT(1);
1702
i = INT_INTOBJ(opR);
1703
while ( i != 0 ) {
1704
if ( i % 2 == 1 ) pow = ProdInt( pow, opL );
1705
if ( i > 1 ) opL = ProdInt( opL, opL );
1706
i = i / 2;
1707
}
1708
}
1709
1710
/* return the power */
1711
return pow;
1712
}
1713
1714
1715
/****************************************************************************
1716
**
1717
*F PowObjInt(<op>,<n>) . . . . . . . . . . power of an object and an integer
1718
*/
1719
static Obj OneAttr;
1720
1721
Obj PowObjInt (
1722
Obj op,
1723
Obj n )
1724
{
1725
Obj res = 0; /* result */
1726
UInt i, k, l; /* loop variables */
1727
1728
/* if the integer is zero, return the neutral element of the operand */
1729
if ( TNUM_OBJ(n) == T_INT && INT_INTOBJ(n) == 0 ) {
1730
return ONE_MUT( op );
1731
}
1732
1733
/* if the integer is one, return a copy of the operand */
1734
else if ( TNUM_OBJ(n) == T_INT && INT_INTOBJ(n) == 1 ) {
1735
res = CopyObj( op, 1 );
1736
}
1737
1738
/* if the integer is minus one, return the inverse of the operand */
1739
else if ( TNUM_OBJ(n) == T_INT && INT_INTOBJ(n) == -1 ) {
1740
res = INV_MUT( op );
1741
}
1742
1743
/* if the integer is negative, invert the operand and the integer */
1744
else if ( TNUM_OBJ(n) == T_INT && INT_INTOBJ(n) < 0 ) {
1745
res = INV_MUT( op );
1746
if ( res == Fail ) {
1747
return ErrorReturnObj(
1748
"Operations: <obj> must have an inverse",
1749
0L, 0L,
1750
"you can supply an inverse <inv> for <obj> via 'return <inv>;'" );
1751
}
1752
res = POW( res, AINV( n ) );
1753
}
1754
1755
/* if the integer is negative, invert the operand and the integer */
1756
else if ( TNUM_OBJ(n) == T_INTNEG ) {
1757
res = INV_MUT( op );
1758
if ( res == Fail ) {
1759
return ErrorReturnObj(
1760
"Operations: <obj> must have an inverse",
1761
0L, 0L,
1762
"you can supply an inverse <inv> for <obj> via 'return <inv>;'" );
1763
}
1764
res = POW( res, AINV( n ) );
1765
}
1766
1767
/* if the integer is small, compute the power by repeated squaring */
1768
/* the loop invariant is <result> = <res>^<k> * <op>^<l>, <l> < <k> */
1769
/* <res> = 0 means that <res> is the neutral element */
1770
else if ( TNUM_OBJ(n) == T_INT && INT_INTOBJ(n) > 0 ) {
1771
res = 0;
1772
k = 1L << (NR_SMALL_INT_BITS+1);
1773
l = INT_INTOBJ(n);
1774
while ( 1 < k ) {
1775
res = (res == 0 ? res : PROD( res, res ));
1776
k = k / 2;
1777
if ( k <= l ) {
1778
res = (res == 0 ? op : PROD( res, op ));
1779
l = l - k;
1780
}
1781
}
1782
}
1783
1784
/* if the integer is large, compute the power by repeated squaring */
1785
else if ( TNUM_OBJ(n) == T_INTPOS ) {
1786
res = 0;
1787
for ( i = SIZE_OBJ(n)/sizeof(TypDigit); 0 < i; i-- ) {
1788
k = 1L << (8*sizeof(TypDigit));
1789
l = ((TypDigit*) ADDR_OBJ(n))[i-1];
1790
while ( 1 < k ) {
1791
res = (res == 0 ? res : PROD( res, res ));
1792
k = k / 2;
1793
if ( k <= l ) {
1794
res = (res == 0 ? op : PROD( res, op ));
1795
l = l - k;
1796
}
1797
}
1798
}
1799
}
1800
1801
/* return the result */
1802
return res;
1803
}
1804
1805
Obj PowObjIntFunc;
1806
1807
Obj FuncPOW_OBJ_INT (
1808
Obj self,
1809
Obj opL,
1810
Obj opR )
1811
{
1812
return PowObjInt( opL, opR );
1813
}
1814
1815
1816
/****************************************************************************
1817
**
1818
*F ModInt( <intL>, <intR> ) . representative of residue class of an integer
1819
**
1820
** 'ModInt' returns the smallest positive representant of the residue class
1821
** of the integer <intL> modulo the integer <intR>. 'ModInt' handles
1822
** operands of type 'T_INT', 'T_INTPOS', 'T_INTNEG'.
1823
**
1824
** It can also be used in the cases that both operands are small integers
1825
** and the result is a small integer too, i.e., that no overflow occurs.
1826
** This case is usually already handled in 'EvalMod' for a better efficiency.
1827
p**
1828
** Is called from the 'EvalMod' binop so both operands are already evaluated.
1829
*/
1830
Obj ModInt (
1831
Obj opL,
1832
Obj opR )
1833
{
1834
Int i; /* loop count, value for small int */
1835
Int k; /* loop count, value for small int */
1836
UInt c; /* product of two digits */
1837
TypDigit d; /* carry into the next digit */
1838
TypDigit * l; /* pointer into the left operand */
1839
TypDigit * r; /* pointer into the right operand */
1840
TypDigit r1; /* leading digit of the right oper */
1841
TypDigit r2; /* next digit of the right operand */
1842
UInt rs; /* size of the right operand */
1843
UInt e; /* we mult r by 2^e so r1 >= 32768 */
1844
Obj mod; /* handle of the remainder bag */
1845
TypDigit * m; /* pointer into the remainder */
1846
UInt m01; /* leading two digits of the rem. */
1847
TypDigit m2; /* next digit of the remainder */
1848
TypDigit qi; /* guessed digit of the quotient */
1849
1850
/* compute the remainder of two small integers */
1851
if ( ARE_INTOBJS( opL, opR ) ) {
1852
1853
/* pathological case first */
1854
if ( opR == INTOBJ_INT(0) ) {
1855
opR = ErrorReturnObj(
1856
"Integer operations: <divisor> must be nonzero",
1857
0L, 0L,
1858
"you can replace the integer <divisor> via 'return <divisor>;'" );
1859
return MOD( opL, opR );
1860
}
1861
1862
/* get the integer values */
1863
i = INT_INTOBJ(opL);
1864
k = INT_INTOBJ(opR);
1865
1866
/* compute the remainder, make sure we divide only positive numbers*/
1867
if ( 0 <= i && 0 <= k ) i = ( i % k );
1868
else if ( 0 <= i && k < 0 ) i = ( i % -k );
1869
else if ( i < 0 && 0 <= k ) i = ( k - ( -i % k )) % k;
1870
else if ( i < 0 && k < 0 ) i = (-k - ( -i % -k )) % k;
1871
mod = INTOBJ_INT( i );
1872
1873
}
1874
1875
/* compute the remainder of a small integer by a large integer */
1876
else if ( IS_INTOBJ(opL) ) {
1877
1878
/* the small int -(1<<28) mod the large int (1<<28) is 0 */
1879
if ( opL == INTOBJ_INT((UInt)-(Int)(1L<<NR_SMALL_INT_BITS))
1880
&& TNUM_OBJ(opR) == T_INTPOS && SIZE_INT(opR) == 4
1881
&& ADDR_INT(opR)[3] == 0
1882
&& ADDR_INT(opR)[2] == 0
1883
&& ADDR_INT(opR)[1] == (NR_SMALL_INT_BITS-NR_DIGIT_BITS)
1884
&& ADDR_INT(opR)[0] == 0 )
1885
mod = INTOBJ_INT(0);
1886
1887
/* in all other cases the remainder is equal the left operand */
1888
else if ( 0 <= INT_INTOBJ(opL) )
1889
mod = opL;
1890
else if ( TNUM_OBJ(opR) == T_INTPOS )
1891
mod = SumInt( opL, opR );
1892
else
1893
mod = DiffInt( opL, opR );
1894
1895
}
1896
1897
/* compute the remainder of a large integer by a small integer */
1898
else if ( IS_INTOBJ(opR)
1899
&& INT_INTOBJ(opR) < INTBASE
1900
&& -(Int)INTBASE <= INT_INTOBJ(opR) ) {
1901
1902
/* pathological case first */
1903
if ( opR == INTOBJ_INT(0) ) {
1904
opR = ErrorReturnObj(
1905
"Integer operations: <divisor> must be nonzero",
1906
0L, 0L,
1907
"you can replace the integer <divisor> via 'return <divisor>;'" );
1908
return MOD( opL, opR );
1909
}
1910
1911
/* get the integer value, make positive */
1912
i = INT_INTOBJ(opR); if ( i < 0 ) i = -i;
1913
1914
/* maybe its trivial */
1915
if ( INTBASE % i == 0 ) {
1916
c = ADDR_INT(opL)[0] % i;
1917
}
1918
1919
/* otherwise run through the left operand and divide digitwise */
1920
else {
1921
l = ADDR_INT(opL) + SIZE_INT(opL) - 1;
1922
c = 0;
1923
for ( ; l >= ADDR_INT(opL); l-- ) {
1924
c = (c<<NR_DIGIT_BITS) + (Int)*l;
1925
c = c % i;
1926
}
1927
}
1928
1929
/* now c is the result, it has the same sign as the left operand */
1930
if ( TNUM_OBJ(opL) == T_INTPOS )
1931
mod = INTOBJ_INT( c );
1932
else if ( c == 0 )
1933
mod = INTOBJ_INT( c );
1934
else if ( 0 <= INT_INTOBJ(opR) )
1935
mod = SumInt( INTOBJ_INT( -(Int)c ), opR );
1936
else
1937
mod = DiffInt( INTOBJ_INT( -(Int)c ), opR );
1938
1939
}
1940
1941
/* compute the remainder of a large integer modulo a large integer */
1942
else {
1943
/* a small divisor larger than one digit isn't handled above */
1944
if ( IS_INTOBJ(opR) ) {
1945
if ( 0 < INT_INTOBJ(opR) ) {
1946
mod = NewBag( T_INTPOS, 4*sizeof(TypDigit) );
1947
ADDR_INT(mod)[0] = (TypDigit)(INT_INTOBJ(opR));
1948
ADDR_INT(mod)[1] = (TypDigit)(INT_INTOBJ(opR)>>NR_DIGIT_BITS);
1949
opR = mod;
1950
}
1951
else {
1952
mod = NewBag( T_INTNEG, 4*sizeof(TypDigit) );
1953
ADDR_INT(mod)[0] = (TypDigit)(-INT_INTOBJ(opR));
1954
ADDR_INT(mod)[1] = (TypDigit)((-INT_INTOBJ(opR))>>NR_DIGIT_BITS);
1955
opR = mod;
1956
}
1957
}
1958
1959
/* trivial case first */
1960
if ( SIZE_INT(opL) < SIZE_INT(opR) ) {
1961
if ( TNUM_OBJ(opL) == T_INTPOS )
1962
return opL;
1963
else if ( TNUM_OBJ(opR) == T_INTPOS )
1964
return SumInt( opL, opR );
1965
else
1966
return DiffInt( opL, opR );
1967
}
1968
1969
/* copy the left operand into a new bag, this holds the remainder */
1970
mod = NewBag( TNUM_OBJ(opL), (SIZE_INT(opL)+4)*sizeof(TypDigit) );
1971
l = ADDR_INT(opL);
1972
m = ADDR_INT(mod);
1973
for ( k = SIZE_INT(opL)-1; k >= 0; k-- )
1974
*m++ = *l++;
1975
1976
/* get the size of the right operand, and get the leading 2 digits */
1977
rs = SIZE_INT(opR);
1978
r = ADDR_INT(opR);
1979
while ( r[rs-1] == 0 ) rs--;
1980
for ( e = 0;
1981
((Int)r[rs-1]<<e) + (e ? r[rs-2]>>(NR_DIGIT_BITS-e) : 0)
1982
< INTBASE/2; e++ ) ;
1983
1984
r1 = ((Int)r[rs-1]<<e) + (e ? r[rs-2]>>(NR_DIGIT_BITS-e) : 0);
1985
r2 = ((Int)r[rs-2]<<e) + ((rs>=3 && e)? r[rs-3]>>(NR_DIGIT_BITS-e) : 0);
1986
1987
/* run through the digits in the quotient */
1988
for ( i = SIZE_INT(mod)-SIZE_INT(opR)-1; i >= 0; i-- ) {
1989
1990
/* guess the factor */
1991
m = ADDR_INT(mod) + rs + i;
1992
m01 = ((INTBASE*m[0]+m[-1])<<e)
1993
+ (e ? m[-2]>>(NR_DIGIT_BITS-e) : 0);
1994
if ( m01 == 0 ) continue;
1995
m2 = ((Int)m[-2]<<e) + ((e && rs+i>=3) ? m[-3]>>(NR_DIGIT_BITS-e) : 0);
1996
if ( ((Int)m[0]<<e)+(e ? m[-1]>>(NR_DIGIT_BITS-e) : 0) < r1 )
1997
qi = m01 / r1;
1998
else qi = INTBASE - 1;
1999
while ( m01-(Int)qi*r1 < (Int)INTBASE
2000
&& INTBASE*(m01-(Int)qi*r1)+m2 < (Int)qi*r2 )
2001
qi--;
2002
2003
/* m = m - qi * r; */
2004
d = 0;
2005
m = ADDR_INT(mod) + i;
2006
r = ADDR_INT(opR);
2007
for ( k = 0; k < (Int)rs; ++k, ++m, ++r ) {
2008
c = (Int)*m - (Int)qi * *r - (Int)d;
2009
*m = (TypDigit)c;
2010
d = -(TypDigit)(c>>NR_DIGIT_BITS);
2011
}
2012
c = (Int)*m - d; *m = (TypDigit)c; d = -(TypDigit)(c>>NR_DIGIT_BITS);
2013
2014
/* if we have a borrow then add back */
2015
if ( d != 0 ) {
2016
d = 0;
2017
m = ADDR_INT(mod) + i;
2018
r = ADDR_INT(opR);
2019
for ( k = 0; k < (Int)rs; ++k, ++m, ++r ) {
2020
c = (Int)*m + (Int)*r + (Int)d;
2021
*m = (TypDigit)c;
2022
d = (TypDigit)(c>>NR_DIGIT_BITS);
2023
}
2024
c = (Int)*m + d; *m = (TypDigit)c; d = (TypDigit)(c>>NR_DIGIT_BITS);
2025
qi--;
2026
}
2027
2028
}
2029
2030
/* remove the leading zeroes */
2031
m = ADDR_INT(mod) + SIZE_INT(mod);
2032
if ( (m[-4] == 0 && m[-3] == 0 && m[-2] == 0 && m[-1] == 0)
2033
|| (SIZE_INT(mod) == 4 && m[-2] == 0 && m[-1] == 0) ) {
2034
2035
/* find the number of significant digits */
2036
m = ADDR_INT(mod);
2037
for ( k = SIZE_INT(mod); k != 0; k-- ) {
2038
if ( m[k-1] != 0 )
2039
break;
2040
}
2041
2042
/* reduce to small integer if possible, otherwise shrink bag */
2043
2044
if ( k <= 2 && TNUM_OBJ(mod) == T_INTPOS
2045
&& (UInt)(INTBASE*m[1]+m[0])<(1L<<NR_SMALL_INT_BITS) )
2046
mod = INTOBJ_INT( INTBASE*m[1]+m[0] );
2047
else if ( k <= 2 && TNUM_OBJ(mod) == T_INTNEG
2048
&& (UInt)(INTBASE*m[1]+m[0])<=(1L<<NR_SMALL_INT_BITS) )
2049
mod = INTOBJ_INT( -(Int)(INTBASE*m[1]+m[0]) );
2050
else
2051
ResizeBag( mod, (((k + 3) / 4) * 4) * sizeof(TypDigit) );
2052
}
2053
2054
/* make the representative positive */
2055
if ( (TNUM_OBJ(mod) == T_INT && INT_INTOBJ(mod) < 0)
2056
|| TNUM_OBJ(mod) == T_INTNEG ) {
2057
if ( TNUM_OBJ(opR) == T_INTPOS )
2058
mod = SumInt( mod, opR );
2059
else
2060
mod = DiffInt( mod, opR );
2061
}
2062
2063
}
2064
2065
/* return the result */
2066
return mod;
2067
}
2068
2069
2070
/****************************************************************************
2071
**
2072
*F FuncIS_INT( <self>, <val> ) . . . . . . . . . . internal function 'IsInt'
2073
**
2074
** 'FuncIS_INT' implements the internal filter 'IsInt'.
2075
**
2076
** 'IsInt( <val> )'
2077
**
2078
** 'IsInt' returns 'true' if the value <val> is an integer and 'false'
2079
** otherwise.
2080
*/
2081
Obj IsIntFilt;
2082
2083
Obj FuncIS_INT (
2084
Obj self,
2085
Obj val )
2086
{
2087
/* return 'true' if <obj> is an integer and 'false' otherwise */
2088
if ( TNUM_OBJ(val) == T_INT
2089
|| TNUM_OBJ(val) == T_INTPOS
2090
|| TNUM_OBJ(val) == T_INTNEG ) {
2091
return True;
2092
}
2093
else if ( TNUM_OBJ(val) < FIRST_EXTERNAL_TNUM ) {
2094
return False;
2095
}
2096
else {
2097
return DoFilter( self, val );
2098
}
2099
}
2100
2101
2102
/****************************************************************************
2103
**
2104
*F QuoInt( <intL>, <intR> ) . . . . . . . . . . . quotient of two integers
2105
**
2106
** 'QuoInt' returns the integer part of the two integers <intL> and <intR>.
2107
** 'QuoInt' handles operands of type 'T_INT', 'T_INTPOS' and 'T_INTNEG'.
2108
**
2109
** It can also be used in the cases that both operands are small integers
2110
** and the result is a small integer too, i.e., that no overflow occurs.
2111
**
2112
** Note that this routine is not called from 'EvalQuo', the division of two
2113
** integers yields a rational and is therefor performed in 'QuoRat'.
2114
** This operation is however available through the internal function 'Quo'.
2115
*/
2116
Obj QuoInt (
2117
Obj opL,
2118
Obj opR )
2119
{
2120
Int i; /* loop count, value for small int */
2121
Int k; /* loop count, value for small int */
2122
UInt c; /* product of two digits */
2123
TypDigit d; /* carry into the next digit */
2124
TypDigit * l; /* pointer into the left operand */
2125
UInt l01; /* leading two digits of the left */
2126
TypDigit l2; /* next digit of the left operand */
2127
TypDigit * r; /* pointer into the right operand */
2128
TypDigit r1; /* leading digit of the right oper */
2129
TypDigit r2; /* next digit of the right operand */
2130
UInt rs; /* size of the right operand */
2131
UInt e; /* we mult r by 2^e so r1 >= 32768 */
2132
TypDigit * q; /* pointer into the quotient */
2133
Obj quo; /* handle of the result bag */
2134
TypDigit qi; /* guessed digit of the quotient */
2135
2136
/* divide to small integers */
2137
if ( ARE_INTOBJS( opL, opR ) ) {
2138
2139
/* pathological case first */
2140
if ( opR == INTOBJ_INT(0) ) {
2141
opR = ErrorReturnObj(
2142
"Integer operations: <divisor> must be nonzero",
2143
0L, 0L,
2144
"you can replace the integer <divisor> via 'return <divisor>;'" );
2145
return QUO( opL, opR );
2146
}
2147
2148
/* the small int -(1<<28) divided by -1 is the large int (1<<28) */
2149
if ( opL == INTOBJ_INT(-(Int)(1L<<NR_SMALL_INT_BITS))
2150
&& opR == INTOBJ_INT(-1) ) {
2151
quo = NewBag( T_INTPOS, 4*sizeof(TypDigit) );
2152
ADDR_INT(quo)[1] = 1L<<(NR_SMALL_INT_BITS-NR_DIGIT_BITS);
2153
ADDR_INT(quo)[0] = 0;
2154
return quo;
2155
}
2156
2157
/* get the integer values */
2158
i = INT_INTOBJ(opL);
2159
k = INT_INTOBJ(opR);
2160
2161
/* divide, make sure we divide only positive numbers */
2162
if ( 0 <= i && 0 <= k ) i = ( i / k );
2163
else if ( 0 <= i && k < 0 ) i = - ( i / -k );
2164
else if ( i < 0 && 0 <= k ) i = - ( -i / k );
2165
else if ( i < 0 && k < 0 ) i = ( -i / -k );
2166
quo = INTOBJ_INT( i );
2167
2168
}
2169
2170
/* divide a small integer by a large one */
2171
else if ( IS_INTOBJ(opL) ) {
2172
2173
/* the small int -(1<<28) divided by the large int (1<<28) is -1 */
2174
2175
if ( opL == INTOBJ_INT(-(Int)(1L<<NR_SMALL_INT_BITS))
2176
&& TNUM_OBJ(opR) == T_INTPOS && SIZE_INT(opR) == 4
2177
&& ADDR_INT(opR)[3] == 0
2178
&& ADDR_INT(opR)[2] == 0
2179
&& ADDR_INT(opR)[1] == 1L<<(NR_SMALL_INT_BITS-NR_DIGIT_BITS)
2180
&& ADDR_INT(opR)[0] == 0 )
2181
quo = INTOBJ_INT(-1);
2182
2183
/* in all other cases the quotient is of course zero */
2184
else
2185
quo = INTOBJ_INT(0);
2186
2187
}
2188
2189
/* divide a large integer by a small integer */
2190
else if ( IS_INTOBJ(opR)
2191
&& INT_INTOBJ(opR) < INTBASE
2192
&& -(Int)INTBASE <= INT_INTOBJ(opR) ) {
2193
2194
/* pathological case first */
2195
if ( opR == INTOBJ_INT(0) ) {
2196
opR = ErrorReturnObj(
2197
"Integer operations: <divisor> must be nonzero",
2198
0L, 0L,
2199
"you can replace the integer <divisor> via 'return <divisor>;'" );
2200
return QUO( opL, opR );
2201
}
2202
2203
/* get the integer value, make positive */
2204
i = INT_INTOBJ(opR); if ( i < 0 ) i = -i;
2205
2206
/* allocate a bag for the result and set up the pointers */
2207
if ( (TNUM_OBJ(opL)==T_INTPOS && 0 < INT_INTOBJ(opR))
2208
|| (TNUM_OBJ(opL)==T_INTNEG && INT_INTOBJ(opR) < 0) )
2209
quo = NewBag( T_INTPOS, SIZE_OBJ(opL) );
2210
else
2211
quo = NewBag( T_INTNEG, SIZE_OBJ(opL) );
2212
l = ADDR_INT(opL) + SIZE_INT(opL) - 1;
2213
q = ADDR_INT(quo) + SIZE_INT(quo) - 1;
2214
2215
/* run through the left operand and divide digitwise */
2216
c = 0;
2217
for ( ; l >= ADDR_INT(opL); l--, q-- ) {
2218
c = (c<<NR_DIGIT_BITS) + (Int)*l;
2219
*q = (TypDigit)(c / i);
2220
c = c - i * *q;
2221
/*N clever compilers may prefer: c = c % i; */
2222
}
2223
2224
/* remove the leading zeroes, note that there can't be more than 5 */
2225
q = ADDR_INT(quo) + SIZE_INT(quo);
2226
if ( q[-4] == 0 && q[-3] == 0 && q[-2] == 0 && q[-1] == 0 ) {
2227
ResizeBag( quo, (SIZE_INT(quo)-4)*sizeof(TypDigit) );
2228
}
2229
2230
/* reduce to small integer if possible */
2231
q = ADDR_INT(quo) + SIZE_INT(quo);
2232
if ( SIZE_INT(quo) == 4 && q[-2] == 0 && q[-1] == 0 ) {
2233
if ( TNUM_OBJ(quo) == T_INTPOS
2234
&&(UInt)(INTBASE*q[-3]+q[-4])
2235
< (1L<<NR_SMALL_INT_BITS))
2236
quo = INTOBJ_INT( INTBASE*q[-3]+q[-4] );
2237
else if ( TNUM_OBJ(quo) == T_INTNEG
2238
&& (UInt)(INTBASE*q[-3]+q[-4])
2239
<= (1L<<NR_SMALL_INT_BITS) )
2240
quo = INTOBJ_INT( -(Int)(INTBASE*q[-3]+q[-4]) );
2241
}
2242
2243
}
2244
2245
/* divide a large integer by a large integer */
2246
else {
2247
2248
/* a small divisor larger than one digit isn't handled above */
2249
if ( IS_INTOBJ(opR) ) {
2250
if ( 0 < INT_INTOBJ(opR) ) {
2251
quo = NewBag( T_INTPOS, 4*sizeof(TypDigit) );
2252
ADDR_INT(quo)[0] = (TypDigit)(INT_INTOBJ(opR));
2253
ADDR_INT(quo)[1] = (TypDigit)(INT_INTOBJ(opR)>>NR_DIGIT_BITS);
2254
opR = quo;
2255
}
2256
else {
2257
quo = NewBag( T_INTNEG, 4*sizeof(TypDigit) );
2258
ADDR_INT(quo)[0] = (TypDigit)(-INT_INTOBJ(opR));
2259
ADDR_INT(quo)[1] = (TypDigit)((-INT_INTOBJ(opR))>>NR_DIGIT_BITS);
2260
opR = quo;
2261
}
2262
}
2263
2264
/* trivial case first */
2265
if ( SIZE_INT(opL) < SIZE_INT(opR) )
2266
return INTOBJ_INT(0);
2267
2268
/* copy the left operand into a new bag, this holds the remainder */
2269
quo = NewBag( TNUM_OBJ(opL), (SIZE_INT(opL)+4)*sizeof(TypDigit) );
2270
l = ADDR_INT(opL);
2271
q = ADDR_INT(quo);
2272
for ( k = SIZE_INT(opL)-1; k >= 0; k-- )
2273
*q++ = *l++;
2274
opL = quo;
2275
2276
/* get the size of the right operand, and get the leading 2 digits */
2277
rs = SIZE_INT(opR);
2278
r = ADDR_INT(opR);
2279
while ( r[rs-1] == 0 ) rs--;
2280
for ( e = 0;
2281
((Int)r[rs-1]<<e) + (e ? r[rs-2]>>(NR_DIGIT_BITS-e) : 0)
2282
< INTBASE/2 ; e++ ) ;
2283
2284
r1 = ((Int)r[rs-1]<<e) + (e ? r[rs-2]>>(NR_DIGIT_BITS-e) : 0);
2285
r2 = ((Int)r[rs-2]<<e) + ((e && rs>=3) ? r[rs-3]>>(NR_DIGIT_BITS-e) : 0);
2286
2287
/* allocate a bag for the quotient */
2288
if ( TNUM_OBJ(opL) == TNUM_OBJ(opR) )
2289
quo = NewBag( T_INTPOS, SIZE_OBJ(opL)-SIZE_OBJ(opR) );
2290
else
2291
quo = NewBag( T_INTNEG, SIZE_OBJ(opL)-SIZE_OBJ(opR) );
2292
2293
/* run through the digits in the quotient */
2294
for ( i = SIZE_INT(opL)-SIZE_INT(opR)-1; i >= 0; i-- ) {
2295
2296
/* guess the factor */
2297
l = ADDR_INT(opL) + rs + i;
2298
l01 = ((INTBASE*l[0]+l[-1])<<e)
2299
+ (e ? l[-2]>>(NR_DIGIT_BITS-e):0);
2300
2301
if ( l01 == 0 ) continue;
2302
l2 = ((Int)l[-2]<<e) + (( e && rs+i>=3) ? l[-3]>>(NR_DIGIT_BITS-e) : 0);
2303
if ( ((Int)l[0]<<e)+(e ? l[-1]>>(NR_DIGIT_BITS-e):0) < r1 )
2304
qi = l01 / r1;
2305
else qi = INTBASE - 1;
2306
while ( l01-(Int)qi*r1 < (Int)INTBASE
2307
&& INTBASE*(l01-(UInt)qi*r1)+l2 < (UInt)qi*r2 )
2308
qi--;
2309
2310
/* l = l - qi * r; */
2311
d = 0;
2312
l = ADDR_INT(opL) + i;
2313
r = ADDR_INT(opR);
2314
for ( k = 0; k < (Int)rs; ++k, ++l, ++r ) {
2315
c = *l - (Int)qi * *r - d; *l = c; d = -(TypDigit)(c>>NR_DIGIT_BITS);
2316
}
2317
c = (Int)*l - d; d = -(TypDigit)(c>>NR_DIGIT_BITS);
2318
2319
/* if we have a borrow then add back */
2320
if ( d != 0 ) {
2321
d = 0;
2322
l = ADDR_INT(opL) + i;
2323
r = ADDR_INT(opR);
2324
for ( k = 0; k < (Int)rs; ++k, ++l, ++r ) {
2325
c = (Int)*l + (Int)*r + (Int)d;
2326
*l = (TypDigit)c;
2327
d = (TypDigit)(c>>NR_DIGIT_BITS);
2328
}
2329
c = *l + d; d = (TypDigit)(c>>NR_DIGIT_BITS);
2330
qi--;
2331
}
2332
2333
/* store the digit in the quotient */
2334
ADDR_INT(quo)[i] = qi;
2335
2336
}
2337
2338
/* remove the leading zeroes, note that there can't be more than 7 */
2339
q = ADDR_INT(quo) + SIZE_INT(quo);
2340
if ( SIZE_INT(quo) > 4
2341
&& q[-4] == 0 && q[-3] == 0 && q[-2] == 0 && q[-1] == 0 ) {
2342
ResizeBag( quo, (SIZE_INT(quo)-4)*sizeof(TypDigit) );
2343
}
2344
2345
/* reduce to small integer if possible */
2346
q = ADDR_INT(quo) + SIZE_INT(quo);
2347
if ( SIZE_INT(quo) == 4 && q[-2] == 0 && q[-1] == 0 ) {
2348
if ( TNUM_OBJ(quo) == T_INTPOS
2349
&& (UInt)(INTBASE*q[-3]+q[-4]) < (1L<<NR_SMALL_INT_BITS) )
2350
quo = INTOBJ_INT( INTBASE*q[-3]+q[-4] );
2351
else if ( TNUM_OBJ(quo) == T_INTNEG
2352
&& (UInt)(INTBASE*q[-3]+q[-4]) <= (1L<<NR_SMALL_INT_BITS) )
2353
quo = INTOBJ_INT( -(Int)(INTBASE*q[-3]+q[-4]) );
2354
}
2355
2356
}
2357
2358
/* return the result */
2359
return quo;
2360
}
2361
2362
2363
/****************************************************************************
2364
**
2365
*F FuncQUO_INT(<self>,<opL>,<opR>) . . . . . . . internal function 'QuoInt'
2366
**
2367
** 'FuncQUO_INT' implements the internal function 'QuoInt'.
2368
**
2369
** 'QuoInt( <i>, <k> )'
2370
**
2371
** 'Quo' returns the integer part of the quotient of its integer operands.
2372
** If <i> and <k> are positive 'Quo( <i>, <k> )' is the largest positive
2373
** integer <q> such that '<q> * <k> \<= <i>'. If <i> or <k> or both are
2374
** negative we define 'Abs( Quo(<i>,<k>) ) = Quo( Abs(<i>), Abs(<k>) )' and
2375
** 'Sign( Quo(<i>,<k>) ) = Sign(<i>) * Sign(<k>)'. Dividing by 0 causes an
2376
** error. 'Rem' (see "Rem") can be used to compute the remainder.
2377
*/
2378
Obj FuncQUO_INT (
2379
Obj self,
2380
Obj opL,
2381
Obj opR )
2382
{
2383
/* check the arguments */
2384
while ( TNUM_OBJ(opL) != T_INT
2385
&& TNUM_OBJ(opL) != T_INTPOS
2386
&& TNUM_OBJ(opL) != T_INTNEG ) {
2387
opL = ErrorReturnObj(
2388
"QuoInt: <left> must be an integer (not a %s)",
2389
(Int)TNAM_OBJ(opL), 0L,
2390
"you can replace <left> via 'return <left>;'" );
2391
}
2392
while ( TNUM_OBJ(opR) != T_INT
2393
&& TNUM_OBJ(opR) != T_INTPOS
2394
&& TNUM_OBJ(opR) != T_INTNEG ) {
2395
opR = ErrorReturnObj(
2396
"QuoInt: <right> must be an integer (not a %s)",
2397
(Int)TNAM_OBJ(opR), 0L,
2398
"you can replace <right> via 'return <right>;'" );
2399
}
2400
2401
/* return the quotient */
2402
return QuoInt( opL, opR );
2403
}
2404
2405
2406
/****************************************************************************
2407
**
2408
*F RemInt( <intL>, <intR> ) . . . . . . . . . . . remainder of two integers
2409
**
2410
** 'RemInt' returns the remainder of the quotient of the integers <intL>
2411
** and <intR>. 'RemInt' handles operands of type 'T_INT', 'T_INTPOS' and
2412
** 'T_INTNEG'.
2413
**
2414
** Note that the remainder is different from the value returned by the 'mod'
2415
** operator which is always positive.
2416
**
2417
** 'RemInt' is called from 'FunRemInt'.
2418
*/
2419
Obj RemInt (
2420
Obj opL,
2421
Obj opR )
2422
{
2423
Int i; /* loop count, value for small int */
2424
Int k; /* loop count, value for small int */
2425
UInt c; /* product of two digits */
2426
TypDigit d; /* carry into the next digit */
2427
TypDigit * l; /* pointer into the left operand */
2428
TypDigit * r; /* pointer into the right operand */
2429
TypDigit r1; /* leading digit of the right oper */
2430
TypDigit r2; /* next digit of the right operand */
2431
UInt rs; /* size of the right operand */
2432
UInt e; /* we mult r by 2^e so r1 >= 32768 */
2433
Obj rem; /* handle of the remainder bag */
2434
TypDigit * m; /* pointer into the remainder */
2435
UInt m01; /* leading two digits of the rem. */
2436
TypDigit m2; /* next digit of the remainder */
2437
TypDigit qi; /* guessed digit of the quotient */
2438
2439
/* compute the remainder of two small integers */
2440
if ( ARE_INTOBJS( opL, opR ) ) {
2441
2442
/* pathological case first */
2443
if ( opR == INTOBJ_INT(0) ) {
2444
opR = ErrorReturnObj(
2445
"Integer operations: <divisor> must be nonzero",
2446
0L, 0L,
2447
"you can replace the integer <divisor> via 'return <divisor>;'" );
2448
return QUO( opL, opR );
2449
}
2450
2451
/* get the integer values */
2452
i = INT_INTOBJ(opL);
2453
k = INT_INTOBJ(opR);
2454
2455
/* compute the remainder, make sure we divide only positive numbers*/
2456
if ( 0 <= i && 0 <= k ) i = ( i % k );
2457
else if ( 0 <= i && k < 0 ) i = ( i % -k );
2458
else if ( i < 0 && 0 <= k ) i = - ( -i % k );
2459
else if ( i < 0 && k < 0 ) i = - ( -i % -k );
2460
rem = INTOBJ_INT( i );
2461
2462
}
2463
2464
/* compute the remainder of a small integer by a large integer */
2465
else if ( IS_INTOBJ(opL) ) {
2466
2467
/* the small int -(1<<28) rem the large int (1<<28) is 0 */
2468
if ( opL == INTOBJ_INT(-(Int)(1L<<NR_SMALL_INT_BITS))
2469
&& TNUM_OBJ(opR) == T_INTPOS && SIZE_INT(opR) == 4
2470
&& ADDR_INT(opR)[3] == 0
2471
&& ADDR_INT(opR)[2] == 0
2472
&& ADDR_INT(opR)[1] == 1L << (NR_SMALL_INT_BITS-NR_DIGIT_BITS)
2473
&& ADDR_INT(opR)[0] == 0 )
2474
rem = INTOBJ_INT(0);
2475
2476
/* in all other cases the remainder is equal the left operand */
2477
else
2478
rem = opL;
2479
2480
}
2481
2482
/* compute the remainder of a large integer by a small integer */
2483
else if ( IS_INTOBJ(opR)
2484
&& INT_INTOBJ(opR) < INTBASE
2485
&& -(Int)INTBASE <= INT_INTOBJ(opR) ) {
2486
2487
/* pathological case first */
2488
if ( opR == INTOBJ_INT(0) ) {
2489
opR = ErrorReturnObj(
2490
"Integer operations: <divisor> must be nonzero",
2491
0L, 0L,
2492
"you can replace the integer <divisor> via 'return <divisor>;'" );
2493
return QUO( opL, opR );
2494
}
2495
2496
/* get the integer value, make positive */
2497
i = INT_INTOBJ(opR); if ( i < 0 ) i = -i;
2498
2499
/* maybe its trivial */
2500
if ( INTBASE % i == 0 ) {
2501
c = ADDR_INT(opL)[0] % i;
2502
}
2503
2504
/* otherwise run through the left operand and divide digitwise */
2505
else {
2506
l = ADDR_INT(opL) + SIZE_INT(opL) - 1;
2507
c = 0;
2508
for ( ; l >= ADDR_INT(opL); l-- ) {
2509
c = (c<<NR_DIGIT_BITS) + (Int)*l;
2510
c = c % i;
2511
}
2512
}
2513
2514
/* now c is the result, it has the same sign as the left operand */
2515
if ( TNUM_OBJ(opL) == T_INTPOS )
2516
rem = INTOBJ_INT( c );
2517
else
2518
rem = INTOBJ_INT( -(Int)c );
2519
2520
}
2521
2522
/* compute the remainder of a large integer modulo a large integer */
2523
else {
2524
2525
/* a small divisor larger than one digit isn't handled above */
2526
if ( IS_INTOBJ(opR) ) {
2527
if ( 0 < INT_INTOBJ(opR) ) {
2528
rem = NewBag( T_INTPOS, 4*sizeof(TypDigit) );
2529
ADDR_INT(rem)[0] = (TypDigit)(INT_INTOBJ(opR));
2530
ADDR_INT(rem)[1] = (TypDigit)(INT_INTOBJ(opR)>>NR_DIGIT_BITS);
2531
opR = rem;
2532
}
2533
else {
2534
rem = NewBag( T_INTNEG, 4*sizeof(TypDigit) );
2535
ADDR_INT(rem)[0] = (TypDigit)(-INT_INTOBJ(opR));
2536
ADDR_INT(rem)[1] = (TypDigit)((-INT_INTOBJ(opR))>>NR_DIGIT_BITS);
2537
opR = rem;
2538
}
2539
}
2540
2541
/* trivial case first */
2542
if ( SIZE_INT(opL) < SIZE_INT(opR) )
2543
return opL;
2544
2545
/* copy the left operand into a new bag, this holds the remainder */
2546
rem = NewBag( TNUM_OBJ(opL), (SIZE_INT(opL)+4)*sizeof(TypDigit) );
2547
l = ADDR_INT(opL);
2548
m = ADDR_INT(rem);
2549
for ( k = SIZE_INT(opL)-1; k >= 0; k-- )
2550
*m++ = *l++;
2551
2552
/* get the size of the right operand, and get the leading 2 digits */
2553
rs = SIZE_INT(opR);
2554
r = ADDR_INT(opR);
2555
while ( r[rs-1] == 0 ) rs--;
2556
for ( e = 0;
2557
((Int)r[rs-1]<<e) + (e ? r[rs-2]>>(NR_DIGIT_BITS-e): 0)
2558
< INTBASE/2; e++ ) ;
2559
2560
r1 = ((Int)r[rs-1]<<e) + (e ? r[rs-2]>>(NR_DIGIT_BITS-e):0);
2561
r2 = ((Int)r[rs-2]<<e) + ((e && rs>=3) ? r[rs-3]>>(NR_DIGIT_BITS-e) : 0);
2562
2563
/* run through the digits in the quotient */
2564
for ( i = SIZE_INT(rem)-SIZE_INT(opR)-1; i >= 0; i-- ) {
2565
2566
/* guess the factor */
2567
m = ADDR_INT(rem) + rs + i;
2568
m01 = ((INTBASE*m[0]+m[-1])<<e) + (e ? m[-2]>>(NR_DIGIT_BITS-e):0);
2569
if ( m01 == 0 ) continue;
2570
m2 = ((Int)m[-2]<<e) + ((e && rs+i>=3) ? m[-3]>>(NR_DIGIT_BITS-e) : 0);
2571
if ( ((Int)m[0]<<e)+(e ? m[-1]>>(NR_DIGIT_BITS-e): 0) < r1 )
2572
qi = m01 / r1;
2573
else qi = INTBASE - 1;
2574
while ( m01-(Int)qi*r1 < (Int)INTBASE
2575
&& INTBASE*(m01-(Int)qi*r1)+m2 < (Int)qi*r2 )
2576
qi--;
2577
2578
/* m = m - qi * r; */
2579
d = 0;
2580
m = ADDR_INT(rem) + i;
2581
r = ADDR_INT(opR);
2582
for ( k = 0; k < (Int)rs; ++k, ++m, ++r ) {
2583
c = (Int)*m - (Int)qi * *r - (Int)d;
2584
*m = (TypDigit)c;
2585
d = -(TypDigit)(c>>NR_DIGIT_BITS);
2586
}
2587
c = (Int)*m - d; *m = (TypDigit)c; d = -(TypDigit)(c>>NR_DIGIT_BITS);
2588
2589
/* if we have a borrow then add back */
2590
if ( d != 0 ) {
2591
d = 0;
2592
m = ADDR_INT(rem) + i;
2593
r = ADDR_INT(opR);
2594
for ( k = 0; k < (Int)rs; ++k, ++m, ++r ) {
2595
c = (Int)*m + (Int)*r + (Int)d;
2596
*m = (TypDigit)c;
2597
d = (TypDigit)(c>>NR_DIGIT_BITS);
2598
}
2599
c = (Int)*m + d; *m = (TypDigit)c; d = (TypDigit)(c>>NR_DIGIT_BITS);
2600
qi--;
2601
}
2602
2603
}
2604
2605
/* remove the leading zeroes */
2606
m = ADDR_INT(rem) + SIZE_INT(rem);
2607
if ( (m[-4] == 0 && m[-3] == 0 && m[-2] == 0 && m[-1] == 0)
2608
|| (SIZE_INT(rem) == 4 && m[-2] == 0 && m[-1] == 0) ) {
2609
2610
/* find the number of significant digits */
2611
m = ADDR_INT(rem);
2612
for ( k = SIZE_INT(rem); k != 0; k-- ) {
2613
if ( m[k-1] != 0 )
2614
break;
2615
}
2616
2617
/* reduce to small integer if possible, otherwise shrink bag */
2618
if ( k <= 2 && TNUM_OBJ(rem) == T_INTPOS
2619
&& (UInt)(INTBASE*m[1]+m[0]) < (1L<<NR_SMALL_INT_BITS) )
2620
rem = INTOBJ_INT( INTBASE*m[1]+m[0] );
2621
else if ( k <= 2 && TNUM_OBJ(rem) == T_INTNEG
2622
&& (UInt)(INTBASE*m[1]+m[0]) <= (1L<<NR_SMALL_INT_BITS) )
2623
rem = INTOBJ_INT( -(Int)(INTBASE*m[1]+m[0]) );
2624
else
2625
ResizeBag( rem, (((k + 3) / 4) * 4) * sizeof(TypDigit) );
2626
}
2627
2628
}
2629
2630
/* return the result */
2631
return rem;
2632
}
2633
2634
2635
/****************************************************************************
2636
**
2637
*F FuncREM_INT(<self>,<opL>,<opR>) . . . . . . . internal function 'RemInt'
2638
**
2639
** 'FuncREM_INT' implements the internal function 'RemInt'.
2640
**
2641
** 'RemInt( <i>, <k> )'
2642
**
2643
** 'Rem' returns the remainder of its two integer operands, i.e., if <k> is
2644
** not equal to zero 'Rem( <i>, <k> ) = <i> - <k> * Quo( <i>, <k> )'. Note
2645
** that the rules given for 'Quo' (see "Quo") imply that 'Rem( <i>, <k> )'
2646
** has the same sign as <i> and its absolute value is strictly less than the
2647
** absolute value of <k>. Dividing by 0 causes an error.
2648
*/
2649
Obj FuncREM_INT (
2650
Obj self,
2651
Obj opL,
2652
Obj opR )
2653
{
2654
/* check the arguments */
2655
while ( TNUM_OBJ(opL) != T_INT
2656
&& TNUM_OBJ(opL) != T_INTPOS
2657
&& TNUM_OBJ(opL) != T_INTNEG ) {
2658
opL = ErrorReturnObj(
2659
"RemInt: <left> must be an integer (not a %s)",
2660
(Int)TNAM_OBJ(opL), 0L,
2661
"you can replace <left> via 'return <left>;'" );
2662
}
2663
while ( TNUM_OBJ(opR) != T_INT
2664
&& TNUM_OBJ(opR) != T_INTPOS
2665
&& TNUM_OBJ(opR) != T_INTNEG ) {
2666
opR = ErrorReturnObj(
2667
"RemInt: <right> must be an integer (not a %s)",
2668
(Int)TNAM_OBJ(opR), 0L,
2669
"you can replace <right> via 'return <right>;'" );
2670
}
2671
2672
/* return the remainder */
2673
return RemInt( opL, opR );
2674
}
2675
2676
2677
/****************************************************************************
2678
**
2679
*F GcdInt( <opL>, <opR> ) . . . . . . . . . . . . . . . gcd of two integers
2680
**
2681
** 'GcdInt' returns the gcd of the two integers <opL> and <opR>.
2682
**
2683
** It is called from 'FunGcdInt' and the rational package.
2684
*/
2685
Obj GcdInt (
2686
Obj opL,
2687
Obj opR )
2688
{
2689
Int i; /* loop count, value for small int */
2690
Int k; /* loop count, value for small int */
2691
UInt c; /* product of two digits */
2692
TypDigit d; /* carry into the next digit */
2693
TypDigit * r; /* pointer into the right operand */
2694
TypDigit r1; /* leading digit of the right oper */
2695
TypDigit r2; /* next digit of the right operand */
2696
UInt rs; /* size of the right operand */
2697
UInt e; /* we mult r by 2^e so r1 >= 32768 */
2698
TypDigit * l; /* pointer into the left operand */
2699
UInt l01; /* leading two digits of the rem. */
2700
TypDigit l2; /* next digit of the remainder */
2701
UInt ls; /* size of the left operand */
2702
TypDigit qi; /* guessed digit of the quotient */
2703
Obj gcd; /* handle of the result */
2704
2705
/* compute the gcd of two small integers */
2706
if ( ARE_INTOBJS( opL, opR ) ) {
2707
2708
/* get the integer values, make them positive */
2709
i = INT_INTOBJ(opL); if ( i < 0 ) i = -i;
2710
k = INT_INTOBJ(opR); if ( k < 0 ) k = -k;
2711
2712
/* compute the gcd using Euclids algorithm */
2713
while ( k != 0 ) {
2714
c = k;
2715
k = i % k;
2716
i = c;
2717
}
2718
2719
/* now i is the result */
2720
if ( i == (1L<<NR_SMALL_INT_BITS) ) {
2721
gcd = NewBag( T_INTPOS, 4*sizeof(TypDigit) );
2722
ADDR_INT(gcd)[0] = (TypDigit)i;
2723
ADDR_INT(gcd)[1] = (TypDigit)(i>>NR_DIGIT_BITS);
2724
}
2725
else {
2726
gcd = INTOBJ_INT( i );
2727
}
2728
2729
}
2730
2731
/* compute the gcd of a small and a large integer */
2732
else if ( (IS_INTOBJ(opL)
2733
&& INT_INTOBJ(opL) < INTBASE && -(Int)INTBASE <= INT_INTOBJ(opL))
2734
|| (IS_INTOBJ(opR)
2735
&& INT_INTOBJ(opR) < INTBASE && -(Int)INTBASE <= INT_INTOBJ(opR)) ) {
2736
2737
/* make the right operand the small one */
2738
if ( IS_INTOBJ(opL) ) {
2739
gcd = opL; opL = opR; opR = gcd;
2740
}
2741
2742
/* maybe it's trivial */
2743
if ( opR == INTOBJ_INT(0) ) {
2744
if( TNUM_OBJ( opL ) == T_INTNEG ) {
2745
/* If opL is negative, change the sign. We do this by */
2746
/* copying opL into a bag of type T_INTPOS. Note that */
2747
/* opL is a large negative number, so it cannot be the */
2748
/* the negative of 1 << NR_SMALL_INT_BITS. */
2749
gcd = NewBag( T_INTPOS, SIZE_OBJ(opL) );
2750
l = ADDR_INT(opL); r = ADDR_INT(gcd);
2751
for ( k = SIZE_INT(opL); k != 0; k-- ) *r++ = *l++;
2752
2753
return gcd;
2754
}
2755
else return opL;
2756
}
2757
2758
/* get the right operand value, make it positive */
2759
i = INT_INTOBJ(opR); if ( i < 0 ) i = -i;
2760
2761
/* do one remainder operation */
2762
l = ADDR_INT(opL) + SIZE_INT(opL) - 1;
2763
c = 0;
2764
for ( ; l >= ADDR_INT(opL); l-- ) {
2765
c = (c<<NR_DIGIT_BITS) + (Int)*l;
2766
c = c % i;
2767
}
2768
k = c;
2769
2770
/* compute the gcd using Euclids algorithm */
2771
while ( k != 0 ) {
2772
c = k;
2773
k = i % k;
2774
i = c;
2775
}
2776
2777
/* now i is the result */
2778
if ( i == (1L<<NR_SMALL_INT_BITS) ) {
2779
gcd = NewBag( T_INTPOS, 4*sizeof(TypDigit) );
2780
ADDR_INT(gcd)[0] = 0;
2781
ADDR_INT(gcd)[1] = 1L<<(NR_SMALL_INT_BITS-NR_DIGIT_BITS);
2782
}
2783
else {
2784
gcd = INTOBJ_INT( i );
2785
}
2786
2787
}
2788
2789
/* compute the gcd of two large integers */
2790
else {
2791
2792
/* a small divisor larger than one digit isn't handled above */
2793
if ( IS_INTOBJ(opL) ) {
2794
if ( 0 < INT_INTOBJ(opL) ) {
2795
gcd = NewBag( T_INTPOS, 4*sizeof(TypDigit) );
2796
ADDR_INT(gcd)[0] = (TypDigit)(INT_INTOBJ(opL));
2797
ADDR_INT(gcd)[1] = (TypDigit)(INT_INTOBJ(opL)>>NR_DIGIT_BITS);
2798
opL = gcd;
2799
}
2800
else {
2801
gcd = NewBag( T_INTNEG, 4*sizeof(TypDigit) );
2802
ADDR_INT(gcd)[0] = (TypDigit)(-INT_INTOBJ(opL));
2803
ADDR_INT(gcd)[1] = (TypDigit)((-INT_INTOBJ(opL))>>NR_DIGIT_BITS);
2804
opL = gcd;
2805
}
2806
}
2807
2808
/* a small dividend larger than one digit isn't handled above */
2809
if ( IS_INTOBJ(opR) ) {
2810
if ( 0 < INT_INTOBJ(opR) ) {
2811
gcd = NewBag( T_INTPOS, 4*sizeof(TypDigit) );
2812
ADDR_INT(gcd)[0] = (TypDigit)(INT_INTOBJ(opR));
2813
ADDR_INT(gcd)[1] = (TypDigit)(INT_INTOBJ(opR)>>NR_DIGIT_BITS);
2814
opR = gcd;
2815
}
2816
else {
2817
gcd = NewBag( T_INTNEG, 4*sizeof(TypDigit) );
2818
ADDR_INT(gcd)[0] = (TypDigit)(-INT_INTOBJ(opR));
2819
ADDR_INT(gcd)[1] = (TypDigit)((-INT_INTOBJ(opR))>>NR_DIGIT_BITS);
2820
opR = gcd;
2821
}
2822
}
2823
2824
/* copy the left operand into a new bag */
2825
gcd = NewBag( T_INTPOS, (SIZE_INT(opL)+4)*sizeof(TypDigit) );
2826
l = ADDR_INT(opL);
2827
r = ADDR_INT(gcd);
2828
for ( k = SIZE_INT(opL)-1; k >= 0; k-- )
2829
*r++ = *l++;
2830
opL = gcd;
2831
2832
/* get the size of the left operand */
2833
ls = SIZE_INT(opL);
2834
l = ADDR_INT(opL);
2835
while ( ls >= 1 && l[ls-1] == 0 ) ls--;
2836
2837
/* copy the right operand into a new bag */
2838
gcd = NewBag( T_INTPOS, (SIZE_INT(opR)+4)*sizeof(TypDigit) );
2839
r = ADDR_INT(opR);
2840
l = ADDR_INT(gcd);
2841
for ( k = SIZE_INT(opR)-1; k >= 0; k-- )
2842
*l++ = *r++;
2843
opR = gcd;
2844
2845
/* get the size of the right operand */
2846
rs = SIZE_INT(opR);
2847
r = ADDR_INT(opR);
2848
while ( rs >= 1 && r[rs-1] == 0 ) rs--;
2849
2850
/* repeat while the right operand is large */
2851
while ( rs >= 2 ) {
2852
2853
/* get the leading two digits */
2854
for ( e = 0;
2855
((Int)r[rs-1]<<e) + (e ? r[rs-2]>>(NR_DIGIT_BITS-e) : 0)
2856
< INTBASE/2; e++ ) ;
2857
r1 = ((Int)r[rs-1]<<e) + (e ? r[rs-2]>>(NR_DIGIT_BITS-e): 0);
2858
r2 = ((Int)r[rs-2]<<e) + ((e && rs>=3) ? r[rs-3]>>(NR_DIGIT_BITS-e) : 0);
2859
2860
/* run through the digits in the quotient */
2861
for ( i = ls - rs; i >= 0; i-- ) {
2862
2863
/* guess the factor */
2864
l = ADDR_INT(opL) + rs + i;
2865
l01 = ((INTBASE*l[0]+l[-1])<<e) + (e ? l[-2]>>(NR_DIGIT_BITS-e):0);
2866
if ( l01 == 0 ) continue;
2867
l2 = ((Int)l[-2]<<e) + ((e && rs+i>=3) ? l[-3]>>(NR_DIGIT_BITS-e):0);
2868
if ( ((Int)l[0]<<e)+(e ? l[-1]>>(NR_DIGIT_BITS-e):0) < r1 )
2869
qi = l01 / r1;
2870
else qi = INTBASE - 1;
2871
while ( l01-(Int)qi*r1 < (Int)INTBASE
2872
&& (INTBASE*(l01-(Int)qi*r1)+l2) < ((Int)qi*r2) )
2873
qi--;
2874
2875
/* l = l - qi * r; */
2876
d = 0;
2877
l = ADDR_INT(opL) + i;
2878
r = ADDR_INT(opR);
2879
for ( k = 0; k < (Int)rs; ++k, ++l, ++r ) {
2880
c = (Int)*l - (Int)qi * *r - (Int)d;
2881
*l = (TypDigit)c;
2882
d = -(TypDigit)(c>>NR_DIGIT_BITS);
2883
}
2884
c = (Int)*l - d; *l = (TypDigit)c; d = -(TypDigit)(c>>NR_DIGIT_BITS);
2885
2886
/* if we have a borrow then add back */
2887
if ( d != 0 ) {
2888
d = 0;
2889
l = ADDR_INT(opL) + i;
2890
r = ADDR_INT(opR);
2891
for ( k = 0; k < (Int)rs; ++k, ++l, ++r ) {
2892
c = (Int)*l + (Int)*r + (Int)d;
2893
*l = (TypDigit)c;
2894
d = (TypDigit)(c>>NR_DIGIT_BITS);
2895
}
2896
c = *l + d; *l = (TypDigit)c; d = (TypDigit)(c>>NR_DIGIT_BITS);
2897
qi--;
2898
}
2899
}
2900
2901
/* exchange the two operands */
2902
gcd = opL; opL = opR; opR = gcd;
2903
ls = rs;
2904
2905
/* get the size of the right operand */
2906
rs = SIZE_INT(opR);
2907
r = ADDR_INT(opR);
2908
while ( rs >= 1 && r[rs-1] == 0 ) rs--;
2909
2910
}
2911
2912
/* if the right operand is zero now, the left is the gcd */
2913
if ( rs == 0 ) {
2914
2915
/* remove the leading zeroes */
2916
l = ADDR_INT(opL) + SIZE_INT(opL);
2917
if ( (l[-4] == 0 && l[-3] == 0 && l[-2] == 0 && l[-1] == 0)
2918
|| (SIZE_INT(opL) == 4 && l[-2] == 0 && l[-1] == 0) ) {
2919
2920
/* find the number of significant digits */
2921
l = ADDR_INT(opL);
2922
for ( k = SIZE_INT(opL); k != 0; k-- ) {
2923
if ( l[k-1] != 0 )
2924
break;
2925
}
2926
2927
/* reduce to small integer if possible, otherwise shrink b */
2928
if ( k <= 2 && TNUM_OBJ(opL) == T_INTPOS
2929
&& (UInt)(INTBASE*l[1]+l[0]) < (1L<<NR_SMALL_INT_BITS) )
2930
opL = INTOBJ_INT( INTBASE*l[1]+l[0] );
2931
else if ( k <= 2 && TNUM_OBJ(opL) == T_INTNEG
2932
&& (UInt)(INTBASE*l[1]+l[0]) <= (1L<<NR_SMALL_INT_BITS) )
2933
opL = INTOBJ_INT( -(Int)(INTBASE*l[1]+l[0]) );
2934
else
2935
ResizeBag( opL, (((k + 3) / 4) * 4) * sizeof(TypDigit) );
2936
}
2937
2938
gcd = opL;
2939
2940
}
2941
2942
/* otherwise handle one large and one small integer as above */
2943
else {
2944
2945
/* get the right operand value, make it positive */
2946
i = r[0];
2947
2948
/* do one remainder operation */
2949
l = ADDR_INT(opL) + SIZE_INT(opL) - 1;
2950
c = 0;
2951
for ( ; l >= ADDR_INT(opL); l-- ) {
2952
c = (c<<NR_DIGIT_BITS) + (Int)*l;
2953
c = c % i;
2954
}
2955
k = c;
2956
2957
/* compute the gcd using Euclids algorithm */
2958
while ( k != 0 ) {
2959
c = k;
2960
k = i % k;
2961
i = c;
2962
}
2963
2964
/* now i is the result */
2965
if ( i == (1L<<NR_SMALL_INT_BITS) ) {
2966
gcd = NewBag( T_INTPOS, 4*sizeof(TypDigit) );
2967
ADDR_INT(gcd)[0] = 0;
2968
ADDR_INT(gcd)[1] = 1L<<(NR_SMALL_INT_BITS-NR_DIGIT_BITS);
2969
}
2970
else {
2971
gcd = INTOBJ_INT( i );
2972
}
2973
2974
}
2975
2976
}
2977
2978
/* return the result */
2979
return gcd;
2980
}
2981
2982
2983
/****************************************************************************
2984
**
2985
*F FuncGCD_INT(<self>,<opL>,<opR>) . . . . . . . internal function 'GcdInt'
2986
**
2987
** 'FuncGCD_INT' implements the internal function 'GcdInt'.
2988
**
2989
** 'GcdInt( <i>, <k> )'
2990
**
2991
** 'Gcd' returns the greatest common divisor of the two integers <m> and
2992
** <n>, i.e., the greatest integer that divides both <m> and <n>. The
2993
** greatest common divisor is never negative, even if the arguments are. We
2994
** define $gcd( m, 0 ) = gcd( 0, m ) = abs( m )$ and $gcd( 0, 0 ) = 0$.
2995
*/
2996
Obj FuncGCD_INT (
2997
Obj self,
2998
Obj opL,
2999
Obj opR )
3000
{
3001
/* check the arguments */
3002
while ( TNUM_OBJ(opL) != T_INT
3003
&& TNUM_OBJ(opL) != T_INTPOS
3004
&& TNUM_OBJ(opL) != T_INTNEG ) {
3005
opL = ErrorReturnObj(
3006
"GcdInt: <left> must be an integer (not a %s)",
3007
(Int)TNAM_OBJ(opL), 0L,
3008
"you can replace <left> via 'return <left>;'" );
3009
}
3010
while ( TNUM_OBJ(opR) != T_INT
3011
&& TNUM_OBJ(opR) != T_INTPOS
3012
&& TNUM_OBJ(opR) != T_INTNEG ) {
3013
opR = ErrorReturnObj(
3014
"GcdInt: <right> must be an integer (not a %s)",
3015
(Int)TNAM_OBJ(opR), 0L,
3016
"you can replace <right> via 'return <right>;'" );
3017
}
3018
3019
/* return the gcd */
3020
return GcdInt( opL, opR );
3021
}
3022
3023
3024
3025
3026
3027
3028
/****************************************************************************
3029
**
3030
*F SaveInt( <int> )
3031
**
3032
** Since the type is saved, we don't need to worry about sign
3033
*/
3034
3035
void SaveInt( Obj bigint)
3036
{
3037
TypDigit *ptr;
3038
UInt i;
3039
ptr = (TypDigit *)ADDR_OBJ(bigint);
3040
for (i = 0; i < SIZE_INT(bigint); i++)
3041
#ifdef SYS_IS_64_BIT
3042
SaveUInt4(*ptr++);
3043
#else
3044
SaveUInt2(*ptr++);
3045
#endif
3046
return;
3047
}
3048
3049
/****************************************************************************
3050
**
3051
*F LoadInt( <int> )
3052
**
3053
** Since the type is loaded, we don't need to worry about sign
3054
*/
3055
3056
void LoadInt( Obj bigint)
3057
{
3058
TypDigit *ptr;
3059
UInt i;
3060
ptr = (TypDigit *)ADDR_OBJ(bigint);
3061
for (i = 0; i < SIZE_INT(bigint); i++)
3062
#ifdef SYS_IS_64_BIT
3063
*ptr++ = LoadUInt4();
3064
#else
3065
*ptr++ = LoadUInt2();
3066
#endif
3067
return;
3068
}
3069
3070
3071
/****************************************************************************
3072
**
3073
** * * * * * * * "Mersenne twister" random numbers * * * * * * * * * * * * *
3074
**
3075
** Part of this code for fast generation of 32 bit pseudo random numbers with
3076
** a period of length 2^19937-1 and a 623-dimensional equidistribution is
3077
** taken from:
3078
** http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/emt.html
3079
** (Also look in Wikipedia for "Mersenne twister".)
3080
** We use the file mt19937ar.c, version 2002/1/26.
3081
*/
3082
3083
3084
3085
/****************************************************************************
3086
**
3087
*F RandomIntegerMT( <mtstr>, <nrbits> )
3088
**
3089
** Returns an integer with at most <nrbits> bits in uniform distribution.
3090
** <nrbits> must be a small integer. <mtstr> is a string as returned by
3091
** InitRandomMT.
3092
**
3093
** Implementation details are a bit tricky to obtain the same random
3094
** integers on 32 bit and 64 bit machines (which have different long
3095
** integer digit lengths and different ranges of small integers).
3096
**
3097
*/
3098
/* for comparison in case result is small int */
3099
Obj SMALLEST_INTPOS = NULL;
3100
3101
Obj FuncRandomIntegerMT(Obj self, Obj mtstr, Obj nrbits)
3102
{
3103
Obj res;
3104
Int i, n, q, r, qoff, len, nlen;
3105
UInt4 *mt, rand;
3106
TypDigit *pt;
3107
while (! IsStringConv(mtstr)) {
3108
mtstr = ErrorReturnObj(
3109
"<mtstr> must be a string, not a %s)",
3110
(Int)TNAM_OBJ(mtstr), 0L,
3111
"you can replace <mtstr> via 'return <mtstr>;'" );
3112
}
3113
while ((! IsStringConv(mtstr)) || GET_LEN_STRING(mtstr) < 2500) {
3114
mtstr = ErrorReturnObj(
3115
"<mtstr> must be a string with at least 2500 characters, ",
3116
0L, 0L,
3117
"you can replace <mtstr> via 'return <mtstr>;'" );
3118
}
3119
while ((! IS_INTOBJ(nrbits)) || INT_INTOBJ(nrbits) < 0) {
3120
nrbits = ErrorReturnObj(
3121
"<nrbits> must be a small non-negative integer, not a %s)",
3122
(Int)TNAM_OBJ(nrbits), 0L,
3123
"you can replace <mtstr> via 'return <mtstr>;'" );
3124
}
3125
n = INT_INTOBJ(nrbits);
3126
3127
/* small int case */
3128
if (n <= NR_SMALL_INT_BITS) {
3129
mt = (UInt4*) CHARS_STRING(mtstr);
3130
#ifdef SYS_IS_64_BIT
3131
if (n <= 32) {
3132
res = INTOBJ_INT((Int)(nextrandMT_int32(mt) & ((UInt4) -1L >> (32-n))));
3133
}
3134
else {
3135
unsigned long rd;
3136
rd = nextrandMT_int32(mt);
3137
rd += (unsigned long) ((UInt4) nextrandMT_int32(mt) &
3138
((UInt4) -1L >> (64-n))) << 32;
3139
res = INTOBJ_INT((Int)rd);
3140
}
3141
#else
3142
res = INTOBJ_INT((Int)(nextrandMT_int32(mt) & ((UInt4) -1L >> (32-n))));
3143
#endif
3144
}
3145
else {
3146
/* number of Digits */
3147
q = n / NR_DIGIT_BITS;
3148
r = n - q*NR_DIGIT_BITS;
3149
qoff = q + (r==0 ? 0:1);
3150
len = qoff;
3151
len = 4*((len+3) / 4);
3152
res = NewBag( T_INTPOS, len*sizeof(TypDigit) );
3153
pt = ADDR_INT(res);
3154
mt = (UInt4*) CHARS_STRING(mtstr);
3155
#ifdef SYS_IS_64_BIT
3156
for (i = 0; i < qoff; i++, pt++) {
3157
rand = (TypDigit) nextrandMT_int32(mt);
3158
*pt = rand;
3159
}
3160
#else
3161
for (i = 0; i < qoff/2; i++, pt++) {
3162
rand = nextrandMT_int32(mt);
3163
*pt = (TypDigit) (rand & 0xFFFFL);
3164
pt++;
3165
*pt = (TypDigit) ((rand & 0xFFFF0000L) >> 16);
3166
}
3167
if (2*i != qoff) {
3168
rand = nextrandMT_int32(mt);
3169
*pt = (TypDigit) (rand & 0xFFFFL);
3170
}
3171
#endif
3172
if (r != 0) {
3173
ADDR_INT(res)[qoff-1] = ADDR_INT(res)[qoff-1] & ((TypDigit)(-1)
3174
>> (NR_DIGIT_BITS-r));
3175
}
3176
/* shrink bag if necessary */
3177
for (nlen = len, i=0; nlen >0 && ADDR_INT(res)[nlen-1] == 0L;
3178
nlen--, i++);
3179
if (i/4 != 0) {
3180
ResizeBag(res, 4*((nlen+3) / 4)*sizeof(TypDigit));
3181
}
3182
/* convert result if small int */
3183
if (LtInt(res, SMALLEST_INTPOS)) {
3184
res = INTOBJ_INT((Int)(ADDR_INT(res)[0] +
3185
((UInt)ADDR_INT(res)[1] << NR_DIGIT_BITS)));
3186
}
3187
}
3188
3189
return res;
3190
}
3191
3192
3193
3194
/****************************************************************************
3195
**
3196
*F * * * * * * * * * * * * * initialize package * * * * * * * * * * * * * * *
3197
*/
3198
3199
/****************************************************************************
3200
**
3201
*V GVarFilts . . . . . . . . . . . . . . . . . . . list of filters to export
3202
*/
3203
static StructGVarFilt GVarFilts [] = {
3204
3205
{ "IS_INT", "obj", &IsIntFilt,
3206
FuncIS_INT, "src/integer.c:IS_INT" },
3207
3208
{ 0 }
3209
3210
};
3211
3212
3213
/****************************************************************************
3214
**
3215
*V GVarFuncs . . . . . . . . . . . . . . . . . . list of functions to export
3216
*/
3217
static StructGVarFunc GVarFuncs [] = {
3218
3219
{ "QUO_INT", 2, "int1, int2",
3220
FuncQUO_INT, "src/integer.c:QUO_INT" },
3221
3222
{ "REM_INT", 2, "int1, int2",
3223
FuncREM_INT, "src/integer.c:REM_INT" },
3224
3225
{ "GCD_INT", 2, "int1, int2",
3226
FuncGCD_INT, "src/integer.c:GCD_INT" },
3227
3228
{ "PROD_INT_OBJ", 2, "int, obj",
3229
FuncPROD_INT_OBJ, "src/integer.c:PROD_INT_OBJ" },
3230
3231
{ "POW_OBJ_INT", 2, "obj, int",
3232
FuncPOW_OBJ_INT, "src/integer.c:POW_OBJ_INT" },
3233
3234
{ "HexStringInt", 1, "integer",
3235
FuncHexStringInt, "src/integer.c:HexStringInt" },
3236
3237
{ "IntHexString", 1, "string",
3238
FuncIntHexString, "src/integer.c:IntHexString" },
3239
3240
{ "Log2Int", 1, "int",
3241
FuncLog2Int, "src/integer.c:Log2Int" },
3242
3243
{ "STRING_INT", 1, "int",
3244
FuncSTRING_INT, "src/integer.c:STRING_INT" },
3245
3246
3247
{ "RandomIntegerMT", 2, "mtstr, nrbits",
3248
FuncRandomIntegerMT, "src/integer.c:RandomIntegerMT" },
3249
3250
3251
{ 0 }
3252
3253
};
3254
3255
3256
/****************************************************************************
3257
**
3258
*F InitKernel( <module> ) . . . . . . . . initialise kernel data structures
3259
*/
3260
static Int InitKernel (
3261
StructInitInfo * module )
3262
{
3263
UInt t1, t2;
3264
3265
/* init filters and functions */
3266
InitHdlrFiltsFromTable( GVarFilts );
3267
InitHdlrFuncsFromTable( GVarFuncs );
3268
3269
/* install the marking functions */
3270
InfoBags[ T_INT ].name = "integer";
3271
#ifdef SYS_IS_64_BIT
3272
InfoBags[ T_INTPOS ].name = "integer (>= 2^60)";
3273
InfoBags[ T_INTNEG ].name = "integer (< -2^60)";
3274
#else
3275
InfoBags[ T_INTPOS ].name = "integer (>= 2^28)";
3276
InfoBags[ T_INTNEG ].name = "integer (< -2^28)";
3277
#endif
3278
InitMarkFuncBags( T_INTPOS, MarkNoSubBags );
3279
InitMarkFuncBags( T_INTNEG, MarkNoSubBags );
3280
3281
/* Install the saving methods */
3282
SaveObjFuncs [ T_INTPOS ] = SaveInt;
3283
SaveObjFuncs [ T_INTNEG ] = SaveInt;
3284
LoadObjFuncs [ T_INTPOS ] = LoadInt;
3285
LoadObjFuncs [ T_INTNEG ] = LoadInt;
3286
3287
/* install the printing function */
3288
PrintObjFuncs[ T_INT ] = PrintInt;
3289
PrintObjFuncs[ T_INTPOS ] = PrintInt;
3290
PrintObjFuncs[ T_INTNEG ] = PrintInt;
3291
3292
/* install the comparison methods */
3293
for ( t1 = T_INT; t1 <= T_INTNEG; t1++ ) {
3294
for ( t2 = T_INT; t2 <= T_INTNEG; t2++ ) {
3295
EqFuncs [ t1 ][ t2 ] = EqInt;
3296
LtFuncs [ t1 ][ t2 ] = LtInt;
3297
}
3298
}
3299
3300
/* install the unary arithmetic methods */
3301
for ( t1 = T_INT; t1 <= T_INTNEG; t1++ ) {
3302
ZeroFuncs[ t1 ] = ZeroInt;
3303
ZeroMutFuncs[ t1 ] = ZeroInt;
3304
AInvFuncs[ t1 ] = AInvInt;
3305
AInvMutFuncs[ t1 ] = AInvInt;
3306
OneFuncs [ t1 ] = OneInt;
3307
OneMutFuncs [ t1 ] = OneInt;
3308
}
3309
3310
/* install the default product and power methods */
3311
for ( t1 = T_INT; t1 <= T_INTNEG; t1++ ) {
3312
for ( t2 = FIRST_CONSTANT_TNUM; t2 <= LAST_CONSTANT_TNUM; t2++ ) {
3313
ProdFuncs[ t1 ][ t2 ] = ProdIntObj;
3314
PowFuncs [ t2 ][ t1 ] = PowObjInt;
3315
}
3316
for ( t2 = FIRST_RECORD_TNUM; t2 <= LAST_RECORD_TNUM; t2++ ) {
3317
ProdFuncs[ t1 ][ t2 ] = ProdIntObj;
3318
PowFuncs [ t2 ][ t1 ] = PowObjInt;
3319
}
3320
for ( t2 = FIRST_LIST_TNUM; t2 <= LAST_LIST_TNUM; t2++ ) {
3321
ProdFuncs[ t1 ][ t2 ] = ProdIntObj;
3322
PowFuncs [ t2 ][ t1 ] = PowObjInt;
3323
}
3324
}
3325
3326
/* install the binary arithmetic methods */
3327
for ( t1 = T_INT; t1 <= T_INTNEG; t1++ ) {
3328
for ( t2 = T_INT; t2 <= T_INTNEG; t2++ ) {
3329
EqFuncs [ t1 ][ t2 ] = EqInt;
3330
LtFuncs [ t1 ][ t2 ] = LtInt;
3331
SumFuncs [ t1 ][ t2 ] = SumInt;
3332
DiffFuncs[ t1 ][ t2 ] = DiffInt;
3333
ProdFuncs[ t1 ][ t2 ] = ProdInt;
3334
PowFuncs [ t1 ][ t2 ] = PowInt;
3335
ModFuncs [ t1 ][ t2 ] = ModInt;
3336
}
3337
}
3338
3339
/* gvars to import from the library */
3340
ImportGVarFromLibrary( "TYPE_INT_SMALL_ZERO", &TYPE_INT_SMALL_ZERO );
3341
ImportGVarFromLibrary( "TYPE_INT_SMALL_POS", &TYPE_INT_SMALL_POS );
3342
ImportGVarFromLibrary( "TYPE_INT_SMALL_NEG", &TYPE_INT_SMALL_NEG );
3343
ImportGVarFromLibrary( "TYPE_INT_LARGE_POS", &TYPE_INT_LARGE_POS );
3344
ImportGVarFromLibrary( "TYPE_INT_LARGE_NEG", &TYPE_INT_LARGE_NEG );
3345
ImportGVarFromLibrary( "SMALLEST_INTPOS", &SMALLEST_INTPOS );
3346
3347
ImportFuncFromLibrary( "String", &String );
3348
ImportFuncFromLibrary( "One", &OneAttr);
3349
3350
/* install the type functions */
3351
TypeObjFuncs[ T_INT ] = TypeIntSmall;
3352
TypeObjFuncs[ T_INTPOS ] = TypeIntLargePos;
3353
TypeObjFuncs[ T_INTNEG ] = TypeIntLargeNeg;
3354
3355
MakeBagTypePublic( T_INTPOS );
3356
MakeBagTypePublic( T_INTNEG );
3357
3358
/* return success */
3359
return 0;
3360
}
3361
3362
3363
/****************************************************************************
3364
**
3365
*F InitLibrary( <module> ) . . . . . . . initialise library data structures
3366
*/
3367
static Int InitLibrary (
3368
StructInitInfo * module )
3369
{
3370
UInt gvar;
3371
/* init filters and functions */
3372
InitGVarFiltsFromTable( GVarFilts );
3373
InitGVarFuncsFromTable( GVarFuncs );
3374
3375
/* hold smallest large integer */
3376
SMALLEST_INTPOS = NewBag( T_INTPOS, 4*sizeof(TypDigit) );
3377
ADDR_INT(SMALLEST_INTPOS)[1] =
3378
(TypDigit) (1L << (NR_SMALL_INT_BITS - NR_DIGIT_BITS));
3379
gvar = GVarName("SMALLEST_INTPOS");
3380
MakeReadWriteGVar( gvar);
3381
AssGVar( gvar, SMALLEST_INTPOS );
3382
MakeReadOnlyGVar(gvar);
3383
3384
/* return success */
3385
return 0;
3386
}
3387
3388
3389
/****************************************************************************
3390
**
3391
*F InitInfoInt() . . . . . . . . . . . . . . . . . . table of init functions
3392
*/
3393
static StructInitInfo module = {
3394
MODULE_BUILTIN, /* type */
3395
"integer", /* name */
3396
0, /* revision entry of c file */
3397
0, /* revision entry of h file */
3398
0, /* version */
3399
0, /* crc */
3400
InitKernel, /* initKernel */
3401
InitLibrary, /* initLibrary */
3402
0, /* checkInit */
3403
0, /* preSave */
3404
0, /* postSave */
3405
0 /* postRestore */
3406
};
3407
3408
StructInitInfo * InitInfoInt ( void )
3409
{
3410
return &module;
3411
}
3412
3413
/* matches the USE_GMP test at the top of the file */
3414
#endif
3415
3416
/****************************************************************************
3417
**
3418
*E integer.c . . . . . . . . . . . . . . . . . . . . . . . . . . . ends here
3419
*/
3420
3421