Mercurial > dropbear

comparison libtommath/bn_s_mp_exptmod.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_S_MP_EXPTMOD_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 #ifdef MP_LOW_MEM
18 #define TAB_SIZE 32
19 #else
20 #define TAB_SIZE 256
21 #endif
22
23 int s_mp_exptmod (mp_int * G, mp_int * X, mp_int * P, mp_int * Y, int redmode)
24 {
25 mp_int M[TAB_SIZE], res, mu;
26 mp_digit buf;
27 int err, bitbuf, bitcpy, bitcnt, mode, digidx, x, y, winsize;
28 int (*redux)(mp_int*,mp_int*,mp_int*);
29
30 /* find window size */
31 x = mp_count_bits (X);
32 if (x <= 7) {
33 winsize = 2;
34 } else if (x <= 36) {
35 winsize = 3;
36 } else if (x <= 140) {
37 winsize = 4;
38 } else if (x <= 450) {
39 winsize = 5;
40 } else if (x <= 1303) {
41 winsize = 6;
42 } else if (x <= 3529) {
43 winsize = 7;
44 } else {
45 winsize = 8;
46 }
47
48 #ifdef MP_LOW_MEM
49 if (winsize > 5) {
50 winsize = 5;
51 }
52 #endif
53
54 /* init M array */
55 /* init first cell */
56 if ((err = mp_init(&M[1])) != MP_OKAY) {
57 return err;
58 }
59
60 /* now init the second half of the array */
61 for (x = 1<<(winsize-1); x < (1 << winsize); x++) {
62 if ((err = mp_init(&M[x])) != MP_OKAY) {
63 for (y = 1<<(winsize-1); y < x; y++) {
64 mp_clear (&M[y]);
65 }
66 mp_clear(&M[1]);
67 return err;
68 }
69 }
70
71 /* create mu, used for Barrett reduction */
72 if ((err = mp_init (&mu)) != MP_OKAY) {
73 goto LBL_M;
74 }
75
76 if (redmode == 0) {
77 if ((err = mp_reduce_setup (&mu, P)) != MP_OKAY) {
78 goto LBL_MU;
79 }
80 redux = mp_reduce;
81 } else {
82 if ((err = mp_reduce_2k_setup_l (P, &mu)) != MP_OKAY) {
83 goto LBL_MU;
84 }
85 redux = mp_reduce_2k_l;
86 }
87
88 /* create M table
89 *
90 * The M table contains powers of the base,
91 * e.g. M[x] = G**x mod P
92 *
93 * The first half of the table is not
94 * computed though accept for M[0] and M[1]
95 */
96 if ((err = mp_mod (G, P, &M[1])) != MP_OKAY) {
97 goto LBL_MU;
98 }
99
100 /* compute the value at M[1<<(winsize-1)] by squaring
101 * M[1] (winsize-1) times
102 */
103 if ((err = mp_copy (&M[1], &M[1 << (winsize - 1)])) != MP_OKAY) {
104 goto LBL_MU;
105 }
106
107 for (x = 0; x < (winsize - 1); x++) {
108 /* square it */
109 if ((err = mp_sqr (&M[1 << (winsize - 1)],
110 &M[1 << (winsize - 1)])) != MP_OKAY) {
111 goto LBL_MU;
112 }
113
114 /* reduce modulo P */
115 if ((err = redux (&M[1 << (winsize - 1)], P, &mu)) != MP_OKAY) {
116 goto LBL_MU;
117 }
118 }
119
120 /* create upper table, that is M[x] = M[x-1] * M[1] (mod P)
121 * for x = (2**(winsize - 1) + 1) to (2**winsize - 1)
122 */
123 for (x = (1 << (winsize - 1)) + 1; x < (1 << winsize); x++) {
124 if ((err = mp_mul (&M[x - 1], &M[1], &M[x])) != MP_OKAY) {
125 goto LBL_MU;
126 }
127 if ((err = redux (&M[x], P, &mu)) != MP_OKAY) {
128 goto LBL_MU;
129 }
130 }
131
132 /* setup result */
133 if ((err = mp_init (&res)) != MP_OKAY) {
134 goto LBL_MU;
135 }
136 mp_set (&res, 1);
137
138 /* set initial mode and bit cnt */
139 mode = 0;
140 bitcnt = 1;
141 buf = 0;
142 digidx = X->used - 1;
143 bitcpy = 0;
144 bitbuf = 0;
145
146 for (;;) {
147 /* grab next digit as required */
148 if (--bitcnt == 0) {
149 /* if digidx == -1 we are out of digits */
150 if (digidx == -1) {
151 break;
152 }
153 /* read next digit and reset the bitcnt */
154 buf = X->dp[digidx--];
155 bitcnt = (int) DIGIT_BIT;
156 }
157
158 /* grab the next msb from the exponent */
159 y = (buf >> (mp_digit)(DIGIT_BIT - 1)) & 1;
160 buf <<= (mp_digit)1;
161
162 /* if the bit is zero and mode == 0 then we ignore it
163 * These represent the leading zero bits before the first 1 bit
164 * in the exponent. Technically this opt is not required but it
165 * does lower the # of trivial squaring/reductions used
166 */
167 if (mode == 0 && y == 0) {
168 continue;
169 }
170
171 /* if the bit is zero and mode == 1 then we square */
172 if (mode == 1 && y == 0) {
173 if ((err = mp_sqr (&res, &res)) != MP_OKAY) {
174 goto LBL_RES;
175 }
176 if ((err = redux (&res, P, &mu)) != MP_OKAY) {
177 goto LBL_RES;
178 }
179 continue;
180 }
181
182 /* else we add it to the window */
183 bitbuf |= (y << (winsize - ++bitcpy));
184 mode = 2;
185
186 if (bitcpy == winsize) {
187 /* ok window is filled so square as required and multiply */
188 /* square first */
189 for (x = 0; x < winsize; x++) {
190 if ((err = mp_sqr (&res, &res)) != MP_OKAY) {
191 goto LBL_RES;
192 }
193 if ((err = redux (&res, P, &mu)) != MP_OKAY) {
194 goto LBL_RES;
195 }
196 }
197
198 /* then multiply */
199 if ((err = mp_mul (&res, &M[bitbuf], &res)) != MP_OKAY) {
200 goto LBL_RES;
201 }
202 if ((err = redux (&res, P, &mu)) != MP_OKAY) {
203 goto LBL_RES;
204 }
205
206 /* empty window and reset */
207 bitcpy = 0;
208 bitbuf = 0;
209 mode = 1;
210 }
211 }
212
213 /* if bits remain then square/multiply */
214 if (mode == 2 && bitcpy > 0) {
215 /* square then multiply if the bit is set */
216 for (x = 0; x < bitcpy; x++) {
217 if ((err = mp_sqr (&res, &res)) != MP_OKAY) {
218 goto LBL_RES;
219 }
220 if ((err = redux (&res, P, &mu)) != MP_OKAY) {
221 goto LBL_RES;
222 }
223
224 bitbuf <<= 1;
225 if ((bitbuf & (1 << winsize)) != 0) {
226 /* then multiply */
227 if ((err = mp_mul (&res, &M[1], &res)) != MP_OKAY) {
228 goto LBL_RES;
229 }
230 if ((err = redux (&res, P, &mu)) != MP_OKAY) {
231 goto LBL_RES;
232 }
233 }
234 }
235 }
236
237 mp_exch (&res, Y);
238 err = MP_OKAY;
239 LBL_RES:mp_clear (&res);
240 LBL_MU:mp_clear (&mu);
241 LBL_M:
242 mp_clear(&M[1]);
243 for (x = 1<<(winsize-1); x < (1 << winsize); x++) {
244 mp_clear (&M[x]);
245 }
246 return err;
247 }
248 #endif
249
250 /* $Source: /cvs/libtom/libtommath/bn_s_mp_exptmod.c,v $ */
251 /* $Revision: 1.4 $ */
252 /* $Date: 2006/03/31 14:18:44 $ */