annotate bn_mp_div.c @ 145:a96ff234ff19 libtommath

propagate of fc94c38452d9fd684a8e1eb9e3a73120aac0d38f and 3e4de4cbef3d9035a7b2f0c25e9f86e297f9f6d1 from branch 'au.asn.ucc.matt.ltm-orig' to 'au.asn.ucc.matt.ltm-db'
author Matt Johnston <matt@ucc.asn.au>
date Sun, 19 Dec 2004 15:57:19 +0000
parents d29b64170cf0
children d8254fc979e9
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
142
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
1 #include <tommath.h>
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
2 #ifdef BN_MP_DIV_C
2
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
3 /* LibTomMath, multiple-precision integer library -- Tom St Denis
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
4 *
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
5 * LibTomMath is a library that provides multiple-precision
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
6 * integer arithmetic as well as number theoretic functionality.
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
7 *
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
8 * The library was designed directly after the MPI library by
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
9 * Michael Fromberger but has been written from scratch with
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
10 * additional optimizations in place.
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
11 *
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
12 * The library is free for all purposes without any express
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
13 * guarantee it works.
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
14 *
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
15 * Tom St Denis, [email protected], http://math.libtomcrypt.org
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
16 */
142
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
17
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
18 #ifdef BN_MP_DIV_SMALL
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
19
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
20 /* slower bit-bang division... also smaller */
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
21 int mp_div(mp_int * a, mp_int * b, mp_int * c, mp_int * d)
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
22 {
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
23 mp_int ta, tb, tq, q;
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
24 int res, n, n2;
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
25
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
26 /* is divisor zero ? */
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
27 if (mp_iszero (b) == 1) {
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
28 return MP_VAL;
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
29 }
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
30
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
31 /* if a < b then q=0, r = a */
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
32 if (mp_cmp_mag (a, b) == MP_LT) {
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
33 if (d != NULL) {
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
34 res = mp_copy (a, d);
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
35 } else {
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
36 res = MP_OKAY;
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
37 }
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
38 if (c != NULL) {
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
39 mp_zero (c);
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
40 }
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
41 return res;
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
42 }
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
43
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
44 /* init our temps */
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
45 if ((res = mp_init_multi(&ta, &tb, &tq, &q, NULL) != MP_OKAY)) {
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
46 return res;
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
47 }
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
48
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
49
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
50 mp_set(&tq, 1);
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
51 n = mp_count_bits(a) - mp_count_bits(b);
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
52 if (((res = mp_copy(a, &ta)) != MP_OKAY) ||
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
53 ((res = mp_copy(b, &tb)) != MP_OKAY) ||
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
54 ((res = mp_mul_2d(&tb, n, &tb)) != MP_OKAY) ||
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
55 ((res = mp_mul_2d(&tq, n, &tq)) != MP_OKAY)) {
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
56 goto __ERR;
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
57 }
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
58
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
59 while (n-- >= 0) {
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
60 if (mp_cmp(&tb, &ta) != MP_GT) {
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
61 if (((res = mp_sub(&ta, &tb, &ta)) != MP_OKAY) ||
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
62 ((res = mp_add(&q, &tq, &q)) != MP_OKAY)) {
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
63 goto __ERR;
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
64 }
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
65 }
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
66 if (((res = mp_div_2d(&tb, 1, &tb, NULL)) != MP_OKAY) ||
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
67 ((res = mp_div_2d(&tq, 1, &tq, NULL)) != MP_OKAY)) {
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
68 goto __ERR;
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
69 }
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
70 }
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
71
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
72 /* now q == quotient and ta == remainder */
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
73 n = a->sign;
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
74 n2 = (a->sign == b->sign ? MP_ZPOS : MP_NEG);
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
75 if (c != NULL) {
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
76 mp_exch(c, &q);
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
77 c->sign = n2;
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
78 }
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
79 if (d != NULL) {
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
80 mp_exch(d, &ta);
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
81 d->sign = n;
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
82 }
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
83 __ERR:
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
84 mp_clear_multi(&ta, &tb, &tq, &q, NULL);
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
85 return res;
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
86 }
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
87
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
88 #else
2
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
89
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
90 /* integer signed division.
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
91 * c*b + d == a [e.g. a/b, c=quotient, d=remainder]
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
92 * HAC pp.598 Algorithm 14.20
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
93 *
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
94 * Note that the description in HAC is horribly
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
95 * incomplete. For example, it doesn't consider
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
96 * the case where digits are removed from 'x' in
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
97 * the inner loop. It also doesn't consider the
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
98 * case that y has fewer than three digits, etc..
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
99 *
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
100 * The overall algorithm is as described as
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
101 * 14.20 from HAC but fixed to treat these cases.
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
102 */
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
103 int mp_div (mp_int * a, mp_int * b, mp_int * c, mp_int * d)
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
104 {
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
105 mp_int q, x, y, t1, t2;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
106 int res, n, t, i, norm, neg;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
107
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
108 /* is divisor zero ? */
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
109 if (mp_iszero (b) == 1) {
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
110 return MP_VAL;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
111 }
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
112
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
113 /* if a < b then q=0, r = a */
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
114 if (mp_cmp_mag (a, b) == MP_LT) {
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
115 if (d != NULL) {
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
116 res = mp_copy (a, d);
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
117 } else {
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
118 res = MP_OKAY;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
119 }
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
120 if (c != NULL) {
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
121 mp_zero (c);
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
122 }
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
123 return res;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
124 }
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
125
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
126 if ((res = mp_init_size (&q, a->used + 2)) != MP_OKAY) {
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
127 return res;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
128 }
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
129 q.used = a->used + 2;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
130
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
131 if ((res = mp_init (&t1)) != MP_OKAY) {
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
132 goto __Q;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
133 }
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
134
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
135 if ((res = mp_init (&t2)) != MP_OKAY) {
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
136 goto __T1;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
137 }
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
138
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
139 if ((res = mp_init_copy (&x, a)) != MP_OKAY) {
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
140 goto __T2;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
141 }
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
142
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
143 if ((res = mp_init_copy (&y, b)) != MP_OKAY) {
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
144 goto __X;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
145 }
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
146
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
147 /* fix the sign */
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
148 neg = (a->sign == b->sign) ? MP_ZPOS : MP_NEG;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
149 x.sign = y.sign = MP_ZPOS;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
150
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
151 /* normalize both x and y, ensure that y >= b/2, [b == 2**DIGIT_BIT] */
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
152 norm = mp_count_bits(&y) % DIGIT_BIT;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
153 if (norm < (int)(DIGIT_BIT-1)) {
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
154 norm = (DIGIT_BIT-1) - norm;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
155 if ((res = mp_mul_2d (&x, norm, &x)) != MP_OKAY) {
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
156 goto __Y;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
157 }
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
158 if ((res = mp_mul_2d (&y, norm, &y)) != MP_OKAY) {
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
159 goto __Y;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
160 }
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
161 } else {
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
162 norm = 0;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
163 }
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
164
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
165 /* note hac does 0 based, so if used==5 then its 0,1,2,3,4, e.g. use 4 */
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
166 n = x.used - 1;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
167 t = y.used - 1;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
168
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
169 /* while (x >= y*b**n-t) do { q[n-t] += 1; x -= y*b**{n-t} } */
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
170 if ((res = mp_lshd (&y, n - t)) != MP_OKAY) { /* y = y*b**{n-t} */
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
171 goto __Y;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
172 }
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
173
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
174 while (mp_cmp (&x, &y) != MP_LT) {
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
175 ++(q.dp[n - t]);
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
176 if ((res = mp_sub (&x, &y, &x)) != MP_OKAY) {
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
177 goto __Y;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
178 }
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
179 }
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
180
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
181 /* reset y by shifting it back down */
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
182 mp_rshd (&y, n - t);
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
183
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
184 /* step 3. for i from n down to (t + 1) */
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
185 for (i = n; i >= (t + 1); i--) {
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
186 if (i > x.used) {
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
187 continue;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
188 }
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
189
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
190 /* step 3.1 if xi == yt then set q{i-t-1} to b-1,
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
191 * otherwise set q{i-t-1} to (xi*b + x{i-1})/yt */
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
192 if (x.dp[i] == y.dp[t]) {
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
193 q.dp[i - t - 1] = ((((mp_digit)1) << DIGIT_BIT) - 1);
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
194 } else {
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
195 mp_word tmp;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
196 tmp = ((mp_word) x.dp[i]) << ((mp_word) DIGIT_BIT);
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
197 tmp |= ((mp_word) x.dp[i - 1]);
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
198 tmp /= ((mp_word) y.dp[t]);
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
199 if (tmp > (mp_word) MP_MASK)
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
200 tmp = MP_MASK;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
201 q.dp[i - t - 1] = (mp_digit) (tmp & (mp_word) (MP_MASK));
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
202 }
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
203
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
204 /* while (q{i-t-1} * (yt * b + y{t-1})) >
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
205 xi * b**2 + xi-1 * b + xi-2
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
206
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
207 do q{i-t-1} -= 1;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
208 */
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
209 q.dp[i - t - 1] = (q.dp[i - t - 1] + 1) & MP_MASK;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
210 do {
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
211 q.dp[i - t - 1] = (q.dp[i - t - 1] - 1) & MP_MASK;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
212
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
213 /* find left hand */
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
214 mp_zero (&t1);
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
215 t1.dp[0] = (t - 1 < 0) ? 0 : y.dp[t - 1];
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
216 t1.dp[1] = y.dp[t];
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
217 t1.used = 2;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
218 if ((res = mp_mul_d (&t1, q.dp[i - t - 1], &t1)) != MP_OKAY) {
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
219 goto __Y;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
220 }
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
221
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
222 /* find right hand */
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
223 t2.dp[0] = (i - 2 < 0) ? 0 : x.dp[i - 2];
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
224 t2.dp[1] = (i - 1 < 0) ? 0 : x.dp[i - 1];
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
225 t2.dp[2] = x.dp[i];
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
226 t2.used = 3;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
227 } while (mp_cmp_mag(&t1, &t2) == MP_GT);
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
228
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
229 /* step 3.3 x = x - q{i-t-1} * y * b**{i-t-1} */
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
230 if ((res = mp_mul_d (&y, q.dp[i - t - 1], &t1)) != MP_OKAY) {
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
231 goto __Y;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
232 }
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
233
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
234 if ((res = mp_lshd (&t1, i - t - 1)) != MP_OKAY) {
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
235 goto __Y;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
236 }
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
237
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
238 if ((res = mp_sub (&x, &t1, &x)) != MP_OKAY) {
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
239 goto __Y;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
240 }
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
241
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
242 /* if x < 0 then { x = x + y*b**{i-t-1}; q{i-t-1} -= 1; } */
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
243 if (x.sign == MP_NEG) {
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
244 if ((res = mp_copy (&y, &t1)) != MP_OKAY) {
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
245 goto __Y;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
246 }
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
247 if ((res = mp_lshd (&t1, i - t - 1)) != MP_OKAY) {
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
248 goto __Y;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
249 }
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
250 if ((res = mp_add (&x, &t1, &x)) != MP_OKAY) {
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
251 goto __Y;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
252 }
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
253
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
254 q.dp[i - t - 1] = (q.dp[i - t - 1] - 1UL) & MP_MASK;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
255 }
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
256 }
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
257
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
258 /* now q is the quotient and x is the remainder
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
259 * [which we have to normalize]
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
260 */
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
261
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
262 /* get sign before writing to c */
142
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
263 x.sign = x.used == 0 ? MP_ZPOS : a->sign;
2
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
264
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
265 if (c != NULL) {
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
266 mp_clamp (&q);
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
267 mp_exch (&q, c);
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
268 c->sign = neg;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
269 }
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
270
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
271 if (d != NULL) {
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
272 mp_div_2d (&x, norm, &x, NULL);
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
273 mp_exch (&x, d);
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
274 }
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
275
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
276 res = MP_OKAY;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
277
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
278 __Y:mp_clear (&y);
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
279 __X:mp_clear (&x);
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
280 __T2:mp_clear (&t2);
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
281 __T1:mp_clear (&t1);
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
282 __Q:mp_clear (&q);
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
283 return res;
86e0b50a9b58 ltm 0.30 orig import
Matt Johnston <matt@ucc.asn.au>
parents:
diff changeset
284 }
142
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
285
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
286 #endif
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
287
d29b64170cf0 import of libtommath 0.32
Matt Johnston <matt@ucc.asn.au>
parents: 2
diff changeset
288 #endif