comparison libtommath/bn_mp_montgomery_reduce.c @ 399:a707e6148060

merge of '5fdf69ca60d1683cdd9f4c2595134bed26394834' and '6b61c50f4cf888bea302ac8fcf5dbb573b443251'
author Matt Johnston <matt@ucc.asn.au>
date Sat, 03 Feb 2007 08:20:34 +0000
parents 5ff8218bcee9
children 60fc6476e044
comparison
equal deleted inserted replaced
394:17d097fc111c 399:a707e6148060
1 #include <tommath.h>
2 #ifdef BN_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 int
20 mp_montgomery_reduce (mp_int * x, mp_int * n, mp_digit rho)
21 {
22 int ix, res, digs;
23 mp_digit mu;
24
25 /* can the fast reduction [comba] method be used?
26 *
27 * Note that unlike in mul you're safely allowed *less*
28 * than the available columns [255 per default] since carries
29 * are fixed up in the inner loop.
30 */
31 digs = n->used * 2 + 1;
32 if ((digs < MP_WARRAY) &&
33 n->used <
34 (1 << ((CHAR_BIT * sizeof (mp_word)) - (2 * DIGIT_BIT)))) {
35 return fast_mp_montgomery_reduce (x, n, rho);
36 }
37
38 /* grow the input as required */
39 if (x->alloc < digs) {
40 if ((res = mp_grow (x, digs)) != MP_OKAY) {
41 return res;
42 }
43 }
44 x->used = digs;
45
46 for (ix = 0; ix < n->used; ix++) {
47 /* mu = ai * rho mod b
48 *
49 * The value of rho must be precalculated via
50 * montgomery_setup() such that
51 * it equals -1/n0 mod b this allows the
52 * following inner loop to reduce the
53 * input one digit at a time
54 */
55 mu = (mp_digit) (((mp_word)x->dp[ix]) * ((mp_word)rho) & MP_MASK);
56
57 /* a = a + mu * m * b**i */
58 {
59 register int iy;
60 register mp_digit *tmpn, *tmpx, u;
61 register mp_word r;
62
63 /* alias for digits of the modulus */
64 tmpn = n->dp;
65
66 /* alias for the digits of x [the input] */
67 tmpx = x->dp + ix;
68
69 /* set the carry to zero */
70 u = 0;
71
72 /* Multiply and add in place */
73 for (iy = 0; iy < n->used; iy++) {
74 /* compute product and sum */
75 r = ((mp_word)mu) * ((mp_word)*tmpn++) +
76 ((mp_word) u) + ((mp_word) * tmpx);
77
78 /* get carry */
79 u = (mp_digit)(r >> ((mp_word) DIGIT_BIT));
80
81 /* fix digit */
82 *tmpx++ = (mp_digit)(r & ((mp_word) MP_MASK));
83 }
84 /* At this point the ix'th digit of x should be zero */
85
86
87 /* propagate carries upwards as required*/
88 while (u) {
89 *tmpx += u;
90 u = *tmpx >> DIGIT_BIT;
91 *tmpx++ &= MP_MASK;
92 }
93 }
94 }
95
96 /* at this point the n.used'th least
97 * significant digits of x are all zero
98 * which means we can shift x to the
99 * right by n.used digits and the
100 * residue is unchanged.
101 */
102
103 /* x = x/b**n.used */
104 mp_clamp(x);
105 mp_rshd (x, n->used);
106
107 /* if x >= n then x = x - n */
108 if (mp_cmp_mag (x, n) != MP_LT) {
109 return s_mp_sub (x, n, x);
110 }
111
112 return MP_OKAY;
113 }
114 #endif
115
116 /* $Source: /cvs/libtom/libtommath/bn_mp_montgomery_reduce.c,v $ */
117 /* $Revision: 1.3 $ */
118 /* $Date: 2006/03/31 14:18:44 $ */