steve
/
factor
1
0
Fork 0
Code Issues Pull Requests Releases Wiki Activity

vm: use tabs instead of spaces in bignum_gcd.

db4
John Benediktsson 2012-08-22 21:31:27 -07:00
parent 54fd95cd81
commit b541986b4e
1 changed files with 188 additions and 188 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

376
vm/bignum.cpp
View File

@ -1720,225 +1720,225 @@ int factor_vm::bignum_unsigned_logbitp(int shift, bignum * bignum)
/* Allocates memory */ /* Allocates memory */
bignum * factor_vm::bignum_gcd(bignum * a, bignum * b) bignum * factor_vm::bignum_gcd(bignum * a, bignum * b)
{ {
GC_BIGNUM(a); GC_BIGNUM(a);
GC_BIGNUM(b); GC_BIGNUM(b);
bignum * d; bignum * d;
bignum_length_type size_a, size_b; bignum_length_type size_a, size_b;
bignum_digit_type * scan_a, * scan_b, * scan_d, * a_end, * b_end; bignum_digit_type * scan_a, * scan_b, * scan_d, * a_end, * b_end;
if (BIGNUM_NEGATIVE_P (a)) { if (BIGNUM_NEGATIVE_P (a)) {
size_a = BIGNUM_LENGTH (a); size_a = BIGNUM_LENGTH (a);
d = allot_bignum (size_a, 0); d = allot_bignum (size_a, 0);
scan_d = BIGNUM_START_PTR (d); scan_d = BIGNUM_START_PTR (d);
scan_a = BIGNUM_START_PTR (a); scan_a = BIGNUM_START_PTR (a);
a_end = scan_a + size_a; a_end = scan_a + size_a;
while (scan_a < a_end) while (scan_a < a_end)
(*scan_d++) = (*scan_a++); (*scan_d++) = (*scan_a++);
a = d; a = d;
} }
if (BIGNUM_NEGATIVE_P (b)) { if (BIGNUM_NEGATIVE_P (b)) {
size_b = BIGNUM_LENGTH (b); size_b = BIGNUM_LENGTH (b);
d = allot_bignum (size_b, 0); d = allot_bignum (size_b, 0);
scan_d = BIGNUM_START_PTR (d); scan_d = BIGNUM_START_PTR (d);
scan_b = BIGNUM_START_PTR (b); scan_b = BIGNUM_START_PTR (b);
b_end = scan_b + size_b; b_end = scan_b + size_b;
while (scan_b < b_end) while (scan_b < b_end)
(*scan_d++) = (*scan_b++); (*scan_d++) = (*scan_b++);
b = d; b = d;
} }
if (bignum_compare(a, b) == bignum_comparison_less) { if (bignum_compare(a, b) == bignum_comparison_less) {
std::swap(a, b); std::swap(a, b);
} }
while (BIGNUM_LENGTH (b) != 0) { while (BIGNUM_LENGTH (b) != 0) {
d = bignum_remainder (a, b); d = bignum_remainder (a, b);
GC_BIGNUM(d); GC_BIGNUM(d);
if (d == BIGNUM_OUT_OF_BAND) { if (d == BIGNUM_OUT_OF_BAND) {
return d; return d;
} }
a = b; a = b;
b = d; b = d;
} }
return a; return a;
} }
#else #else
/* Allocates memory */ /* Allocates memory */
bignum * factor_vm::bignum_gcd(bignum * a, bignum * b) bignum * factor_vm::bignum_gcd(bignum * a, bignum * b)
{ {
GC_BIGNUM(a); GC_BIGNUM(a);
GC_BIGNUM(b); GC_BIGNUM(b);
bignum *c, *d, *e, *f; bignum *c, *d, *e, *f;
bignum_twodigit_type x, y, q, s, t, A, B, C, D; bignum_twodigit_type x, y, q, s, t, A, B, C, D;
int nbits, k; int nbits, k;
bignum_length_type size_a, size_b, size_c; bignum_length_type size_a, size_b, size_c;
bignum_digit_type *scan_a, *scan_b, *scan_c, *scan_d; bignum_digit_type *scan_a, *scan_b, *scan_c, *scan_d;
bignum_digit_type *a_end, *b_end, *c_end; bignum_digit_type *a_end, *b_end, *c_end;
/* clone the bignums so we can modify them in-place */ /* clone the bignums so we can modify them in-place */
size_a = BIGNUM_LENGTH (a); size_a = BIGNUM_LENGTH (a);
c = allot_bignum (size_a, 0); c = allot_bignum (size_a, 0);
scan_a = BIGNUM_START_PTR (a); scan_a = BIGNUM_START_PTR (a);
a_end = scan_a + size_a; a_end = scan_a + size_a;
GC_BIGNUM(c); GC_BIGNUM(c);
scan_c = BIGNUM_START_PTR (c); scan_c = BIGNUM_START_PTR (c);
while (scan_a < a_end) while (scan_a < a_end)
(*scan_c++) = (*scan_a++); (*scan_c++) = (*scan_a++);
a = c; a = c;
size_b = BIGNUM_LENGTH (b); size_b = BIGNUM_LENGTH (b);
d = allot_bignum (size_b, 0); d = allot_bignum (size_b, 0);
scan_b = BIGNUM_START_PTR (b); scan_b = BIGNUM_START_PTR (b);
b_end = scan_b + size_b; b_end = scan_b + size_b;
GC_BIGNUM(d); GC_BIGNUM(d);
scan_d = BIGNUM_START_PTR (d); scan_d = BIGNUM_START_PTR (d);
while (scan_b < b_end) while (scan_b < b_end)
(*scan_d++) = (*scan_b++); (*scan_d++) = (*scan_b++);
b = d; b = d;
/* Initial reduction: make sure that 0 <= b <= a. */ /* Initial reduction: make sure that 0 <= b <= a. */
if (bignum_compare(a, b) == bignum_comparison_less) { if (bignum_compare(a, b) == bignum_comparison_less) {
std::swap(a, b); std::swap(a, b);
std::swap(size_a, size_b); std::swap(size_a, size_b);
} }
while (size_a > 1) { while (size_a > 1) {
nbits = log2 (BIGNUM_REF (a, size_a-1)); nbits = log2 (BIGNUM_REF (a, size_a-1));
x = ((BIGNUM_REF (a, size_a-1) << (BIGNUM_DIGIT_LENGTH - nbits)) | x = ((BIGNUM_REF (a, size_a-1) << (BIGNUM_DIGIT_LENGTH - nbits)) |
(BIGNUM_REF (a, size_a-2) >> nbits)); (BIGNUM_REF (a, size_a-2) >> nbits));
y = ((size_b >= size_a - 1 ? BIGNUM_REF (b, size_a-2) >> nbits : 0) | 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)); (size_b >= size_a ? BIGNUM_REF (b, size_a-1) << (BIGNUM_DIGIT_LENGTH - nbits) : 0));
/* inner loop of Lehmer's algorithm; */ /* inner loop of Lehmer's algorithm; */
A = 1; B = 0; C = 0; D = 1; A = 1; B = 0; C = 0; D = 1;
for (k=0 ;; k++) { for (k=0 ;; k++) {
if (y - C == 0) if (y - C == 0)
break; break;
q = (x + (A - 1)) / (y - C); q = (x + (A - 1)) / (y - C);
s = B + (q * D); s = B + (q * D);
t = x - (q * y); t = x - (q * y);
if (s > t) if (s > t)
break; break;
x = y; x = y;
y = t; y = t;
t = A + (q * C); t = A + (q * C);
A = D; B = C; C = s; D = t; A = D; B = C; C = s; D = t;
} }
if (k == 0) { if (k == 0) {
/* no progress; do a Euclidean step */ /* no progress; do a Euclidean step */
if (size_b == 0) { if (size_b == 0) {
return bignum_trim (a); return bignum_trim (a);
} }
e = bignum_trim (a); e = bignum_trim (a);
GC_BIGNUM(e); GC_BIGNUM(e);
f = bignum_trim (b); f = bignum_trim (b);
GC_BIGNUM(f); GC_BIGNUM(f);
c = bignum_remainder (e, f); c = bignum_remainder (e, f);
GC_BIGNUM (c); GC_BIGNUM (c);
if (c == BIGNUM_OUT_OF_BAND) { if (c == BIGNUM_OUT_OF_BAND) {
return c; return c;
} }
// copy 'b' to 'a' // copy 'b' to 'a'
scan_a = BIGNUM_START_PTR (a); scan_a = BIGNUM_START_PTR (a);
scan_b = BIGNUM_START_PTR (b); scan_b = BIGNUM_START_PTR (b);
a_end = scan_a + size_a; a_end = scan_a + size_a;
b_end = scan_b + size_b; b_end = scan_b + size_b;
while (scan_b < b_end) *(scan_a++) = *(scan_b++); while (scan_b < b_end) *(scan_a++) = *(scan_b++);
while (scan_a < a_end) *(scan_a++) = 0; while (scan_a < a_end) *(scan_a++) = 0;
size_a = size_b; size_a = size_b;
// copy 'c' to 'b' // copy 'c' to 'b'
scan_b = BIGNUM_START_PTR (b); scan_b = BIGNUM_START_PTR (b);
scan_c = BIGNUM_START_PTR (c); scan_c = BIGNUM_START_PTR (c);
size_c = BIGNUM_LENGTH (c); size_c = BIGNUM_LENGTH (c);
c_end = scan_c + size_c; c_end = scan_c + size_c;
while (scan_c < c_end) *(scan_b++) = *(scan_c++); while (scan_c < c_end) *(scan_b++) = *(scan_c++);
while (scan_b < b_end) *(scan_b++) = 0; while (scan_b < b_end) *(scan_b++) = 0;
size_b = size_c; size_b = size_c;
continue; continue;
} }
/* /*
a, b = A*b - B*a, D*a - C*b if k is odd a, b = A*b - B*a, D*a - C*b if k is odd
a, b = A*a - B*b, D*b - C*a if k is even a, b = A*a - B*b, D*b - C*a if k is even
*/ */
scan_a = BIGNUM_START_PTR (a); scan_a = BIGNUM_START_PTR (a);
scan_b = BIGNUM_START_PTR (b); scan_b = BIGNUM_START_PTR (b);
scan_c = scan_a; scan_c = scan_a;
scan_d = scan_b; scan_d = scan_b;
a_end = scan_a + size_a; a_end = scan_a + size_a;
b_end = scan_b + size_b; b_end = scan_b + size_b;
s = 0; s = 0;
t = 0; t = 0;
if (k & 1) { if (k & 1) {
while (scan_b < b_end) { while (scan_b < b_end) {
s += (A * *scan_b) - (B * *scan_a); s += (A * *scan_b) - (B * *scan_a);
t += (D * *scan_a++) - (C * *scan_b++); t += (D * *scan_a++) - (C * *scan_b++);
*scan_c++ = (bignum_digit_type) (s & BIGNUM_DIGIT_MASK); *scan_c++ = (bignum_digit_type) (s & BIGNUM_DIGIT_MASK);
*scan_d++ = (bignum_digit_type) (t & BIGNUM_DIGIT_MASK); *scan_d++ = (bignum_digit_type) (t & BIGNUM_DIGIT_MASK);
s >>= BIGNUM_DIGIT_LENGTH; s >>= BIGNUM_DIGIT_LENGTH;
t >>= BIGNUM_DIGIT_LENGTH; t >>= BIGNUM_DIGIT_LENGTH;
} }
while (scan_a < a_end) { while (scan_a < a_end) {
s -= (B * *scan_a); s -= (B * *scan_a);
t += (D * *scan_a++); t += (D * *scan_a++);
*scan_c++ = (bignum_digit_type) (s & BIGNUM_DIGIT_MASK); *scan_c++ = (bignum_digit_type) (s & BIGNUM_DIGIT_MASK);
//*scan_d++ = (bignum_digit_type) (t & BIGNUM_DIGIT_MASK); //*scan_d++ = (bignum_digit_type) (t & BIGNUM_DIGIT_MASK);
s >>= BIGNUM_DIGIT_LENGTH; s >>= BIGNUM_DIGIT_LENGTH;
t >>= BIGNUM_DIGIT_LENGTH; t >>= BIGNUM_DIGIT_LENGTH;
} }
} }
else { else {
while (scan_b < b_end) { while (scan_b < b_end) {
s += (A * *scan_a) - (B * *scan_b); s += (A * *scan_a) - (B * *scan_b);
t += (D * *scan_b++) - (C * *scan_a++); t += (D * *scan_b++) - (C * *scan_a++);
*scan_c++ = (bignum_digit_type) (s & BIGNUM_DIGIT_MASK); *scan_c++ = (bignum_digit_type) (s & BIGNUM_DIGIT_MASK);
*scan_d++ = (bignum_digit_type) (t & BIGNUM_DIGIT_MASK); *scan_d++ = (bignum_digit_type) (t & BIGNUM_DIGIT_MASK);
s >>= BIGNUM_DIGIT_LENGTH; s >>= BIGNUM_DIGIT_LENGTH;
t >>= BIGNUM_DIGIT_LENGTH; t >>= BIGNUM_DIGIT_LENGTH;
} }
while (scan_a < a_end) { while (scan_a < a_end) {
s += (A * *scan_a); s += (A * *scan_a);
t -= (C * *scan_a++); t -= (C * *scan_a++);
*scan_c++ = (bignum_digit_type) (s & BIGNUM_DIGIT_MASK); *scan_c++ = (bignum_digit_type) (s & BIGNUM_DIGIT_MASK);
//*scan_d++ = (bignum_digit_type) (t & BIGNUM_DIGIT_MASK); //*scan_d++ = (bignum_digit_type) (t & BIGNUM_DIGIT_MASK);
s >>= BIGNUM_DIGIT_LENGTH; s >>= BIGNUM_DIGIT_LENGTH;
t >>= BIGNUM_DIGIT_LENGTH; t >>= BIGNUM_DIGIT_LENGTH;
} }
} }
BIGNUM_ASSERT (s == 0); BIGNUM_ASSERT (s == 0);
BIGNUM_ASSERT (t == 0); BIGNUM_ASSERT (t == 0);
// update size_a and size_b to remove any zeros at end // update size_a and size_b to remove any zeros at end
while (size_a > 0 && *(--scan_a) == 0) size_a--; while (size_a > 0 && *(--scan_a) == 0) size_a--;
while (size_b > 0 && *(--scan_b) == 0) size_b--; while (size_b > 0 && *(--scan_b) == 0) size_b--;
BIGNUM_ASSERT (size_a >= size_b); BIGNUM_ASSERT (size_a >= size_b);
} }
/* a fits into a fixnum, so b must too */ /* a fits into a fixnum, so b must too */
fixnum xx = bignum_to_fixnum (a); fixnum xx = bignum_to_fixnum (a);
fixnum yy = bignum_to_fixnum (b); fixnum yy = bignum_to_fixnum (b);
fixnum tt; fixnum tt;
/* usual Euclidean algorithm for longs */ /* usual Euclidean algorithm for longs */
while (yy != 0) { while (yy != 0) {
tt = yy; tt = yy;
yy = xx % yy; yy = xx % yy;
xx = tt; xx = tt;
} }
return fixnum_to_bignum (xx); return fixnum_to_bignum (xx);
} }
#endif #endif