← Back to branch summary

~ubuntu-branches/ubuntu/precise/kompozer/precise

« back to all changes in this revision

Viewing changes to mozilla/security/nss/lib/freebl/mpi/mp_gf2m.c

  • Committer: Bazaar Package Importer
  • Author(s): Anthony Yarusso
  • Date: 2007-08-27 01:11:03 UTC
  • Revision ID: james.westby@ubuntu.com-20070827011103-2jgf4s6532gqu2ka
Tags: upstream-0.7.10
Import upstream version 0.7.10

Show diffs side-by-side

added added

removed removed

Lines of Context:
mozilla/security/nss/lib/freebl/mpi/mp_gf2m.c
 
1
/*
 
2
 * Version: MPL 1.1/GPL 2.0/LGPL 2.1
 
3
 *
 
4
 * The contents of this file are subject to the Mozilla Public License Version
 
5
 * 1.1 (the "License"); you may not use this file except in compliance with
 
6
 * the License. You may obtain a copy of the License at
 
7
 * http://www.mozilla.org/MPL/
 
8
 *
 
9
 * Software distributed under the License is distributed on an "AS IS" basis,
 
10
 * WITHOUT WARRANTY OF ANY KIND, either express or implied. See the License
 
11
 * for the specific language governing rights and limitations under the
 
12
 * License.
 
13
 *
 
14
 * The Original Code is the Multi-precision Binary Polynomial Arithmetic 
 
15
 * Library.
 
16
 *
 
17
 * The Initial Developer of the Original Code is Sun Microsystems, Inc.
 
18
 * Portions created by Sun Microsystems, Inc. are Copyright (C) 2003
 
19
 * Sun Microsystems, Inc. All Rights Reserved.
 
20
 *
 
21
 * Contributor(s):
 
22
 *      Sheueling Chang Shantz <sheueling.chang@sun.com> and
 
23
 *      Douglas Stebila <douglas@stebila.ca> of Sun Laboratories.
 
24
 *
 
25
 * Alternatively, the contents of this file may be used under the terms of
 
26
 * either the GNU General Public License Version 2 or later (the "GPL"), or
 
27
 * the GNU Lesser General Public License Version 2.1 or later (the "LGPL"),
 
28
 * in which case the provisions of the GPL or the LGPL are applicable instead
 
29
 * of those above. If you wish to allow use of your version of this file only
 
30
 * under the terms of either the GPL or the LGPL, and not to allow others to
 
31
 * use your version of this file under the terms of the MPL, indicate your
 
32
 * decision by deleting the provisions above and replace them with the notice
 
33
 * and other provisions required by the GPL or the LGPL. If you do not delete
 
34
 * the provisions above, a recipient may use your version of this file under
 
35
 * the terms of any one of the MPL, the GPL or the LGPL.
 
36
 *
 
37
 */
 
38
 
 
39
#include "mp_gf2m.h"
 
40
#include "mp_gf2m-priv.h"
 
41
#include "mplogic.h"
 
42
#include "mpi-priv.h"
 
43
 
 
44
const mp_digit mp_gf2m_sqr_tb[16] =
 
45
{
 
46
      0,     1,     4,     5,    16,    17,    20,    21,
 
47
     64,    65,    68,    69,    80,    81,    84,    85
 
48
};
 
49
 
 
50
/* Multiply two binary polynomials mp_digits a, b.
 
51
 * Result is a polynomial with degree < 2 * MP_DIGIT_BITS - 1.
 
52
 * Output in two mp_digits rh, rl.
 
53
 */
 
54
#if MP_DIGIT_BITS == 32
 
55
void 
 
56
s_bmul_1x1(mp_digit *rh, mp_digit *rl, const mp_digit a, const mp_digit b)
 
