Mercurial > dropbear

comparison libtommath/bn_fast_mp_montgomery_reduce.c @ 389:5ff8218bcee9

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
propagate from branch 'au.asn.ucc.matt.ltm.dropbear' (head 2af95f00ebd5bb7a28b3817db1218442c935388e) to branch 'au.asn.ucc.matt.dropbear' (head ecd779509ef23a8cdf64888904fc9b31d78aa933)
author Matt Johnston <matt@ucc.asn.au>
date Thu, 11 Jan 2007 03:14:55 +0000
parents eed26cff980b
children 60fc6476e044
comparison
equal deleted inserted replaced
388:fb54020f78e1 389:5ff8218bcee9
1 #include <tommath.h>
2 #ifdef BN_FAST_MP_MONTGOMERY_REDUCE_C
3 /* LibTomMath, multiple-precision integer library -- Tom St Denis
4 *
5 * LibTomMath is a library that provides multiple-precision
6 * integer arithmetic as well as number theoretic functionality.
7 *
8 * The library was designed directly after the MPI library by
9 * Michael Fromberger but has been written from scratch with
10 * additional optimizations in place.
11 *
12 * The library is free for all purposes without any express
13 * guarantee it works.
14 *
15 * Tom St Denis, [email protected], http://math.libtomcrypt.com
16 */
17
18 /* computes xR**-1 == x (mod N) via Montgomery Reduction
19 *
20 * This is an optimized implementation of montgomery_reduce
21 * which uses the comba method to quickly calculate the columns of the
22 * reduction.
23 *
24 * Based on Algorithm 14.32 on pp.601 of HAC.
25 */
26 int fast_mp_montgomery_reduce (mp_int * x, mp_int * n, mp_digit rho)
27 {
28 int ix, res, olduse;
29 mp_word W[MP_WARRAY];
30
31 /* get old used count */
32 olduse = x->used;
33
34 /* grow a as required */
35 if (x->alloc < n->used + 1) {
36 if ((res = mp_grow (x, n->used + 1)) != MP_OKAY) {
37 return res;
38 }
39 }
40
41 /* first we have to get the digits of the input into
42 * an array of double precision words W[...]
43 */
44 {
45 register mp_word *_W;
46 register mp_digit *tmpx;
47
48 /* alias for the W[] array */
49 _W = W;
50
51 /* alias for the digits of x*/
52 tmpx = x->dp;
53
54 /* copy the digits of a into W[0..a->used-1] */
55 for (ix = 0; ix < x->used; ix++) {
56 *_W++ = *tmpx++;
57 }
58
59 /* zero the high words of W[a->used..m->used*2] */
60 for (; ix < n->used * 2 + 1; ix++) {
61 *_W++ = 0;
62 }
63 }
64
65 /* now we proceed to zero successive digits
66 * from the least significant upwards
67 */
68 for (ix = 0; ix < n->used; ix++) {
69 /* mu = ai * m' mod b
70 *
71 * We avoid a double precision multiplication (which isn't required)
72 * by casting the value down to a mp_digit. Note this requires
73 * that W[ix-1] have the carry cleared (see after the inner loop)
74 */
75 register mp_digit mu;
76 mu = (mp_digit) (((W[ix] & MP_MASK) * rho) & MP_MASK);
77
78 /* a = a + mu * m * b**i
79 *
80 * This is computed in place and on the fly. The multiplication
81 * by b**i is handled by offseting which columns the results
82 * are added to.
83 *
84 * Note the comba method normally doesn't handle carries in the
85 * inner loop In this case we fix the carry from the previous
86 * column since the Montgomery reduction requires digits of the
87 * result (so far) [see above] to work. This is
88 * handled by fixing up one carry after the inner loop. The
89 * carry fixups are done in order so after these loops the
90 * first m->used words of W[] have the carries fixed
91 */
92 {
93 register int iy;
94 register mp_digit *tmpn;
95 register mp_word *_W;
96
97 /* alias for the digits of the modulus */
98 tmpn = n->dp;
99
100 /* Alias for the columns set by an offset of ix */
101 _W = W + ix;
102
103 /* inner loop */
104 for (iy = 0; iy < n->used; iy++) {
105 *_W++ += ((mp_word)mu) * ((mp_word)*tmpn++);
106 }
107 }
108
109 /* now fix carry for next digit, W[ix+1] */
110 W[ix + 1] += W[ix] >> ((mp_word) DIGIT_BIT);
111 }
112
113 /* now we have to propagate the carries and
114 * shift the words downward [all those least
115 * significant digits we zeroed].
116 */
117 {
118 register mp_digit *tmpx;
119 register mp_word *_W, *_W1;
120
121 /* nox fix rest of carries */
122
123 /* alias for current word */
124 _W1 = W + ix;
125
126 /* alias for next word, where the carry goes */
127 _W = W + ++ix;
128
129 for (; ix <= n->used * 2 + 1; ix++) {
130 *_W++ += *_W1++ >> ((mp_word) DIGIT_BIT);
131 }
132
133 /* copy out, A = A/b**n
134 *
135 * The result is A/b**n but instead of converting from an
136 * array of mp_word to mp_digit than calling mp_rshd
137 * we just copy them in the right order
138 */
139
140 /* alias for destination word */
141 tmpx = x->dp;
142
143 /* alias for shifted double precision result */
144 _W = W + n->used;
145
146 for (ix = 0; ix < n->used + 1; ix++) {
147 *tmpx++ = (mp_digit)(*_W++ & ((mp_word) MP_MASK));
148 }
149
150 /* zero oldused digits, if the input a was larger than
151 * m->used+1 we'll have to clear the digits
152 */
153 for (; ix < olduse; ix++) {
154 *tmpx++ = 0;
155 }
156 }
157
158 /* set the max used and clamp */
159 x->used = n->used + 1;
160 mp_clamp (x);
161
162 /* if A >= m then A = A - m */
163 if (mp_cmp_mag (x, n) != MP_LT) {
164 return s_mp_sub (x, n, x);
165 }
166 return MP_OKAY;
167 }
168 #endif
169
170 /* $Source: /cvs/libtom/libtommath/bn_fast_mp_montgomery_reduce.c,v $ */
171 /* $Revision: 1.3 $ */
172 /* $Date: 2006/03/31 14:18:44 $ */