Mercurial > dropbear

comparison bn_mp_toom_sqr.c @ 282:91fbc376f010 libtommath-orig libtommath-0.35

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Import of libtommath 0.35 From ltm-0.35.tar.bz2 SHA1 of 3f193dbae9351e92d02530994fa18236f7fde01c
author Matt Johnston <matt@ucc.asn.au>
date Wed, 08 Mar 2006 13:16:18 +0000
parents
children 97db060d0ef5
comparison
equal deleted inserted replaced
-1:000000000000 282:91fbc376f010
1 #include <tommath.h>
2 #ifdef BN_MP_TOOM_SQR_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.org
16 */
17
18 /* squaring using Toom-Cook 3-way algorithm */
19 int
20 mp_toom_sqr(mp_int *a, mp_int *b)
21 {
22 mp_int w0, w1, w2, w3, w4, tmp1, a0, a1, a2;
23 int res, B;
24
25 /* init temps */
26 if ((res = mp_init_multi(&w0, &w1, &w2, &w3, &w4, &a0, &a1, &a2, &tmp1, NULL)) != MP_OKAY) {
27 return res;
28 }
29
30 /* B */
31 B = a->used / 3;
32
33 /* a = a2 * B**2 + a1 * B + a0 */
34 if ((res = mp_mod_2d(a, DIGIT_BIT * B, &a0)) != MP_OKAY) {
35 goto ERR;
36 }
37
38 if ((res = mp_copy(a, &a1)) != MP_OKAY) {
39 goto ERR;
40 }
41 mp_rshd(&a1, B);
42 mp_mod_2d(&a1, DIGIT_BIT * B, &a1);
43
44 if ((res = mp_copy(a, &a2)) != MP_OKAY) {
45 goto ERR;
46 }
47 mp_rshd(&a2, B*2);
48
49 /* w0 = a0*a0 */
50 if ((res = mp_sqr(&a0, &w0)) != MP_OKAY) {
51 goto ERR;
52 }
53
54 /* w4 = a2 * a2 */
55 if ((res = mp_sqr(&a2, &w4)) != MP_OKAY) {
56 goto ERR;
57 }
58
59 /* w1 = (a2 + 2(a1 + 2a0))**2 */
60 if ((res = mp_mul_2(&a0, &tmp1)) != MP_OKAY) {
61 goto ERR;
62 }
63 if ((res = mp_add(&tmp1, &a1, &tmp1)) != MP_OKAY) {
64 goto ERR;
65 }
66 if ((res = mp_mul_2(&tmp1, &tmp1)) != MP_OKAY) {
67 goto ERR;
68 }
69 if ((res = mp_add(&tmp1, &a2, &tmp1)) != MP_OKAY) {
70 goto ERR;
71 }
72
73 if ((res = mp_sqr(&tmp1, &w1)) != MP_OKAY) {
74 goto ERR;
75 }
76
77 /* w3 = (a0 + 2(a1 + 2a2))**2 */
78 if ((res = mp_mul_2(&a2, &tmp1)) != MP_OKAY) {
79 goto ERR;
80 }
81 if ((res = mp_add(&tmp1, &a1, &tmp1)) != MP_OKAY) {
82 goto ERR;
83 }
84 if ((res = mp_mul_2(&tmp1, &tmp1)) != MP_OKAY) {
85 goto ERR;
86 }
87 if ((res = mp_add(&tmp1, &a0, &tmp1)) != MP_OKAY) {
88 goto ERR;
89 }
90
91 if ((res = mp_sqr(&tmp1, &w3)) != MP_OKAY) {
92 goto ERR;
93 }
94
95
96 /* w2 = (a2 + a1 + a0)**2 */
97 if ((res = mp_add(&a2, &a1, &tmp1)) != MP_OKAY) {
98 goto ERR;
99 }
100 if ((res = mp_add(&tmp1, &a0, &tmp1)) != MP_OKAY) {
101 goto ERR;
102 }
103 if ((res = mp_sqr(&tmp1, &w2)) != MP_OKAY) {
104 goto ERR;
105 }
106
107 /* now solve the matrix
108
109 0 0 0 0 1
110 1 2 4 8 16
111 1 1 1 1 1
112 16 8 4 2 1
113 1 0 0 0 0
114
115 using 12 subtractions, 4 shifts, 2 small divisions and 1 small multiplication.
116 */
117
118 /* r1 - r4 */
119 if ((res = mp_sub(&w1, &w4, &w1)) != MP_OKAY) {
120 goto ERR;
121 }
122 /* r3 - r0 */
123 if ((res = mp_sub(&w3, &w0, &w3)) != MP_OKAY) {
124 goto ERR;
125 }
126 /* r1/2 */
127 if ((res = mp_div_2(&w1, &w1)) != MP_OKAY) {
128 goto ERR;
129 }
130 /* r3/2 */
131 if ((res = mp_div_2(&w3, &w3)) != MP_OKAY) {
132 goto ERR;
133 }
134 /* r2 - r0 - r4 */
135 if ((res = mp_sub(&w2, &w0, &w2)) != MP_OKAY) {
136 goto ERR;
137 }
138 if ((res = mp_sub(&w2, &w4, &w2)) != MP_OKAY) {
139 goto ERR;
140 }
141 /* r1 - r2 */
142 if ((res = mp_sub(&w1, &w2, &w1)) != MP_OKAY) {
143 goto ERR;
144 }
145 /* r3 - r2 */
146 if ((res = mp_sub(&w3, &w2, &w3)) != MP_OKAY) {
147 goto ERR;
148 }
149 /* r1 - 8r0 */
150 if ((res = mp_mul_2d(&w0, 3, &tmp1)) != MP_OKAY) {
151 goto ERR;
152 }
153 if ((res = mp_sub(&w1, &tmp1, &w1)) != MP_OKAY) {
154 goto ERR;
155 }
156 /* r3 - 8r4 */
157 if ((res = mp_mul_2d(&w4, 3, &tmp1)) != MP_OKAY) {
158 goto ERR;
159 }
160 if ((res = mp_sub(&w3, &tmp1, &w3)) != MP_OKAY) {
161 goto ERR;
162 }
163 /* 3r2 - r1 - r3 */
164 if ((res = mp_mul_d(&w2, 3, &w2)) != MP_OKAY) {
165 goto ERR;
166 }
167 if ((res = mp_sub(&w2, &w1, &w2)) != MP_OKAY) {
168 goto ERR;
169 }
170 if ((res = mp_sub(&w2, &w3, &w2)) != MP_OKAY) {
171 goto ERR;
172 }
173 /* r1 - r2 */
174 if ((res = mp_sub(&w1, &w2, &w1)) != MP_OKAY) {
175 goto ERR;
176 }
177 /* r3 - r2 */
178 if ((res = mp_sub(&w3, &w2, &w3)) != MP_OKAY) {
179 goto ERR;
180 }
181 /* r1/3 */
182 if ((res = mp_div_3(&w1, &w1, NULL)) != MP_OKAY) {
183 goto ERR;
184 }
185 /* r3/3 */
186 if ((res = mp_div_3(&w3, &w3, NULL)) != MP_OKAY) {
187 goto ERR;
188 }
189
190 /* at this point shift W[n] by B*n */
191 if ((res = mp_lshd(&w1, 1*B)) != MP_OKAY) {
192 goto ERR;
193 }
194 if ((res = mp_lshd(&w2, 2*B)) != MP_OKAY) {
195 goto ERR;
196 }
197 if ((res = mp_lshd(&w3, 3*B)) != MP_OKAY) {
198 goto ERR;
199 }
200 if ((res = mp_lshd(&w4, 4*B)) != MP_OKAY) {
201 goto ERR;
202 }
203
204 if ((res = mp_add(&w0, &w1, b)) != MP_OKAY) {
205 goto ERR;
206 }
207 if ((res = mp_add(&w2, &w3, &tmp1)) != MP_OKAY) {
208 goto ERR;
209 }
210 if ((res = mp_add(&w4, &tmp1, &tmp1)) != MP_OKAY) {
211 goto ERR;
212 }
213 if ((res = mp_add(&tmp1, b, b)) != MP_OKAY) {
214 goto ERR;
215 }
216
217 ERR:
218 mp_clear_multi(&w0, &w1, &w2, &w3, &w4, &a0, &a1, &a2, &tmp1, NULL);
219 return res;
220 }
221
222 #endif