57
{
 
58
    register mp_digit h, l, s;
 
59
    mp_digit tab[8], top2b = a >> 30; 
 
60
    register mp_digit a1, a2, a4;
 
61
 
 
62
    a1 = a & (0x3FFFFFFF); a2 = a1 << 1; a4 = a2 << 1;
 
63
 
 
64
    tab[0] =  0; tab[1] = a1;    tab[2] = a2;    tab[3] = a1^a2;
 
65
    tab[4] = a4; tab[5] = a1^a4; tab[6] = a2^a4; tab[7] = a1^a2^a4;
 
66
 
 
67
    s = tab[b       & 0x7]; l  = s;
 
68
    s = tab[b >>  3 & 0x7]; l ^= s <<  3; h  = s >> 29;
 
69
    s = tab[b >>  6 & 0x7]; l ^= s <<  6; h ^= s >> 26;
 
70
    s = tab[b >>  9 & 0x7]; l ^= s <<  9; h ^= s >> 23;
 
71
    s = tab[b >> 12 & 0x7]; l ^= s << 12; h ^= s >> 20;
 
72
    s = tab[b >> 15 & 0x7]; l ^= s << 15; h ^= s >> 17;
 
73
    s = tab[b >> 18 & 0x7]; l ^= s << 18; h ^= s >> 14;
 
74
    s = tab[b >> 21 & 0x7]; l ^= s << 21; h ^= s >> 11;
 
75
    s = tab[b >> 24 & 0x7]; l ^= s << 24; h ^= s >>  8;
 
76
    s = tab[b >> 27 & 0x7]; l ^= s << 27; h ^= s >>  5;
 
77
    s = tab[b >> 30      ]; l ^= s << 30; h ^= s >>  2;
 
78
 
 
79
    /* compensate for the top two bits of a */
 
80
 
 
81
    if (top2b & 01) { l ^= b << 30; h ^= b >> 2; } 
 
82
    if (top2b & 02) { l ^= b << 31; h ^= b >> 1; } 
 
83
 
 
84
    *rh = h; *rl = l;
 
85
} 
 
86
#else
 
87
void 
 
88
s_bmul_1x1(mp_digit *rh, mp_digit *rl, const mp_digit a, const mp_digit b)
 
89
{
 
90
    register mp_digit h, l, s;
 
91
    mp_digit tab[16], top3b = a >> 61;
 
92
    register mp_digit a1, a2, a4, a8;
 
93
 
 
94
    a1 = a & (0x1FFFFFFFFFFFFFFF); a2 = a1 << 1; 
 
95
    a4 = a2 << 1; a8 = a4 << 1;
 
96
    tab[ 0] = 0;     tab[ 1] = a1;       tab[ 2] = a2;       tab[ 3] = a1^a2;
 
97
    tab[ 4] = a4;    tab[ 5] = a1^a4;    tab[ 6] = a2^a4;    tab[ 7] = a1^a2^a4;
 
98
    tab[ 8] = a8;    tab[ 9] = a1^a8;    tab[10] = a2^a8;    tab[11] = a1^a2^a8;
 
99
    tab[12] = a4^a8; tab[13] = a1^a4^a8; tab[14] = a2^a4^a8; tab[15] = a1^a2^a4^a8;
 
100
 
 
101
    s = tab[b       & 0xF]; l  = s;
 
102
    s = tab[b >>  4 & 0xF]; l ^= s <<  4; h  = s >> 60;
 
103
    s = tab[b >>  8 & 0xF]; l ^= s <<  8; h ^= s >> 56;
 
104
    s = tab[b >> 12 & 0xF]; l ^= s << 12; h ^= s >> 52;
 
105
    s = tab[b >> 16 & 0xF]; l ^= s << 16; h ^= s >> 48;
 
106
    s = tab[b >> 20 & 0xF]; l ^= s << 20; h ^= s >> 44;
 
107
    s = tab[b >> 24 & 0xF]; l ^= s << 24; h ^= s >> 40;
 
108
    s = tab[b >> 28 & 0xF]; l ^= s << 28; h ^= s >> 36;
 
109
    s = tab[b >> 32 & 0xF]; l ^= s << 32; h ^= s >> 32;
 
110
    s = tab[b >> 36 & 0xF]; l ^= s << 36; h ^= s >> 28;
 
111
    s = tab[b >> 40 & 0xF]; l ^= s << 40; h ^= s >> 24;
 
112
    s = tab[b >> 44 & 0xF]; l ^= s << 44; h ^= s >> 20;
 
113
    s = tab[b >> 48 & 0xF]; l ^= s << 48; h ^= s >> 16;
 
114
    s = tab[b >> 52 & 0xF]; l ^= s << 52; h ^= s >> 12;
 
115
    s = tab[b >> 56 & 0xF]; l ^= s << 56; h ^= s >>  8;
 
116
    s = tab[b >> 60      ]; l ^= s << 60; h ^= s >>  4;
 
117
 
 
118
    /* compensate for the top three bits of a */
 
119
 
 
120
    if (top3b & 01) { l ^= b << 61; h ^= b >> 3; } 
 
121
    if (top3b & 02) { l ^= b << 62; h ^= b >> 2; } 
 
122
    if (top3b & 04) { l ^= b << 63; h ^= b >> 1; } 
 
123
 
 
124
    *rh = h; *rl = l;
 
125
} 
 
