-bignum * factor_vm::bignum_gcd(bignum * a, bignum * b)
-{
- GC_BIGNUM(a);
- GC_BIGNUM(b);
- bignum *c, *d, *e, *f;
- bignum_twodigit_type x, y, q, s, t, A, B, C, D;
- int nbits, k;
- bignum_length_type size_a, size_b, size_c;
- bignum_digit_type *scan_a, *scan_b, *scan_c, *scan_d;
- bignum_digit_type *a_end, *b_end, *c_end;
-
- /* clone the bignums so we can modify them in-place */
- size_a = BIGNUM_LENGTH (a);
- c = allot_bignum (size_a, 0);
- scan_a = BIGNUM_START_PTR (a);
- a_end = scan_a + size_a;
- GC_BIGNUM(c);
- scan_c = BIGNUM_START_PTR (c);
- while (scan_a < a_end)
- (*scan_c++) = (*scan_a++);
- a = c;
- size_b = BIGNUM_LENGTH (b);
- d = allot_bignum (size_b, 0);
- scan_b = BIGNUM_START_PTR (b);
- b_end = scan_b + size_b;
- GC_BIGNUM(d);
- scan_d = BIGNUM_START_PTR (d);
- while (scan_b < b_end)
- (*scan_d++) = (*scan_b++);
- b = d;
-
- /* Initial reduction: make sure that 0 <= b <= a. */
- if (bignum_compare(a, b) == bignum_comparison_less) {
- std::swap(a, b);
- std::swap(size_a, size_b);
- }
-
- while (size_a > 1) {
- nbits = log2 (BIGNUM_REF (a, size_a-1));
- x = ((BIGNUM_REF (a, size_a-1) << (BIGNUM_DIGIT_LENGTH - nbits)) |
- (BIGNUM_REF (a, size_a-2) >> nbits));
- y = ((size_b >= size_a - 1 ? BIGNUM_REF (b, size_a-2) >> nbits : 0) |
- (size_b >= size_a ? BIGNUM_REF (b, size_a-1) << (BIGNUM_DIGIT_LENGTH - nbits) : 0));
-
- /* inner loop of Lehmer's algorithm; */
- A = 1; B = 0; C = 0; D = 1;
- for (k=0 ;; k++) {
- if (y - C == 0)
- break;
-
- q = (x + (A - 1)) / (y - C);
-
- s = B + (q * D);
- t = x - (q * y);
-
- if (s > t)
- break;
-
- x = y;
- y = t;
-
- t = A + (q * C);
-
- A = D; B = C; C = s; D = t;
- }
-
- if (k == 0) {
- /* no progress; do a Euclidean step */
- if (size_b == 0) {
- return bignum_trim (a);
- }
- e = bignum_trim (a);
- GC_BIGNUM(e);
- f = bignum_trim (b);
- GC_BIGNUM(f);
- c = bignum_remainder (e, f);
- GC_BIGNUM (c);
- if (c == BIGNUM_OUT_OF_BAND) {
- return c;
- }
-
- // copy 'b' to 'a'
- scan_a = BIGNUM_START_PTR (a);
- scan_b = BIGNUM_START_PTR (b);
- a_end = scan_a + size_a;
- b_end = scan_b + size_b;
- while (scan_b < b_end) *(scan_a++) = *(scan_b++);
- while (scan_a < a_end) *(scan_a++) = 0;
- size_a = size_b;
-
- // copy 'c' to 'b'
- scan_b = BIGNUM_START_PTR (b);
- scan_c = BIGNUM_START_PTR (c);
- size_c = BIGNUM_LENGTH (c);
- c_end = scan_c + size_c;
- while (scan_c < c_end) *(scan_b++) = *(scan_c++);
- while (scan_b < b_end) *(scan_b++) = 0;
- size_b = size_c;
-
- continue;
- }
-
- /*
+bignum* factor_vm::bignum_gcd(bignum* a, bignum* b) {
+ GC_BIGNUM(a);
+ GC_BIGNUM(b);
+ bignum* c, *d, *e, *f;
+ bignum_twodigit_type x, y, q, s, t, A, B, C, D;
+ int nbits, k;
+ bignum_length_type size_a, size_b, size_c;
+ bignum_digit_type* scan_a, *scan_b, *scan_c, *scan_d;
+ bignum_digit_type* a_end, *b_end, *c_end;
+
+ /* clone the bignums so we can modify them in-place */
+ size_a = BIGNUM_LENGTH(a);
+ c = allot_bignum(size_a, 0);
+ scan_a = BIGNUM_START_PTR(a);
+ a_end = scan_a + size_a;
+ GC_BIGNUM(c);
+ scan_c = BIGNUM_START_PTR(c);
+ while (scan_a < a_end)
+ (*scan_c++) = (*scan_a++);
+ a = c;
+ size_b = BIGNUM_LENGTH(b);
+ d = allot_bignum(size_b, 0);
+ scan_b = BIGNUM_START_PTR(b);
+ b_end = scan_b + size_b;
+ GC_BIGNUM(d);
+ scan_d = BIGNUM_START_PTR(d);
+ while (scan_b < b_end)
+ (*scan_d++) = (*scan_b++);
+ b = d;
+
+ /* Initial reduction: make sure that 0 <= b <= a. */
+ if (bignum_compare(a, b) == bignum_comparison_less) {
+ std::swap(a, b);
+ std::swap(size_a, size_b);
+ }
+
+ while (size_a > 1) {
+ nbits = log2(BIGNUM_REF(a, size_a - 1));
+ x = ((BIGNUM_REF(a, size_a - 1) << (BIGNUM_DIGIT_LENGTH - nbits)) |
+ (BIGNUM_REF(a, size_a - 2) >> nbits));
+ y = ((size_b >= size_a - 1 ? BIGNUM_REF(b, size_a - 2) >> nbits : 0) |
+ (size_b >= size_a
+ ? BIGNUM_REF(b, size_a - 1) << (BIGNUM_DIGIT_LENGTH - nbits)
+ : 0));
+
+ /* inner loop of Lehmer's algorithm; */
+ A = 1;
+ B = 0;
+ C = 0;
+ D = 1;
+ for (k = 0;; k++) {
+ if (y - C == 0)
+ break;
+
+ q = (x + (A - 1)) / (y - C);
+
+ s = B + (q * D);
+ t = x - (q * y);
+
+ if (s > t)
+ break;
+
+ x = y;
+ y = t;
+
+ t = A + (q * C);
+
+ A = D;
+ B = C;
+ C = s;
+ D = t;
+ }
+
+ if (k == 0) {
+ /* no progress; do a Euclidean step */
+ if (size_b == 0) {
+ return bignum_trim(a);
+ }
+ e = bignum_trim(a);
+ GC_BIGNUM(e);
+ f = bignum_trim(b);
+ GC_BIGNUM(f);
+ c = bignum_remainder(e, f);
+ GC_BIGNUM(c);
+ if (c == BIGNUM_OUT_OF_BAND) {
+ return c;
+ }
+
+ // copy 'b' to 'a'
+ scan_a = BIGNUM_START_PTR(a);
+ scan_b = BIGNUM_START_PTR(b);
+ a_end = scan_a + size_a;
+ b_end = scan_b + size_b;
+ while (scan_b < b_end)
+ *(scan_a++) = *(scan_b++);
+ while (scan_a < a_end)
+ *(scan_a++) = 0;
+ size_a = size_b;
+
+ // copy 'c' to 'b'
+ scan_b = BIGNUM_START_PTR(b);
+ scan_c = BIGNUM_START_PTR(c);
+ size_c = BIGNUM_LENGTH(c);
+ c_end = scan_c + size_c;
+ while (scan_c < c_end)
+ *(scan_b++) = *(scan_c++);
+ while (scan_b < b_end)
+ *(scan_b++) = 0;
+ size_b = size_c;
+
+ continue;
+ }
+
+ /*