2 Copyright (C) 1989-94 Massachusetts Institute of Technology
3 Portions copyright (C) 2004-2008 Slava Pestov
5 This material was developed by the Scheme project at the Massachusetts
6 Institute of Technology, Department of Electrical Engineering and
7 Computer Science. Permission to copy and modify this software, to
8 redistribute either the original software or a modified version, and
9 to use this software for any purpose is granted, subject to the
10 following restrictions and understandings.
12 1. Any copy made of this software must include this copyright notice
15 2. Users of this software agree to make their best efforts (a) to
16 return to the MIT Scheme project any improvements or extensions that
17 they make, so that these may be included in future releases; and (b)
18 to inform MIT of noteworthy uses of this software.
20 3. All materials developed as a consequence of the use of this
21 software shall duly acknowledge such use, in accordance with the usual
22 standards of acknowledging credit in academic research.
24 4. MIT has made no warrantee or representation that the operation of
25 this software will be error-free, and MIT is under no obligation to
26 provide any services, by way of maintenance, update, or otherwise.
28 5. In conjunction with products arising from the use of this material,
29 there shall be no use of the name of the Massachusetts Institute of
30 Technology nor of any adaptation thereof in any advertising,
31 promotional, or sales literature without prior written consent from
34 /* Changes for Scheme 48:
35 * - Converted to ANSI.
36 * - Added bitwise operations.
37 * - Added s48 to the beginning of all externally visible names.
38 * - Cached the bignum representations of -1, 0, and 1.
41 /* Changes for Factor:
42 * - Adapt bignumint.h for Factor memory manager
43 * - Add more bignum <-> C type conversions
44 * - Remove unused functions
45 * - Add local variable GC root recording
46 * - Remove s48 prefix from function names
47 * - Various fixes for Win64
63 int factorvm::bignum_equal_p(bignum * x, bignum * y)
68 : ((! (BIGNUM_ZERO_P (y)))
69 && ((BIGNUM_NEGATIVE_P (x))
70 ? (BIGNUM_NEGATIVE_P (y))
71 : (! (BIGNUM_NEGATIVE_P (y))))
72 && (bignum_equal_p_unsigned (x, y))));
75 enum bignum_comparison factorvm::bignum_compare(bignum * x, bignum * y)
79 ? ((BIGNUM_ZERO_P (y))
80 ? bignum_comparison_equal
81 : (BIGNUM_NEGATIVE_P (y))
82 ? bignum_comparison_greater
83 : bignum_comparison_less)
85 ? ((BIGNUM_NEGATIVE_P (x))
86 ? bignum_comparison_less
87 : bignum_comparison_greater)
88 : (BIGNUM_NEGATIVE_P (x))
89 ? ((BIGNUM_NEGATIVE_P (y))
90 ? (bignum_compare_unsigned (y, x))
91 : (bignum_comparison_less))
92 : ((BIGNUM_NEGATIVE_P (y))
93 ? (bignum_comparison_greater)
94 : (bignum_compare_unsigned (x, y))));
97 /* allocates memory */
98 bignum *factorvm::bignum_add(bignum * x, bignum * y)
103 : (BIGNUM_ZERO_P (y))
105 : ((BIGNUM_NEGATIVE_P (x))
106 ? ((BIGNUM_NEGATIVE_P (y))
107 ? (bignum_add_unsigned (x, y, 1))
108 : (bignum_subtract_unsigned (y, x)))
109 : ((BIGNUM_NEGATIVE_P (y))
110 ? (bignum_subtract_unsigned (x, y))
111 : (bignum_add_unsigned (x, y, 0)))));
114 /* allocates memory */
115 bignum *factorvm::bignum_subtract(bignum * x, bignum * y)
119 ? ((BIGNUM_ZERO_P (y))
121 : (bignum_new_sign (y, (! (BIGNUM_NEGATIVE_P (y))))))
122 : ((BIGNUM_ZERO_P (y))
124 : ((BIGNUM_NEGATIVE_P (x))
125 ? ((BIGNUM_NEGATIVE_P (y))
126 ? (bignum_subtract_unsigned (y, x))
127 : (bignum_add_unsigned (x, y, 1)))
128 : ((BIGNUM_NEGATIVE_P (y))
129 ? (bignum_add_unsigned (x, y, 0))
130 : (bignum_subtract_unsigned (x, y))))));
133 /* allocates memory */
134 bignum *factorvm::bignum_multiply(bignum * x, bignum * y)
136 bignum_length_type x_length = (BIGNUM_LENGTH (x));
137 bignum_length_type y_length = (BIGNUM_LENGTH (y));
139 ((BIGNUM_NEGATIVE_P (x))
140 ? (! (BIGNUM_NEGATIVE_P (y)))
141 : (BIGNUM_NEGATIVE_P (y)));
142 if (BIGNUM_ZERO_P (x))
144 if (BIGNUM_ZERO_P (y))
148 bignum_digit_type digit = (BIGNUM_REF (x, 0));
150 return (bignum_maybe_new_sign (y, negative_p));
151 if (digit < BIGNUM_RADIX_ROOT)
152 return (bignum_multiply_unsigned_small_factor (y, digit, negative_p));
156 bignum_digit_type digit = (BIGNUM_REF (y, 0));
158 return (bignum_maybe_new_sign (x, negative_p));
159 if (digit < BIGNUM_RADIX_ROOT)
160 return (bignum_multiply_unsigned_small_factor (x, digit, negative_p));
162 return (bignum_multiply_unsigned (x, y, negative_p));
165 /* allocates memory */
166 void factorvm::bignum_divide(bignum * numerator, bignum * denominator, bignum * * quotient, bignum * * remainder)
168 if (BIGNUM_ZERO_P (denominator))
170 divide_by_zero_error();
173 if (BIGNUM_ZERO_P (numerator))
175 (*quotient) = numerator;
176 (*remainder) = numerator;
180 int r_negative_p = (BIGNUM_NEGATIVE_P (numerator));
182 ((BIGNUM_NEGATIVE_P (denominator)) ? (! r_negative_p) : r_negative_p);
183 switch (bignum_compare_unsigned (numerator, denominator))
185 case bignum_comparison_equal:
187 (*quotient) = (BIGNUM_ONE (q_negative_p));
188 (*remainder) = (BIGNUM_ZERO ());
191 case bignum_comparison_less:
193 (*quotient) = (BIGNUM_ZERO ());
194 (*remainder) = numerator;
197 case bignum_comparison_greater:
199 if ((BIGNUM_LENGTH (denominator)) == 1)
201 bignum_digit_type digit = (BIGNUM_REF (denominator, 0));
205 (bignum_maybe_new_sign (numerator, q_negative_p));
206 (*remainder) = (BIGNUM_ZERO ());
209 else if (digit < BIGNUM_RADIX_ROOT)
211 bignum_divide_unsigned_small_denominator
214 q_negative_p, r_negative_p);
219 bignum_divide_unsigned_medium_denominator
222 q_negative_p, r_negative_p);
226 bignum_divide_unsigned_large_denominator
227 (numerator, denominator,
229 q_negative_p, r_negative_p);
236 /* allocates memory */
237 bignum *factorvm::bignum_quotient(bignum * numerator, bignum * denominator)
239 if (BIGNUM_ZERO_P (denominator))
241 divide_by_zero_error();
242 return (BIGNUM_OUT_OF_BAND);
244 if (BIGNUM_ZERO_P (numerator))
248 ((BIGNUM_NEGATIVE_P (denominator))
249 ? (! (BIGNUM_NEGATIVE_P (numerator)))
250 : (BIGNUM_NEGATIVE_P (numerator)));
251 switch (bignum_compare_unsigned (numerator, denominator))
253 case bignum_comparison_equal:
254 return (BIGNUM_ONE (q_negative_p));
255 case bignum_comparison_less:
256 return (BIGNUM_ZERO ());
257 case bignum_comparison_greater:
258 default: /* to appease gcc -Wall */
261 if ((BIGNUM_LENGTH (denominator)) == 1)
263 bignum_digit_type digit = (BIGNUM_REF (denominator, 0));
265 return (bignum_maybe_new_sign (numerator, q_negative_p));
266 if (digit < BIGNUM_RADIX_ROOT)
267 bignum_divide_unsigned_small_denominator
269 ("ient), ((bignum * *) 0),
272 bignum_divide_unsigned_medium_denominator
274 ("ient), ((bignum * *) 0),
278 bignum_divide_unsigned_large_denominator
279 (numerator, denominator,
280 ("ient), ((bignum * *) 0),
288 /* allocates memory */
289 bignum *factorvm::bignum_remainder(bignum * numerator, bignum * denominator)
291 if (BIGNUM_ZERO_P (denominator))
293 divide_by_zero_error();
294 return (BIGNUM_OUT_OF_BAND);
296 if (BIGNUM_ZERO_P (numerator))
298 switch (bignum_compare_unsigned (numerator, denominator))
300 case bignum_comparison_equal:
301 return (BIGNUM_ZERO ());
302 case bignum_comparison_less:
304 case bignum_comparison_greater:
305 default: /* to appease gcc -Wall */
308 if ((BIGNUM_LENGTH (denominator)) == 1)
310 bignum_digit_type digit = (BIGNUM_REF (denominator, 0));
312 return (BIGNUM_ZERO ());
313 if (digit < BIGNUM_RADIX_ROOT)
315 (bignum_remainder_unsigned_small_denominator
316 (numerator, digit, (BIGNUM_NEGATIVE_P (numerator))));
317 bignum_divide_unsigned_medium_denominator
319 ((bignum * *) 0), (&remainder),
320 0, (BIGNUM_NEGATIVE_P (numerator)));
323 bignum_divide_unsigned_large_denominator
324 (numerator, denominator,
325 ((bignum * *) 0), (&remainder),
326 0, (BIGNUM_NEGATIVE_P (numerator)));
332 #define FOO_TO_BIGNUM(name,type,utype) \
333 bignum * factorvm::name##_to_bignum(type n) \
336 bignum_digit_type result_digits [BIGNUM_DIGITS_FOR(type)]; \
337 bignum_digit_type * end_digits = result_digits; \
338 /* Special cases win when these small constants are cached. */ \
339 if (n == 0) return (BIGNUM_ZERO ()); \
340 if (n == 1) return (BIGNUM_ONE (0)); \
341 if (n < (type)0 && n == (type)-1) return (BIGNUM_ONE (1)); \
343 utype accumulator = ((negative_p = (n < (type)0)) ? (-n) : n); \
346 (*end_digits++) = (accumulator & BIGNUM_DIGIT_MASK); \
347 accumulator >>= BIGNUM_DIGIT_LENGTH; \
349 while (accumulator != 0); \
353 (allot_bignum ((end_digits - result_digits), negative_p)); \
354 bignum_digit_type * scan_digits = result_digits; \
355 bignum_digit_type * scan_result = (BIGNUM_START_PTR (result)); \
356 while (scan_digits < end_digits) \
357 (*scan_result++) = (*scan_digits++); \
362 /* all below allocate memory */
363 FOO_TO_BIGNUM(cell,cell,cell)
364 FOO_TO_BIGNUM(fixnum,fixnum,cell)
365 FOO_TO_BIGNUM(long_long,s64,u64)
366 FOO_TO_BIGNUM(ulong_long,u64,u64)
368 #define BIGNUM_TO_FOO(name,type,utype) \
369 type factorvm::bignum_to_##name(bignum * bignum) \
371 if (BIGNUM_ZERO_P (bignum)) \
374 utype accumulator = 0; \
375 bignum_digit_type * start = (BIGNUM_START_PTR (bignum)); \
376 bignum_digit_type * scan = (start + (BIGNUM_LENGTH (bignum))); \
377 while (start < scan) \
378 accumulator = ((accumulator << BIGNUM_DIGIT_LENGTH) + (*--scan)); \
379 return ((BIGNUM_NEGATIVE_P (bignum)) ? (-((type)accumulator)) : accumulator); \
383 /* all of the below allocate memory */
384 BIGNUM_TO_FOO(cell,cell,cell);
385 BIGNUM_TO_FOO(fixnum,fixnum,cell);
386 BIGNUM_TO_FOO(long_long,s64,u64)
387 BIGNUM_TO_FOO(ulong_long,u64,u64)
389 double factorvm::bignum_to_double(bignum * bignum)
391 if (BIGNUM_ZERO_P (bignum))
394 double accumulator = 0;
395 bignum_digit_type * start = (BIGNUM_START_PTR (bignum));
396 bignum_digit_type * scan = (start + (BIGNUM_LENGTH (bignum)));
398 accumulator = ((accumulator * BIGNUM_RADIX) + (*--scan));
399 return ((BIGNUM_NEGATIVE_P (bignum)) ? (-accumulator) : accumulator);
403 #define DTB_WRITE_DIGIT(factor) \
405 significand *= (factor); \
406 digit = ((bignum_digit_type) significand); \
408 significand -= ((double) digit); \
411 /* allocates memory */
412 #define inf std::numeric_limits<double>::infinity()
414 bignum *factorvm::double_to_bignum(double x)
416 if (x == inf || x == -inf || x != x) return (BIGNUM_ZERO ());
418 double significand = (frexp (x, (&exponent)));
419 if (exponent <= 0) return (BIGNUM_ZERO ());
420 if (exponent == 1) return (BIGNUM_ONE (x < 0));
421 if (significand < 0) significand = (-significand);
423 bignum_length_type length = (BIGNUM_BITS_TO_DIGITS (exponent));
424 bignum * result = (allot_bignum (length, (x < 0)));
425 bignum_digit_type * start = (BIGNUM_START_PTR (result));
426 bignum_digit_type * scan = (start + length);
427 bignum_digit_type digit;
428 int odd_bits = (exponent % BIGNUM_DIGIT_LENGTH);
430 DTB_WRITE_DIGIT ((fixnum)1 << odd_bits);
433 if (significand == 0)
439 DTB_WRITE_DIGIT (BIGNUM_RADIX);
445 #undef DTB_WRITE_DIGIT
449 int factorvm::bignum_equal_p_unsigned(bignum * x, bignum * y)
451 bignum_length_type length = (BIGNUM_LENGTH (x));
452 if (length != (BIGNUM_LENGTH (y)))
456 bignum_digit_type * scan_x = (BIGNUM_START_PTR (x));
457 bignum_digit_type * scan_y = (BIGNUM_START_PTR (y));
458 bignum_digit_type * end_x = (scan_x + length);
459 while (scan_x < end_x)
460 if ((*scan_x++) != (*scan_y++))
466 enum bignum_comparison factorvm::bignum_compare_unsigned(bignum * x, bignum * y)
468 bignum_length_type x_length = (BIGNUM_LENGTH (x));
469 bignum_length_type y_length = (BIGNUM_LENGTH (y));
470 if (x_length < y_length)
471 return (bignum_comparison_less);
472 if (x_length > y_length)
473 return (bignum_comparison_greater);
475 bignum_digit_type * start_x = (BIGNUM_START_PTR (x));
476 bignum_digit_type * scan_x = (start_x + x_length);
477 bignum_digit_type * scan_y = ((BIGNUM_START_PTR (y)) + y_length);
478 while (start_x < scan_x)
480 bignum_digit_type digit_x = (*--scan_x);
481 bignum_digit_type digit_y = (*--scan_y);
482 if (digit_x < digit_y)
483 return (bignum_comparison_less);
484 if (digit_x > digit_y)
485 return (bignum_comparison_greater);
488 return (bignum_comparison_equal);
493 /* allocates memory */
494 bignum *factorvm::bignum_add_unsigned(bignum * x, bignum * y, int negative_p)
496 GC_BIGNUM(x); GC_BIGNUM(y);
498 if ((BIGNUM_LENGTH (y)) > (BIGNUM_LENGTH (x)))
505 bignum_length_type x_length = (BIGNUM_LENGTH (x));
507 bignum * r = (allot_bignum ((x_length + 1), negative_p));
509 bignum_digit_type sum;
510 bignum_digit_type carry = 0;
511 bignum_digit_type * scan_x = (BIGNUM_START_PTR (x));
512 bignum_digit_type * scan_r = (BIGNUM_START_PTR (r));
514 bignum_digit_type * scan_y = (BIGNUM_START_PTR (y));
515 bignum_digit_type * end_y = (scan_y + (BIGNUM_LENGTH (y)));
516 while (scan_y < end_y)
518 sum = ((*scan_x++) + (*scan_y++) + carry);
519 if (sum < BIGNUM_RADIX)
526 (*scan_r++) = (sum - BIGNUM_RADIX);
532 bignum_digit_type * end_x = ((BIGNUM_START_PTR (x)) + x_length);
534 while (scan_x < end_x)
536 sum = ((*scan_x++) + 1);
537 if (sum < BIGNUM_RADIX)
544 (*scan_r++) = (sum - BIGNUM_RADIX);
546 while (scan_x < end_x)
547 (*scan_r++) = (*scan_x++);
554 return (bignum_shorten_length (r, x_length));
560 /* allocates memory */
561 bignum *factorvm::bignum_subtract_unsigned(bignum * x, bignum * y)
563 GC_BIGNUM(x); GC_BIGNUM(y);
566 switch (bignum_compare_unsigned (x, y))
568 case bignum_comparison_equal:
569 return (BIGNUM_ZERO ());
570 case bignum_comparison_less:
578 case bignum_comparison_greater:
583 bignum_length_type x_length = (BIGNUM_LENGTH (x));
585 bignum * r = (allot_bignum (x_length, negative_p));
587 bignum_digit_type difference;
588 bignum_digit_type borrow = 0;
589 bignum_digit_type * scan_x = (BIGNUM_START_PTR (x));
590 bignum_digit_type * scan_r = (BIGNUM_START_PTR (r));
592 bignum_digit_type * scan_y = (BIGNUM_START_PTR (y));
593 bignum_digit_type * end_y = (scan_y + (BIGNUM_LENGTH (y)));
594 while (scan_y < end_y)
596 difference = (((*scan_x++) - (*scan_y++)) - borrow);
599 (*scan_r++) = (difference + BIGNUM_RADIX);
604 (*scan_r++) = difference;
610 bignum_digit_type * end_x = ((BIGNUM_START_PTR (x)) + x_length);
612 while (scan_x < end_x)
614 difference = ((*scan_x++) - borrow);
616 (*scan_r++) = (difference + BIGNUM_RADIX);
619 (*scan_r++) = difference;
624 BIGNUM_ASSERT (borrow == 0);
625 while (scan_x < end_x)
626 (*scan_r++) = (*scan_x++);
628 return (bignum_trim (r));
633 Maximum value for product_low or product_high:
634 ((R * R) + (R * (R - 2)) + (R - 1))
635 Maximum value for carry: ((R * (R - 1)) + (R - 1))
636 where R == BIGNUM_RADIX_ROOT */
638 /* allocates memory */
639 bignum *factorvm::bignum_multiply_unsigned(bignum * x, bignum * y, int negative_p)
641 GC_BIGNUM(x); GC_BIGNUM(y);
643 if ((BIGNUM_LENGTH (y)) > (BIGNUM_LENGTH (x)))
650 bignum_digit_type carry;
651 bignum_digit_type y_digit_low;
652 bignum_digit_type y_digit_high;
653 bignum_digit_type x_digit_low;
654 bignum_digit_type x_digit_high;
655 bignum_digit_type product_low;
656 bignum_digit_type * scan_r;
657 bignum_digit_type * scan_y;
658 bignum_length_type x_length = (BIGNUM_LENGTH (x));
659 bignum_length_type y_length = (BIGNUM_LENGTH (y));
662 (allot_bignum_zeroed ((x_length + y_length), negative_p));
664 bignum_digit_type * scan_x = (BIGNUM_START_PTR (x));
665 bignum_digit_type * end_x = (scan_x + x_length);
666 bignum_digit_type * start_y = (BIGNUM_START_PTR (y));
667 bignum_digit_type * end_y = (start_y + y_length);
668 bignum_digit_type * start_r = (BIGNUM_START_PTR (r));
669 #define x_digit x_digit_high
670 #define y_digit y_digit_high
671 #define product_high carry
672 while (scan_x < end_x)
674 x_digit = (*scan_x++);
675 x_digit_low = (HD_LOW (x_digit));
676 x_digit_high = (HD_HIGH (x_digit));
679 scan_r = (start_r++);
680 while (scan_y < end_y)
682 y_digit = (*scan_y++);
683 y_digit_low = (HD_LOW (y_digit));
684 y_digit_high = (HD_HIGH (y_digit));
687 (x_digit_low * y_digit_low) +
690 ((x_digit_high * y_digit_low) +
691 (x_digit_low * y_digit_high) +
692 (HD_HIGH (product_low)) +
695 (HD_CONS ((HD_LOW (product_high)), (HD_LOW (product_low))));
697 ((x_digit_high * y_digit_high) +
698 (HD_HIGH (product_high)));
702 return (bignum_trim (r));
709 /* allocates memory */
710 bignum *factorvm::bignum_multiply_unsigned_small_factor(bignum * x, bignum_digit_type y, int negative_p)
714 bignum_length_type length_x = (BIGNUM_LENGTH (x));
716 bignum * p = (allot_bignum ((length_x + 1), negative_p));
718 bignum_destructive_copy (x, p);
719 (BIGNUM_REF (p, length_x)) = 0;
720 bignum_destructive_scale_up (p, y);
721 return (bignum_trim (p));
724 void factorvm::bignum_destructive_add(bignum * bignum, bignum_digit_type n)
726 bignum_digit_type * scan = (BIGNUM_START_PTR (bignum));
727 bignum_digit_type digit;
728 digit = ((*scan) + n);
729 if (digit < BIGNUM_RADIX)
734 (*scan++) = (digit - BIGNUM_RADIX);
737 digit = ((*scan) + 1);
738 if (digit < BIGNUM_RADIX)
743 (*scan++) = (digit - BIGNUM_RADIX);
747 void factorvm::bignum_destructive_scale_up(bignum * bignum, bignum_digit_type factor)
749 bignum_digit_type carry = 0;
750 bignum_digit_type * scan = (BIGNUM_START_PTR (bignum));
751 bignum_digit_type two_digits;
752 bignum_digit_type product_low;
753 #define product_high carry
754 bignum_digit_type * end = (scan + (BIGNUM_LENGTH (bignum)));
755 BIGNUM_ASSERT ((factor > 1) && (factor < BIGNUM_RADIX_ROOT));
758 two_digits = (*scan);
759 product_low = ((factor * (HD_LOW (two_digits))) + (HD_LOW (carry)));
761 ((factor * (HD_HIGH (two_digits))) +
762 (HD_HIGH (product_low)) +
764 (*scan++) = (HD_CONS ((HD_LOW (product_high)), (HD_LOW (product_low))));
765 carry = (HD_HIGH (product_high));
767 /* A carry here would be an overflow, i.e. it would not fit.
768 Hopefully the callers allocate enough space that this will
771 BIGNUM_ASSERT (carry == 0);
778 /* For help understanding this algorithm, see:
779 Knuth, Donald E., "The Art of Computer Programming",
780 volume 2, "Seminumerical Algorithms"
781 section 4.3.1, "Multiple-Precision Arithmetic". */
783 /* allocates memory */
784 void factorvm::bignum_divide_unsigned_large_denominator(bignum * numerator, bignum * denominator, bignum * * quotient, bignum * * remainder, int q_negative_p, int r_negative_p)
786 GC_BIGNUM(numerator); GC_BIGNUM(denominator);
788 bignum_length_type length_n = ((BIGNUM_LENGTH (numerator)) + 1);
789 bignum_length_type length_d = (BIGNUM_LENGTH (denominator));
792 ((quotient != ((bignum * *) 0))
793 ? (allot_bignum ((length_n - length_d), q_negative_p))
794 : BIGNUM_OUT_OF_BAND);
797 bignum * u = (allot_bignum (length_n, r_negative_p));
801 BIGNUM_ASSERT (length_d > 1);
803 bignum_digit_type v1 = (BIGNUM_REF ((denominator), (length_d - 1)));
804 while (v1 < (BIGNUM_RADIX / 2))
812 bignum_destructive_copy (numerator, u);
813 (BIGNUM_REF (u, (length_n - 1))) = 0;
814 bignum_divide_unsigned_normalized (u, denominator, q);
818 bignum * v = (allot_bignum (length_d, 0));
820 bignum_destructive_normalization (numerator, u, shift);
821 bignum_destructive_normalization (denominator, v, shift);
822 bignum_divide_unsigned_normalized (u, v, q);
823 if (remainder != ((bignum * *) 0))
824 bignum_destructive_unnormalization (u, shift);
832 if (quotient != ((bignum * *) 0))
835 if (remainder != ((bignum * *) 0))
841 void factorvm::bignum_divide_unsigned_normalized(bignum * u, bignum * v, bignum * q)
843 bignum_length_type u_length = (BIGNUM_LENGTH (u));
844 bignum_length_type v_length = (BIGNUM_LENGTH (v));
845 bignum_digit_type * u_start = (BIGNUM_START_PTR (u));
846 bignum_digit_type * u_scan = (u_start + u_length);
847 bignum_digit_type * u_scan_limit = (u_start + v_length);
848 bignum_digit_type * u_scan_start = (u_scan - v_length);
849 bignum_digit_type * v_start = (BIGNUM_START_PTR (v));
850 bignum_digit_type * v_end = (v_start + v_length);
851 bignum_digit_type * q_scan = NULL;
852 bignum_digit_type v1 = (v_end[-1]);
853 bignum_digit_type v2 = (v_end[-2]);
854 bignum_digit_type ph; /* high half of double-digit product */
855 bignum_digit_type pl; /* low half of double-digit product */
856 bignum_digit_type guess;
857 bignum_digit_type gh; /* high half-digit of guess */
858 bignum_digit_type ch; /* high half of double-digit comparand */
859 bignum_digit_type v2l = (HD_LOW (v2));
860 bignum_digit_type v2h = (HD_HIGH (v2));
861 bignum_digit_type cl; /* low half of double-digit comparand */
862 #define gl ph /* low half-digit of guess */
865 bignum_digit_type gm; /* memory loc for reference parameter */
866 if (q != BIGNUM_OUT_OF_BAND)
867 q_scan = ((BIGNUM_START_PTR (q)) + (BIGNUM_LENGTH (q)));
868 while (u_scan_limit < u_scan)
874 (((((uj * BIGNUM_RADIX) + uj1) % v1) * BIGNUM_RADIX) + uj2);
875 guess = (((uj * BIGNUM_RADIX) + uj1) / v1); */
877 ch = (bignum_digit_divide (uj, (u_scan[-1]), v1, (&gm)));
883 ch = ((u_scan[-1]) + v1);
884 guess = (BIGNUM_RADIX - 1);
888 /* product = (guess * v2); */
889 gl = (HD_LOW (guess));
890 gh = (HD_HIGH (guess));
892 ph = ((v2l * gh) + (v2h * gl) + (HD_HIGH (pl)));
893 pl = (HD_CONS ((HD_LOW (ph)), (HD_LOW (pl))));
894 ph = ((v2h * gh) + (HD_HIGH (ph)));
895 /* if (comparand >= product) */
896 if ((ch > ph) || ((ch == ph) && (cl >= pl)))
899 /* comparand += (v1 << BIGNUM_DIGIT_LENGTH) */
901 /* if (comparand >= (BIGNUM_RADIX * BIGNUM_RADIX)) */
902 if (ch >= BIGNUM_RADIX)
905 qj = (bignum_divide_subtract (v_start, v_end, guess, (--u_scan_start)));
906 if (q != BIGNUM_OUT_OF_BAND)
915 bignum_digit_type factorvm::bignum_divide_subtract(bignum_digit_type * v_start, bignum_digit_type * v_end, bignum_digit_type guess, bignum_digit_type * u_start)
917 bignum_digit_type * v_scan = v_start;
918 bignum_digit_type * u_scan = u_start;
919 bignum_digit_type carry = 0;
920 if (guess == 0) return (0);
922 bignum_digit_type gl = (HD_LOW (guess));
923 bignum_digit_type gh = (HD_HIGH (guess));
925 bignum_digit_type pl;
926 bignum_digit_type vl;
930 while (v_scan < v_end)
935 pl = ((vl * gl) + (HD_LOW (carry)));
936 ph = ((vl * gh) + (vh * gl) + (HD_HIGH (pl)) + (HD_HIGH (carry)));
937 diff = ((*u_scan) - (HD_CONS ((HD_LOW (ph)), (HD_LOW (pl)))));
940 (*u_scan++) = (diff + BIGNUM_RADIX);
941 carry = ((vh * gh) + (HD_HIGH (ph)) + 1);
946 carry = ((vh * gh) + (HD_HIGH (ph)));
951 diff = ((*u_scan) - carry);
953 (*u_scan) = (diff + BIGNUM_RADIX);
963 /* Subtraction generated carry, implying guess is one too large.
964 Add v back in to bring it back down. */
968 while (v_scan < v_end)
970 bignum_digit_type sum = ((*v_scan++) + (*u_scan) + carry);
971 if (sum < BIGNUM_RADIX)
978 (*u_scan++) = (sum - BIGNUM_RADIX);
984 bignum_digit_type sum = ((*u_scan) + carry);
985 (*u_scan) = ((sum < BIGNUM_RADIX) ? sum : (sum - BIGNUM_RADIX));
990 /* allocates memory */
991 void factorvm::bignum_divide_unsigned_medium_denominator(bignum * numerator,bignum_digit_type denominator, bignum * * quotient, bignum * * remainder,int q_negative_p, int r_negative_p)
993 GC_BIGNUM(numerator);
995 bignum_length_type length_n = (BIGNUM_LENGTH (numerator));
996 bignum_length_type length_q;
1001 /* Because `bignum_digit_divide' requires a normalized denominator. */
1002 while (denominator < (BIGNUM_RADIX / 2))
1009 length_q = length_n;
1011 q = (allot_bignum (length_q, q_negative_p));
1012 bignum_destructive_copy (numerator, q);
1016 length_q = (length_n + 1);
1018 q = (allot_bignum (length_q, q_negative_p));
1019 bignum_destructive_normalization (numerator, q, shift);
1022 bignum_digit_type r = 0;
1023 bignum_digit_type * start = (BIGNUM_START_PTR (q));
1024 bignum_digit_type * scan = (start + length_q);
1025 bignum_digit_type qj;
1027 while (start < scan)
1029 r = (bignum_digit_divide (r, (*--scan), denominator, (&qj)));
1033 q = bignum_trim (q);
1035 if (remainder != ((bignum * *) 0))
1040 (*remainder) = (bignum_digit_to_bignum (r, r_negative_p));
1043 if (quotient != ((bignum * *) 0))
1049 void factorvm::bignum_destructive_normalization(bignum * source, bignum * target, int shift_left)
1051 bignum_digit_type digit;
1052 bignum_digit_type * scan_source = (BIGNUM_START_PTR (source));
1053 bignum_digit_type carry = 0;
1054 bignum_digit_type * scan_target = (BIGNUM_START_PTR (target));
1055 bignum_digit_type * end_source = (scan_source + (BIGNUM_LENGTH (source)));
1056 bignum_digit_type * end_target = (scan_target + (BIGNUM_LENGTH (target)));
1057 int shift_right = (BIGNUM_DIGIT_LENGTH - shift_left);
1058 bignum_digit_type mask = (((cell)1 << shift_right) - 1);
1059 while (scan_source < end_source)
1061 digit = (*scan_source++);
1062 (*scan_target++) = (((digit & mask) << shift_left) | carry);
1063 carry = (digit >> shift_right);
1065 if (scan_target < end_target)
1066 (*scan_target) = carry;
1068 BIGNUM_ASSERT (carry == 0);
1072 void factorvm::bignum_destructive_unnormalization(bignum * bignum, int shift_right)
1074 bignum_digit_type * start = (BIGNUM_START_PTR (bignum));
1075 bignum_digit_type * scan = (start + (BIGNUM_LENGTH (bignum)));
1076 bignum_digit_type digit;
1077 bignum_digit_type carry = 0;
1078 int shift_left = (BIGNUM_DIGIT_LENGTH - shift_right);
1079 bignum_digit_type mask = (((fixnum)1 << shift_right) - 1);
1080 while (start < scan)
1083 (*scan) = ((digit >> shift_right) | carry);
1084 carry = ((digit & mask) << shift_left);
1086 BIGNUM_ASSERT (carry == 0);
1090 /* This is a reduced version of the division algorithm, applied to the
1091 case of dividing two bignum digits by one bignum digit. It is
1092 assumed that the numerator, denominator are normalized. */
1094 #define BDD_STEP(qn, j) \
1099 uj_uj1 = (HD_CONS (uj, (u[j + 1]))); \
1100 guess = (uj_uj1 / v1); \
1101 comparand = (HD_CONS ((uj_uj1 % v1), (u[j + 2]))); \
1105 guess = (BIGNUM_RADIX_ROOT - 1); \
1106 comparand = (HD_CONS (((u[j + 1]) + v1), (u[j + 2]))); \
1108 while ((guess * v2) > comparand) \
1111 comparand += (v1 << BIGNUM_HALF_DIGIT_LENGTH); \
1112 if (comparand >= BIGNUM_RADIX) \
1115 qn = (bignum_digit_divide_subtract (v1, v2, guess, (&u[j]))); \
1118 bignum_digit_type factorvm::bignum_digit_divide(bignum_digit_type uh, bignum_digit_type ul, bignum_digit_type v, bignum_digit_type * q) /* return value */
1120 bignum_digit_type guess;
1121 bignum_digit_type comparand;
1122 bignum_digit_type v1 = (HD_HIGH (v));
1123 bignum_digit_type v2 = (HD_LOW (v));
1124 bignum_digit_type uj;
1125 bignum_digit_type uj_uj1;
1126 bignum_digit_type q1;
1127 bignum_digit_type q2;
1128 bignum_digit_type u [4];
1142 (u[0]) = (HD_HIGH (uh));
1143 (u[1]) = (HD_LOW (uh));
1144 (u[2]) = (HD_HIGH (ul));
1145 (u[3]) = (HD_LOW (ul));
1150 (*q) = (HD_CONS (q1, q2));
1151 return (HD_CONS ((u[2]), (u[3])));
1156 #define BDDS_MULSUB(vn, un, carry_in) \
1158 product = ((vn * guess) + carry_in); \
1159 diff = (un - (HD_LOW (product))); \
1162 un = (diff + BIGNUM_RADIX_ROOT); \
1163 carry = ((HD_HIGH (product)) + 1); \
1168 carry = (HD_HIGH (product)); \
1172 #define BDDS_ADD(vn, un, carry_in) \
1174 sum = (vn + un + carry_in); \
1175 if (sum < BIGNUM_RADIX_ROOT) \
1182 un = (sum - BIGNUM_RADIX_ROOT); \
1187 bignum_digit_type factorvm::bignum_digit_divide_subtract(bignum_digit_type v1, bignum_digit_type v2, bignum_digit_type guess, bignum_digit_type * u)
1190 bignum_digit_type product;
1191 bignum_digit_type diff;
1192 bignum_digit_type carry;
1193 BDDS_MULSUB (v2, (u[2]), 0);
1194 BDDS_MULSUB (v1, (u[1]), carry);
1197 diff = ((u[0]) - carry);
1199 (u[0]) = (diff + BIGNUM_RADIX);
1207 bignum_digit_type sum;
1208 bignum_digit_type carry;
1209 BDDS_ADD(v2, (u[2]), 0);
1210 BDDS_ADD(v1, (u[1]), carry);
1220 /* allocates memory */
1221 void factorvm::bignum_divide_unsigned_small_denominator(bignum * numerator, bignum_digit_type denominator, bignum * * quotient, bignum * * remainder,int q_negative_p, int r_negative_p)
1223 GC_BIGNUM(numerator);
1225 bignum * q = (bignum_new_sign (numerator, q_negative_p));
1228 bignum_digit_type r = (bignum_destructive_scale_down (q, denominator));
1230 q = (bignum_trim (q));
1232 if (remainder != ((bignum * *) 0))
1233 (*remainder) = (bignum_digit_to_bignum (r, r_negative_p));
1240 /* Given (denominator > 1), it is fairly easy to show that
1241 (quotient_high < BIGNUM_RADIX_ROOT), after which it is easy to see
1242 that all digits are < BIGNUM_RADIX. */
1244 bignum_digit_type factorvm::bignum_destructive_scale_down(bignum * bignum, bignum_digit_type denominator)
1246 bignum_digit_type numerator;
1247 bignum_digit_type remainder = 0;
1248 bignum_digit_type two_digits;
1249 #define quotient_high remainder
1250 bignum_digit_type * start = (BIGNUM_START_PTR (bignum));
1251 bignum_digit_type * scan = (start + (BIGNUM_LENGTH (bignum)));
1252 BIGNUM_ASSERT ((denominator > 1) && (denominator < BIGNUM_RADIX_ROOT));
1253 while (start < scan)
1255 two_digits = (*--scan);
1256 numerator = (HD_CONS (remainder, (HD_HIGH (two_digits))));
1257 quotient_high = (numerator / denominator);
1258 numerator = (HD_CONS ((numerator % denominator), (HD_LOW (two_digits))));
1259 (*scan) = (HD_CONS (quotient_high, (numerator / denominator)));
1260 remainder = (numerator % denominator);
1263 #undef quotient_high
1266 /* allocates memory */
1267 bignum * factorvm::bignum_remainder_unsigned_small_denominator(bignum * n, bignum_digit_type d, int negative_p)
1269 bignum_digit_type two_digits;
1270 bignum_digit_type * start = (BIGNUM_START_PTR (n));
1271 bignum_digit_type * scan = (start + (BIGNUM_LENGTH (n)));
1272 bignum_digit_type r = 0;
1273 BIGNUM_ASSERT ((d > 1) && (d < BIGNUM_RADIX_ROOT));
1274 while (start < scan)
1276 two_digits = (*--scan);
1278 ((HD_CONS (((HD_CONS (r, (HD_HIGH (two_digits)))) % d),
1279 (HD_LOW (two_digits))))
1282 return (bignum_digit_to_bignum (r, negative_p));
1285 /* allocates memory */
1286 bignum *factorvm::bignum_digit_to_bignum(bignum_digit_type digit, int negative_p)
1289 return (BIGNUM_ZERO ());
1292 bignum * result = (allot_bignum (1, negative_p));
1293 (BIGNUM_REF (result, 0)) = digit;
1298 /* allocates memory */
1299 bignum *factorvm::allot_bignum(bignum_length_type length, int negative_p)
1301 BIGNUM_ASSERT ((length >= 0) || (length < BIGNUM_RADIX));
1302 bignum * result = allot_array_internal<bignum>(length + 1);
1303 BIGNUM_SET_NEGATIVE_P (result, negative_p);
1307 /* allocates memory */
1308 bignum * factorvm::allot_bignum_zeroed(bignum_length_type length, int negative_p)
1310 bignum * result = allot_bignum(length,negative_p);
1311 bignum_digit_type * scan = (BIGNUM_START_PTR (result));
1312 bignum_digit_type * end = (scan + length);
1318 #define BIGNUM_REDUCE_LENGTH(source, length) \
1319 source = reallot_array(source,length + 1)
1321 /* allocates memory */
1322 bignum *factorvm::bignum_shorten_length(bignum * bignum, bignum_length_type length)
1324 bignum_length_type current_length = (BIGNUM_LENGTH (bignum));
1325 BIGNUM_ASSERT ((length >= 0) || (length <= current_length));
1326 if (length < current_length)
1328 BIGNUM_REDUCE_LENGTH (bignum, length);
1329 BIGNUM_SET_NEGATIVE_P (bignum, (length != 0) && (BIGNUM_NEGATIVE_P (bignum)));
1334 /* allocates memory */
1335 bignum *factorvm::bignum_trim(bignum * bignum)
1337 bignum_digit_type * start = (BIGNUM_START_PTR (bignum));
1338 bignum_digit_type * end = (start + (BIGNUM_LENGTH (bignum)));
1339 bignum_digit_type * scan = end;
1340 while ((start <= scan) && ((*--scan) == 0))
1345 bignum_length_type length = (scan - start);
1346 BIGNUM_REDUCE_LENGTH (bignum, length);
1347 BIGNUM_SET_NEGATIVE_P (bignum, (length != 0) && (BIGNUM_NEGATIVE_P (bignum)));
1354 /* allocates memory */
1355 bignum *factorvm::bignum_new_sign(bignum * x, int negative_p)
1358 bignum * result = (allot_bignum ((BIGNUM_LENGTH (x)), negative_p));
1360 bignum_destructive_copy (x, result);
1364 /* allocates memory */
1365 bignum *factorvm::bignum_maybe_new_sign(bignum * x, int negative_p)
1367 if ((BIGNUM_NEGATIVE_P (x)) ? negative_p : (! negative_p))
1372 (allot_bignum ((BIGNUM_LENGTH (x)), negative_p));
1373 bignum_destructive_copy (x, result);
1378 void factorvm::bignum_destructive_copy(bignum * source, bignum * target)
1380 bignum_digit_type * scan_source = (BIGNUM_START_PTR (source));
1381 bignum_digit_type * end_source =
1382 (scan_source + (BIGNUM_LENGTH (source)));
1383 bignum_digit_type * scan_target = (BIGNUM_START_PTR (target));
1384 while (scan_source < end_source)
1385 (*scan_target++) = (*scan_source++);
1390 * Added bitwise operations (and oddp).
1393 /* allocates memory */
1394 bignum *factorvm::bignum_bitwise_not(bignum * x)
1396 return bignum_subtract(BIGNUM_ONE(1), x);
1399 /* allocates memory */
1400 bignum *factorvm::bignum_arithmetic_shift(bignum * arg1, fixnum n)
1402 if (BIGNUM_NEGATIVE_P(arg1) && n < 0)
1403 return bignum_bitwise_not(bignum_magnitude_ash(bignum_bitwise_not(arg1), n));
1405 return bignum_magnitude_ash(arg1, n);
1412 /* allocates memory */
1413 bignum *factorvm::bignum_bitwise_and(bignum * arg1, bignum * arg2)
1416 (BIGNUM_NEGATIVE_P (arg1))
1417 ? (BIGNUM_NEGATIVE_P (arg2))
1418 ? bignum_negneg_bitwise_op(AND_OP, arg1, arg2)
1419 : bignum_posneg_bitwise_op(AND_OP, arg2, arg1)
1420 : (BIGNUM_NEGATIVE_P (arg2))
1421 ? bignum_posneg_bitwise_op(AND_OP, arg1, arg2)
1422 : bignum_pospos_bitwise_op(AND_OP, arg1, arg2)
1426 /* allocates memory */
1427 bignum *factorvm::bignum_bitwise_ior(bignum * arg1, bignum * arg2)
1430 (BIGNUM_NEGATIVE_P (arg1))
1431 ? (BIGNUM_NEGATIVE_P (arg2))
1432 ? bignum_negneg_bitwise_op(IOR_OP, arg1, arg2)
1433 : bignum_posneg_bitwise_op(IOR_OP, arg2, arg1)
1434 : (BIGNUM_NEGATIVE_P (arg2))
1435 ? bignum_posneg_bitwise_op(IOR_OP, arg1, arg2)
1436 : bignum_pospos_bitwise_op(IOR_OP, arg1, arg2)
1440 /* allocates memory */
1441 bignum *factorvm::bignum_bitwise_xor(bignum * arg1, bignum * arg2)
1444 (BIGNUM_NEGATIVE_P (arg1))
1445 ? (BIGNUM_NEGATIVE_P (arg2))
1446 ? bignum_negneg_bitwise_op(XOR_OP, arg1, arg2)
1447 : bignum_posneg_bitwise_op(XOR_OP, arg2, arg1)
1448 : (BIGNUM_NEGATIVE_P (arg2))
1449 ? bignum_posneg_bitwise_op(XOR_OP, arg1, arg2)
1450 : bignum_pospos_bitwise_op(XOR_OP, arg1, arg2)
1454 /* allocates memory */
1455 /* ash for the magnitude */
1456 /* assume arg1 is a big number, n is a long */
1457 bignum *factorvm::bignum_magnitude_ash(bignum * arg1, fixnum n)
1461 bignum * result = NULL;
1462 bignum_digit_type *scan1;
1463 bignum_digit_type *scanr;
1464 bignum_digit_type *end;
1466 fixnum digit_offset,bit_offset;
1468 if (BIGNUM_ZERO_P (arg1)) return (arg1);
1471 digit_offset = n / BIGNUM_DIGIT_LENGTH;
1472 bit_offset = n % BIGNUM_DIGIT_LENGTH;
1474 result = allot_bignum_zeroed (BIGNUM_LENGTH (arg1) + digit_offset + 1,
1475 BIGNUM_NEGATIVE_P(arg1));
1477 scanr = BIGNUM_START_PTR (result) + digit_offset;
1478 scan1 = BIGNUM_START_PTR (arg1);
1479 end = scan1 + BIGNUM_LENGTH (arg1);
1481 while (scan1 < end) {
1482 *scanr = *scanr | (*scan1 & BIGNUM_DIGIT_MASK) << bit_offset;
1483 *scanr = *scanr & BIGNUM_DIGIT_MASK;
1485 *scanr = *scan1++ >> (BIGNUM_DIGIT_LENGTH - bit_offset);
1486 *scanr = *scanr & BIGNUM_DIGIT_MASK;
1490 && (-n >= (BIGNUM_LENGTH (arg1) * (bignum_length_type) BIGNUM_DIGIT_LENGTH)))
1491 result = BIGNUM_ZERO ();
1494 digit_offset = -n / BIGNUM_DIGIT_LENGTH;
1495 bit_offset = -n % BIGNUM_DIGIT_LENGTH;
1497 result = allot_bignum_zeroed (BIGNUM_LENGTH (arg1) - digit_offset,
1498 BIGNUM_NEGATIVE_P(arg1));
1500 scanr = BIGNUM_START_PTR (result);
1501 scan1 = BIGNUM_START_PTR (arg1) + digit_offset;
1502 end = scanr + BIGNUM_LENGTH (result) - 1;
1504 while (scanr < end) {
1505 *scanr = (*scan1++ & BIGNUM_DIGIT_MASK) >> bit_offset ;
1507 *scan1 << (BIGNUM_DIGIT_LENGTH - bit_offset)) & BIGNUM_DIGIT_MASK;
1510 *scanr = (*scan1++ & BIGNUM_DIGIT_MASK) >> bit_offset ;
1512 else if (n == 0) result = arg1;
1514 return (bignum_trim (result));
1517 /* allocates memory */
1518 bignum *factorvm::bignum_pospos_bitwise_op(int op, bignum * arg1, bignum * arg2)
1520 GC_BIGNUM(arg1); GC_BIGNUM(arg2);
1523 bignum_length_type max_length;
1525 bignum_digit_type *scan1, *end1, digit1;
1526 bignum_digit_type *scan2, *end2, digit2;
1527 bignum_digit_type *scanr, *endr;
1529 max_length = (BIGNUM_LENGTH(arg1) > BIGNUM_LENGTH(arg2))
1530 ? BIGNUM_LENGTH(arg1) : BIGNUM_LENGTH(arg2);
1532 result = allot_bignum(max_length, 0);
1534 scanr = BIGNUM_START_PTR(result);
1535 scan1 = BIGNUM_START_PTR(arg1);
1536 scan2 = BIGNUM_START_PTR(arg2);
1537 endr = scanr + max_length;
1538 end1 = scan1 + BIGNUM_LENGTH(arg1);
1539 end2 = scan2 + BIGNUM_LENGTH(arg2);
1541 while (scanr < endr) {
1542 digit1 = (scan1 < end1) ? *scan1++ : 0;
1543 digit2 = (scan2 < end2) ? *scan2++ : 0;
1544 *scanr++ = (op == AND_OP) ? digit1 & digit2 :
1545 (op == IOR_OP) ? digit1 | digit2 :
1548 return bignum_trim(result);
1551 /* allocates memory */
1552 bignum *factorvm::bignum_posneg_bitwise_op(int op, bignum * arg1, bignum * arg2)
1554 GC_BIGNUM(arg1); GC_BIGNUM(arg2);
1557 bignum_length_type max_length;
1559 bignum_digit_type *scan1, *end1, digit1;
1560 bignum_digit_type *scan2, *end2, digit2, carry2;
1561 bignum_digit_type *scanr, *endr;
1563 char neg_p = op == IOR_OP || op == XOR_OP;
1565 max_length = (BIGNUM_LENGTH(arg1) > BIGNUM_LENGTH(arg2) + 1)
1566 ? BIGNUM_LENGTH(arg1) : BIGNUM_LENGTH(arg2) + 1;
1568 result = allot_bignum(max_length, neg_p);
1570 scanr = BIGNUM_START_PTR(result);
1571 scan1 = BIGNUM_START_PTR(arg1);
1572 scan2 = BIGNUM_START_PTR(arg2);
1573 endr = scanr + max_length;
1574 end1 = scan1 + BIGNUM_LENGTH(arg1);
1575 end2 = scan2 + BIGNUM_LENGTH(arg2);
1579 while (scanr < endr) {
1580 digit1 = (scan1 < end1) ? *scan1++ : 0;
1581 digit2 = (~((scan2 < end2) ? *scan2++ : 0) & BIGNUM_DIGIT_MASK)
1584 if (digit2 < BIGNUM_RADIX)
1588 digit2 = (digit2 - BIGNUM_RADIX);
1592 *scanr++ = (op == AND_OP) ? digit1 & digit2 :
1593 (op == IOR_OP) ? digit1 | digit2 :
1598 bignum_negate_magnitude(result);
1600 return bignum_trim(result);
1603 /* allocates memory */
1604 bignum *factorvm::bignum_negneg_bitwise_op(int op, bignum * arg1, bignum * arg2)
1606 GC_BIGNUM(arg1); GC_BIGNUM(arg2);
1609 bignum_length_type max_length;
1611 bignum_digit_type *scan1, *end1, digit1, carry1;
1612 bignum_digit_type *scan2, *end2, digit2, carry2;
1613 bignum_digit_type *scanr, *endr;
1615 char neg_p = op == AND_OP || op == IOR_OP;
1617 max_length = (BIGNUM_LENGTH(arg1) > BIGNUM_LENGTH(arg2))
1618 ? BIGNUM_LENGTH(arg1) + 1 : BIGNUM_LENGTH(arg2) + 1;
1620 result = allot_bignum(max_length, neg_p);
1622 scanr = BIGNUM_START_PTR(result);
1623 scan1 = BIGNUM_START_PTR(arg1);
1624 scan2 = BIGNUM_START_PTR(arg2);
1625 endr = scanr + max_length;
1626 end1 = scan1 + BIGNUM_LENGTH(arg1);
1627 end2 = scan2 + BIGNUM_LENGTH(arg2);
1632 while (scanr < endr) {
1633 digit1 = (~((scan1 < end1) ? *scan1++ : 0) & BIGNUM_DIGIT_MASK) + carry1;
1634 digit2 = (~((scan2 < end2) ? *scan2++ : 0) & BIGNUM_DIGIT_MASK) + carry2;
1636 if (digit1 < BIGNUM_RADIX)
1640 digit1 = (digit1 - BIGNUM_RADIX);
1644 if (digit2 < BIGNUM_RADIX)
1648 digit2 = (digit2 - BIGNUM_RADIX);
1652 *scanr++ = (op == AND_OP) ? digit1 & digit2 :
1653 (op == IOR_OP) ? digit1 | digit2 :
1658 bignum_negate_magnitude(result);
1660 return bignum_trim(result);
1663 void factorvm::bignum_negate_magnitude(bignum * arg)
1665 bignum_digit_type *scan;
1666 bignum_digit_type *end;
1667 bignum_digit_type digit;
1668 bignum_digit_type carry;
1670 scan = BIGNUM_START_PTR(arg);
1671 end = scan + BIGNUM_LENGTH(arg);
1675 while (scan < end) {
1676 digit = (~*scan & BIGNUM_DIGIT_MASK) + carry;
1678 if (digit < BIGNUM_RADIX)
1682 digit = (digit - BIGNUM_RADIX);
1690 /* Allocates memory */
1691 bignum *factorvm::bignum_integer_length(bignum * x)
1695 bignum_length_type index = ((BIGNUM_LENGTH (x)) - 1);
1696 bignum_digit_type digit = (BIGNUM_REF (x, index));
1698 bignum * result = (allot_bignum (2, 0));
1700 (BIGNUM_REF (result, 0)) = index;
1701 (BIGNUM_REF (result, 1)) = 0;
1702 bignum_destructive_scale_up (result, BIGNUM_DIGIT_LENGTH);
1705 bignum_destructive_add (result, ((bignum_digit_type) 1));
1708 return (bignum_trim (result));
1711 /* Allocates memory */
1712 int factorvm::bignum_logbitp(int shift, bignum * arg)
1714 return((BIGNUM_NEGATIVE_P (arg))
1715 ? !bignum_unsigned_logbitp (shift, bignum_bitwise_not (arg))
1716 : bignum_unsigned_logbitp (shift,arg));
1719 int factorvm::bignum_unsigned_logbitp(int shift, bignum * bignum)
1721 bignum_length_type len = (BIGNUM_LENGTH (bignum));
1722 int index = shift / BIGNUM_DIGIT_LENGTH;
1725 bignum_digit_type digit = (BIGNUM_REF (bignum, index));
1726 int p = shift % BIGNUM_DIGIT_LENGTH;
1727 bignum_digit_type mask = ((fixnum)1) << p;
1728 return (digit & mask) ? 1 : 0;
1731 /* Allocates memory */
1732 bignum *factorvm::digit_stream_to_bignum(unsigned int n_digits, unsigned int (*producer)(unsigned int, factorvm*), unsigned int radix, int negative_p)
1734 BIGNUM_ASSERT ((radix > 1) && (radix <= BIGNUM_RADIX_ROOT));
1736 return (BIGNUM_ZERO ());
1739 fixnum digit = ((fixnum) ((*producer) (0,this)));
1740 return (fixnum_to_bignum (negative_p ? (- digit) : digit));
1743 bignum_length_type length;
1745 unsigned int radix_copy = radix;
1746 unsigned int log_radix = 0;
1747 while (radix_copy > 0)
1752 /* This length will be at least as large as needed. */
1753 length = (BIGNUM_BITS_TO_DIGITS (n_digits * log_radix));
1756 bignum * result = (allot_bignum_zeroed (length, negative_p));
1757 while ((n_digits--) > 0)
1759 bignum_destructive_scale_up (result, ((bignum_digit_type) radix));
1760 bignum_destructive_add
1761 (result, ((bignum_digit_type) ((*producer) (n_digits,this))));
1763 return (bignum_trim (result));