126
#endif
 
127
 
 
128
/* Compute xor-multiply of two binary polynomials  (a1, a0) x (b1, b0)  
 
129
 * result is a binary polynomial in 4 mp_digits r[4].
 
130
 * The caller MUST ensure that r has the right amount of space allocated.
 
131
 */
 
132
void 
 
133
s_bmul_2x2(mp_digit *r, const mp_digit a1, const mp_digit a0, const mp_digit b1,
 
134
           const mp_digit b0)
 
135
{
 
136
    mp_digit m1, m0;
 
137
    /* r[3] = h1, r[2] = h0; r[1] = l1; r[0] = l0 */
 
138
    s_bmul_1x1(r+3, r+2, a1, b1);
 
139
    s_bmul_1x1(r+1, r, a0, b0);
 
140
    s_bmul_1x1(&m1, &m0, a0 ^ a1, b0 ^ b1);
 
141
    /* Correction on m1 ^= l1 ^ h1; m0 ^= l0 ^ h0; */
 
142
    r[2] ^= m1 ^ r[1] ^ r[3];  /* h0 ^= m1 ^ l1 ^ h1; */
 
143
    r[1]  = r[3] ^ r[2] ^ r[0] ^ m1 ^ m0;  /* l1 ^= l0 ^ h0 ^ m0; */
 
144
}
 
145
 
 
146
/* Compute xor-multiply of two binary polynomials  (a2, a1, a0) x (b2, b1, b0)  
 
147
 * result is a binary polynomial in 6 mp_digits r[6].
 
148
 * The caller MUST ensure that r has the right amount of space allocated.
 
149
 */
 
150
void 
 
151
s_bmul_3x3(mp_digit *r, const mp_digit a2, const mp_digit a1, const mp_digit a0, 
 
152
        const mp_digit b2, const mp_digit b1, const mp_digit b0)
 
153
{
 
154
        mp_digit zm[4];
 
155
 
 
156
        s_bmul_1x1(r+5, r+4, a2, b2);         /* fill top 2 words */
 
157
        s_bmul_2x2(zm, a1, a2^a0, b1, b2^b0); /* fill middle 4 words */
 
158
        s_bmul_2x2(r, a1, a0, b1, b0);        /* fill bottom 4 words */
 
159
 
 
160
        zm[3] ^= r[3];
 
161
        zm[2] ^= r[2]; 
 
162
        zm[1] ^= r[1] ^ r[5];
 
163
        zm[0] ^= r[0] ^ r[4];
 
164
 
 
165
        r[5]  ^= zm[3];
 
166
        r[4]  ^= zm[2];
 
167
        r[3]  ^= zm[1];
 
168
        r[2]  ^= zm[0];
 
169
}
 
170
 
 
171
/* Compute xor-multiply of two binary polynomials  (a3, a2, a1, a0) x (b3, b2, b1, b0)  
 
172
 * result is a binary polynomial in 8 mp_digits r[8].
 
173
 * The caller MUST ensure that r has the right amount of space allocated.
 
174
 */
 
175
void s_bmul_4x4(mp_digit *r, const mp_digit a3, const mp_digit a2, const mp_digit a1, 
 
176
        const mp_digit a0, const mp_digit b3, const mp_digit b2, const mp_digit b1, 
 
177
        const mp_digit b0)
 
178
{
 
179
        mp_digit zm[4];
 
180
 
 
181
        s_bmul_2x2(r+4, a3, a2, b3, b2);            /* fill top 4 words */
 
182
        s_bmul_2x2(zm, a3^a1, a2^a0, b3^b1, b2^b0); /* fill middle 4 words */
 
183
        s_bmul_2x2(r, a1, a0, b1, b0);              /* fill bottom 4 words */
 
184
 
 
185
        zm[3] ^= r[3] ^ r[7]; 
 
186
        zm[2] ^= r[2] ^ r[6]; 
 
187
        zm[1] ^= r[1] ^ r[5]; 
 
188
        zm[0] ^= r[0] ^ r[4]; 
 
189
 
 
190
        r[5]  ^= zm[3];    
 
191
        r[4]  ^= zm[2];
 
192
        r[3]  ^= zm[1];    
 
193
        r[2]  ^= zm[0];
 
194
}
 
195
 
 
196
/* Compute addition of two binary polynomials a and b,
 
197
 * store result in c; c could be a or b, a and b could be equal; 
 
198
 * c is the bitwise XOR of a and b.
 
199
 */
 
200
mp_err
 
201
mp_badd(const mp_int *a, const mp_int *b, mp_int *c)
 
202
{
 
203
    mp_digit *pa, *pb, *pc;
 
204
    mp_size ix;
 
205
    mp_size used_pa, used_pb;
 
206
    mp_err res = MP_OKAY;
 
207
 
 
208
    /* Add all digits up to the precision of b.  If b had more
 
209
     * precision than a initially, swap a, b first
 
210
     */
 
211
    if (MP_USED(a) >= MP_USED(b)) {
 
212
        pa = MP_DIGITS(a);
 
213
        pb = MP_DIGITS(b);
 
214
        used_pa = MP_USED(a);
 
215
        used_pb = MP_USED(b);
 
216
    } else {
 
217
        pa = MP_DIGITS(b);
 
218
        pb = MP_DIGITS(a);
 
219
        used_pa = MP_USED(b);
 
220
        used_pb = MP_USED(a);
 
221
    }
 
222
 
 
223
    /* Make sure c has enough precision for the output value */
 
224
    MP_CHECKOK( s_mp_pad(c, used_pa) );
 
225
 
 
226
    /* Do word-by-word xor */
 
227
    pc = MP_DIGITS(c);
 
228
    for (ix = 0; ix < used_pb; ix++) {
 
229
        (*pc++) = (*pa++) ^ (*pb++);
 
230
    }
 
231
 
 
232
    /* Finish the rest of digits until we're actually done */
 
233
    for (; ix < used_pa; ++ix) {
 
234
        *pc++ = *pa++;
 
235
    }
 
236
 
 
237
    MP_USED(c) = used_pa;
 
238
    MP_SIGN(c) = ZPOS;
 
239
    s_mp_clamp(c);
 
240
 
 
241
CLEANUP:
 
242
    return res;
 
243
} 
 
244
 
 
245
#define s_mp_div2(a) MP_CHECKOK( mpl_rsh((a), (a), 1) );
 
246
 
 
247
/* Compute binary polynomial multiply d = a * b */
 
248
static void 
 
249
s_bmul_d(const mp_digit *a, mp_size a_len, mp_digit b, mp_digit *d)
 
250
{
 
251
    mp_digit a_i, a0b0, a1b1, carry = 0;
 
252
    while (a_len--) {
 
253
        a_i = *a++;
 
254
        s_bmul_1x1(&a1b1, &a0b0, a_i, b);
 
255
        *d++ = a0b0 ^ carry;
 
256
        carry = a1b1;
 
257
    }
 
258
    *d = carry;
 
259
}
 
260
 
 
261
/* Compute binary polynomial xor multiply accumulate d ^= a * b */
 
262
static void 
 
263
s_bmul_d_add(const mp_digit *a, mp_size a_len, mp_digit b, mp_digit *d)
 
264
{
 
265
    mp_digit a_i, a0b0, a1b1, carry = 0;
 
266
    while (a_len--) {
 
267
        a_i = *a++;
 
268
        s_bmul_1x1(&a1b1, &a0b0, a_i, b);
 
269
        *d++ ^= a0b0 ^ carry;
 
270
        carry = a1b1;
 
271
    }
 
272
    *d ^= carry;
 
273
}
 
274
 
 
275
/* Compute binary polynomial xor multiply c = a * b.  
 
276
 * All parameters may be identical.
 
277
 */
 
278
mp_err 
 
279
mp_bmul(const mp_int *a, const mp_int *b, mp_int *c)
 
280
{
 
281
    mp_digit *pb, b_i;
 
282
    mp_int tmp;
 
283
    mp_size ib, a_used, b_used;
 
284
    mp_err res = MP_OKAY;
 
285
 
 
286
    MP_DIGITS(&tmp) = 0;
 
287
 
 
288
    ARGCHK(a != NULL && b != NULL && c != NULL, MP_BADARG);
 
289
 
 
290
    if (a == c) {
 
291
        MP_CHECKOK( mp_init_copy(&tmp, a) );
 
292
        if (a == b)
 
293
            b = &tmp;
 
294
        a = &tmp;
 
295
    } else if (b == c) {
 
296
        MP_CHECKOK( mp_init_copy(&tmp, b) );
 
297
        b = &tmp;
 
298
    }
 
299
 
 
300
    if (MP_USED(a) < MP_USED(b)) {
 
301
        const mp_int *xch = b;      /* switch a and b if b longer */
 
302
        b = a;
 
303
        a = xch;
 
304
    }
 
305
 
 
306
    MP_USED(c) = 1; MP_DIGIT(c, 0) = 0;
 
307
    MP_CHECKOK( s_mp_pad(c, USED(a) + USED(b)) );
 
308
 
 
309
    pb = MP_DIGITS(b);
 
310
    s_bmul_d(MP_DIGITS(a), MP_USED(a), *pb++, MP_DIGITS(c));
 
311
 
 
312
    /* Outer loop:  Digits of b */
 
313
    a_used = MP_USED(a);
 
314
    b_used = MP_USED(b);
 
315
        MP_USED(c) = a_used + b_used;
 
316
    for (ib = 1; ib < b_used; ib++) {
 
317
        b_i = *pb++;
 
318
 
 
319
        /* Inner product:  Digits of a */
 
320
        if (b_i)
 
321
            s_bmul_d_add(MP_DIGITS(a), a_used, b_i, MP_DIGITS(c) + ib);
 
322
        else
 
323
            MP_DIGIT(c, ib + a_used) = b_i;
 
324
    }
 
325
 
 
326
    s_mp_clamp(c);
 
327
 
 
328
    SIGN(c) = ZPOS;
 
329
 
 
330
CLEANUP:
 
331
    mp_clear(&tmp);
 
332
    return res;
 
333
}
 
334
 
 
335
 
 
336
/* Compute modular reduction of a and store result in r.  
 
337
 * r could be a. 
 
338
 * For modular arithmetic, the irreducible polynomial f(t) is represented 
 
339
 * as an array of int[], where f(t) is of the form: 
 
340
 *     f(t) = t^p[0] + t^p[1] + ... + t^p[k]
 
341
 * where m = p[0] > p[1] > ... > p[k] = 0.
 
342
 */
 
343
mp_err
 
344
mp_bmod(const mp_int *a, const unsigned int p[], mp_int *r)
 
345
{
 
346
    int j, k;
 
347
    int n, dN, d0, d1;
 
348
    mp_digit zz, *z, tmp;
 
349
    mp_size used;
 
350
    mp_err res = MP_OKAY;
 
351
 
 
352
    /* The algorithm does the reduction in place in r, 
 
353
     * if a != r, copy a into r first so reduction can be done in r
 
354
     */
 
355
    if (a != r) {
 
356
        MP_CHECKOK( mp_copy(a, r) );
 
357
    }
 
358
    z = MP_DIGITS(r);
 
359
 
 
360
    /* start reduction */
 
361
    dN = p[0] / MP_DIGIT_BITS;
 
362
    used = MP_USED(r);
 
363
 
 
364
    for (j = used - 1; j > dN;) {
 
365
 
 
366
        zz = z[j];
 
367
        if (zz == 0) {
 
368
            j--; continue;
 
369
        }
 
370
        z[j] = 0;
 
371
 
 
372
        for (k = 1; p[k] > 0; k++) {
 
373
            /* reducing component t^p[k] */
 
374
            n = p[0] - p[k];
 
375
            d0 = n % MP_DIGIT_BITS;  
 
376
            d1 = MP_DIGIT_BITS - d0;
 
377
            n /= MP_DIGIT_BITS;
 
378
            z[j-n] ^= (zz>>d0);
 
379
            if (d0) 
 
380
                z[j-n-1] ^= (zz<<d1);
 
381
        }
 
382
 
 
383
        /* reducing component t^0 */
 
384
        n = dN;  
 
385
        d0 = p[0] % MP_DIGIT_BITS;
 
386
        d1 = MP_DIGIT_BITS - d0;
 
387
        z[j-n] ^= (zz >> d0);
 
388
        if (d0) 
 
389
            z[j-n-1] ^= (zz << d1);
 
390
 
 
391
    }
 
392
 
 
393
    /* final round of reduction */
 
394
    while (j == dN) {
 
395
 
 
396
        d0 = p[0] % MP_DIGIT_BITS;
 
397
        zz = z[dN] >> d0;  
 
398
        if (zz == 0) break;
 
399
        d1 = MP_DIGIT_BITS - d0;
 
400
 
 
401
        /* clear up the top d1 bits */
 
402
        if (d0) z[dN] = (z[dN] << d1) >> d1; 
 
403
        *z ^= zz; /* reduction t^0 component */
 
404
 
 
405
        for (k = 1; p[k] > 0; k++) {
 
406
            /* reducing component t^p[k]*/
 
407
            n = p[k] / MP_DIGIT_BITS;
 
408
            d0 = p[k] % MP_DIGIT_BITS;
 
409
            d1 = MP_DIGIT_BITS - d0;
 
410
            z[n] ^= (zz << d0);
 
411
            tmp = zz >> d1;
 
412
            if (d0 && tmp)
 
413
                z[n+1] ^= tmp;
 
414
        }
 
415
    }
 
416
 
 
417
    s_mp_clamp(r);
 
418
CLEANUP:
 
419
    return res;
 
420
}
 
421
 
 
422
/* Compute the product of two polynomials a and b, reduce modulo p, 
 
423
 * Store the result in r.  r could be a or b; a could be b.
 
424
 */
 
425
mp_err 
 
426
mp_bmulmod(const mp_int *a, const mp_int *b, const unsigned int p[], mp_int *r)
 
427
{
 
428
    mp_err res;
 
429
    
 
430
    if (a == b) return mp_bsqrmod(a, p, r);
 
431
    if ((res = mp_bmul(a, b, r) ) != MP_OKAY)
 
432
        return res;
 
433
    return mp_bmod(r, p, r);
 
434
}
 
435
 
 
436
/* Compute binary polynomial squaring c = a*a mod p .  
 
437
 * Parameter r and a can be identical.
 
438
 */
 
439
 
 
440
mp_err 
 
441
mp_bsqrmod(const mp_int *a, const unsigned int p[], mp_int *r)
 
442
{
 
443
    mp_digit *pa, *pr, a_i;
 
444
    mp_int tmp;
 
445
    mp_size ia, a_used;
 
446
    mp_err res;
 
447
 
 
448
    ARGCHK(a != NULL && r != NULL, MP_BADARG);
 
449
    MP_DIGITS(&tmp) = 0;
 
450
 
 
451
    if (a == r) {
 
452
        MP_CHECKOK( mp_init_copy(&tmp, a) );
 
453
        a = &tmp;
 
454
    }
 
455
 
 
456
    MP_USED(r) = 1; MP_DIGIT(r, 0) = 0;
 
457
    MP_CHECKOK( s_mp_pad(r, 2*USED(a)) );
 
458
 
 
459
    pa = MP_DIGITS(a);
 
460
    pr = MP_DIGITS(r);
 
461
    a_used = MP_USED(a);
 
462
        MP_USED(r) = 2 * a_used;
 
463
 
 
464
    for (ia = 0; ia < a_used; ia++) {
 
465
        a_i = *pa++;
 
466
        *pr++ = gf2m_SQR0(a_i);
 
467
        *pr++ = gf2m_SQR1(a_i);
 
468
    }
 
469
 
 
470
    MP_CHECKOK( mp_bmod(r, p, r) );
 
471
    s_mp_clamp(r);
 
472
    SIGN(r) = ZPOS;
 
473
 
 
474
CLEANUP:
 
475
    mp_clear(&tmp);
 
476
    return res;
 
477
}
 
478
 
 
479
/* Compute binary polynomial y/x mod p, y divided by x, reduce modulo p.
 
480
 * Store the result in r. r could be x or y, and x could equal y.
 
481
 * Uses algorithm Modular_Division_GF(2^m) from 
 
482
 *     Chang-Shantz, S.  "From Euclid's GCD to Montgomery Multiplication to 
 
483
 *     the Great Divide".
 
484
 */
 
485
int 
 
486
mp_bdivmod(const mp_int *y, const mp_int *x, const mp_int *pp, 
 
487
    const unsigned int p[], mp_int *r)
 
488
{
 
489
    mp_int aa, bb, uu;
 
490
    mp_int *a, *b, *u, *v;
 
491
    mp_err res = MP_OKAY;
 
492
 
 
493
    MP_DIGITS(&aa) = 0;
 
494
    MP_DIGITS(&bb) = 0;
 
495
    MP_DIGITS(&uu) = 0;
 
496
 
 
497
    MP_CHECKOK( mp_init_copy(&aa, x) );
 
498
    MP_CHECKOK( mp_init_copy(&uu, y) );
 
499
    MP_CHECKOK( mp_init_copy(&bb, pp) );
 
500
    MP_CHECKOK( s_mp_pad(r, USED(pp)) );
 
501
    MP_USED(r) = 1; MP_DIGIT(r, 0) = 0;
 
502
 
 
503
    a = &aa; b= &bb; u=&uu; v=r;
 
504
    /* reduce x and y mod p */
 
505
    MP_CHECKOK( mp_bmod(a, p, a) );
 
506
    MP_CHECKOK( mp_bmod(u, p, u) );
 
507
 
 
508
    while (!mp_isodd(a)) {
 
509
        s_mp_div2(a);
 
510
        if (mp_isodd(u)) {
 
511
            MP_CHECKOK( mp_badd(u, pp, u) );
 
512
        }
 
513
        s_mp_div2(u);
 
514
    }
 
515
 
 
516
    do {
 
517
        if (mp_cmp_mag(b, a) > 0) {
 
518
            MP_CHECKOK( mp_badd(b, a, b) );
 
519
            MP_CHECKOK( mp_badd(v, u, v) );
 
520
            do {
 
521
                s_mp_div2(b);
 
522
                if (mp_isodd(v)) {
 
523
                    MP_CHECKOK( mp_badd(v, pp, v) );
 
524
                }
 
525
                s_mp_div2(v);
 
526
            } while (!mp_isodd(b));
 
527
        }
 
528
        else if ((MP_DIGIT(a,0) == 1) && (MP_USED(a) == 1))
 
529
            break;
 
530
        else {
 
531
            MP_CHECKOK( mp_badd(a, b, a) );
 
532
            MP_CHECKOK( mp_badd(u, v, u) );
 
533
            do {
 
534
                s_mp_div2(a);
 
535
                if (mp_isodd(u)) {
 
536
                    MP_CHECKOK( mp_badd(u, pp, u) );
 
537
                }
 
538
                s_mp_div2(u);
 
539
            } while (!mp_isodd(a));
 
540
        }
 
541
    } while (1);
 
542
 
 
543
    MP_CHECKOK( mp_copy(u, r) );
 
544
 
 
545
CLEANUP:
 
546
    return res;
 
547
 
 
548
}
 
549
 
 
550
/* Convert the bit-string representation of a polynomial a into an array
 
551
 * of integers corresponding to the bits with non-zero coefficient.
 
552
 * Up to max elements of the array will be filled.  Return value is total
 
553
 * number of coefficients that would be extracted if array was large enough.
 
554
 */
 
555
int
 
556
mp_bpoly2arr(const mp_int *a, unsigned int p[], int max)
 
557
{
 
558
    int i, j, k;
 
559
    mp_digit top_bit, mask;
 
560
 
 
561
    top_bit = 1;
 
562
    top_bit <<= MP_DIGIT_BIT - 1;
 
563
 
 
564
    for (k = 0; k < max; k++) p[k] = 0;
 
565
    k = 0;
 
566
 
 
567
    for (i = MP_USED(a) - 1; i >= 0; i--) {
 
568
        mask = top_bit;
 
569
        for (j = MP_DIGIT_BIT - 1; j >= 0; j--) {
 
570
            if (MP_DIGITS(a)[i] & mask) {
 
571
                if (k < max) p[k] = MP_DIGIT_BIT * i + j;
 
572
                k++;
 
573
            }
 
574
            mask >>= 1;
 
575
        }
 
576
    }
 
577
 
 
578
    return k;
 
579
}
 
580
 
 
581
/* Convert the coefficient array representation of a polynomial to a 
 
582
 * bit-string.  The array must be terminated by 0.
 
583
 */
 
584
mp_err
 
585
mp_barr2poly(const unsigned int p[], mp_int *a)
 
586
{
 
587
 
 
588
    mp_err res = MP_OKAY;
 
589
    int i;
 
590
 
 
591
    mp_zero(a);
 
592
    for (i = 0; p[i] > 0; i++) {
 
593
        MP_CHECKOK( mpl_set_bit(a, p[i], 1) );
 
594
    }
 
595
    MP_CHECKOK( mpl_set_bit(a, 0, 1) );
 
596
        
 
597
CLEANUP:
 
598
    return res;
 
599
}