Import Curve 448 support
[openssl.git] / crypto / ec / curve448 / p448 / arch_arm_32 / f_impl.c
1 /* Copyright (c) 2014 Cryptography Research, Inc.
2  * Released under the MIT License.  See LICENSE.txt for license information.
3  */
4
5 #include "f_field.h"
6
7 static inline void __attribute__((gnu_inline,always_inline))
8 smlal (
9     uint64_t *acc,
10     const uint32_t a,
11     const uint32_t b
12 ) {
13
14 #ifdef  __ARMEL__
15     uint32_t lo = *acc, hi = (*acc)>>32;
16     
17     __asm__ __volatile__ ("smlal %[lo], %[hi], %[a], %[b]"
18         : [lo]"+&r"(lo), [hi]"+&r"(hi)
19         : [a]"r"(a), [b]"r"(b));
20     
21     *acc = lo + (((uint64_t)hi)<<32);
22 #else
23     *acc += (int64_t)(int32_t)a * (int64_t)(int32_t)b;
24 #endif
25 }
26
27 static inline void __attribute__((gnu_inline,always_inline))
28 smlal2 (
29     uint64_t *acc,
30     const uint32_t a,
31     const uint32_t b
32 ) {
33 #ifdef __ARMEL__
34     uint32_t lo = *acc, hi = (*acc)>>32;
35     
36     __asm__ __volatile__ ("smlal %[lo], %[hi], %[a], %[b]"
37         : [lo]"+&r"(lo), [hi]"+&r"(hi)
38         : [a]"r"(a), [b]"r"(2*b));
39     
40     *acc = lo + (((uint64_t)hi)<<32);
41 #else
42     *acc += (int64_t)(int32_t)a * (int64_t)(int32_t)(b * 2);
43 #endif
44 }
45
46 static inline void __attribute__((gnu_inline,always_inline))
47 smull (
48     uint64_t *acc,
49     const uint32_t a,
50     const uint32_t b
51 ) {
52 #ifdef __ARMEL__
53     uint32_t lo, hi;
54     
55     __asm__ __volatile__ ("smull %[lo], %[hi], %[a], %[b]"
56         : [lo]"=&r"(lo), [hi]"=&r"(hi)
57         : [a]"r"(a), [b]"r"(b));
58     
59     *acc = lo + (((uint64_t)hi)<<32);
60 #else
61     *acc = (int64_t)(int32_t)a * (int64_t)(int32_t)b;
62 #endif
63 }
64
65 static inline void __attribute__((gnu_inline,always_inline))
66 smull2 (
67     uint64_t *acc,
68     const uint32_t a,
69     const uint32_t b
70 ) {
71 #ifdef __ARMEL__
72     uint32_t lo, hi;
73     
74     __asm__ /*__volatile__*/ ("smull %[lo], %[hi], %[a], %[b]"
75         : [lo]"=&r"(lo), [hi]"=&r"(hi)
76         : [a]"r"(a), [b]"r"(2*b));
77     
78     *acc = lo + (((uint64_t)hi)<<32);
79 #else
80     *acc = (int64_t)(int32_t)a * (int64_t)(int32_t)(b * 2);
81 #endif
82 }
83
84 void gf_mul (gf_s *__restrict__ cs, const gf as, const gf bs) {
85     
86     const uint32_t *a = as->limb, *b = bs->limb;
87     uint32_t *c = cs->limb;
88
89     uint64_t accum0 = 0, accum1 = 0, accum2, accum3, accumC0, accumC1;
90     uint32_t mask = (1<<28) - 1;  
91
92     uint32_t aa[8], bm[8];
93
94     int i;
95     for (i=0; i<8; i++) {
96         aa[i] = a[i] + a[i+8];
97         bm[i] = b[i] - b[i+8];
98     }
99
100     uint32_t ax,bx;
101     {
102         /* t^3 terms */
103         smull(&accum1, ax = aa[1], bx = b[15]);
104         smull(&accum3, ax = aa[2], bx);
105         smlal(&accum1, ax, bx = b[14]);
106         smlal(&accum3, ax = aa[3], bx);
107         smlal(&accum1, ax, bx = b[13]);
108         smlal(&accum3, ax = aa[4], bx);
109         smlal(&accum1, ax, bx = b[12]);
110         smlal(&accum3, ax = aa[5], bx);
111         smlal(&accum1, ax, bx = b[11]);
112         smlal(&accum3, ax = aa[6], bx);
113         smlal(&accum1, ax, bx = b[10]);
114         smlal(&accum3, ax = aa[7], bx);
115         smlal(&accum1, ax, bx = b[9]);
116         
117         accum0 = accum1;
118         accum2 = accum3;
119         
120         /* t^2 terms */
121         smlal(&accum2, ax = aa[0], bx);
122         smlal(&accum0, ax, bx = b[8]);
123         smlal(&accum2, ax = aa[1], bx);
124         
125         smlal(&accum0, ax = a[9], bx = b[7]);
126         smlal(&accum2, ax = a[10], bx);
127         smlal(&accum0, ax, bx = b[6]);
128         smlal(&accum2, ax = a[11], bx);
129         smlal(&accum0, ax, bx = b[5]);
130         smlal(&accum2, ax = a[12], bx);
131         smlal(&accum0, ax, bx = b[4]);
132         smlal(&accum2, ax = a[13], bx);
133         smlal(&accum0, ax, bx = b[3]);
134         smlal(&accum2, ax = a[14], bx);
135         smlal(&accum0, ax, bx = b[2]);
136         smlal(&accum2, ax = a[15], bx);
137         smlal(&accum0, ax, bx = b[1]);
138         
139         /* t terms */
140         accum1 += accum0;
141         accum3 += accum2;
142         smlal(&accum3, ax = a[8], bx);
143         smlal(&accum1, ax, bx = b[0]);
144         smlal(&accum3, ax = a[9], bx);
145         
146         smlal(&accum1, ax = a[1], bx = bm[7]);
147         smlal(&accum3, ax = a[2], bx);
148         smlal(&accum1, ax, bx = bm[6]);
149         smlal(&accum3, ax = a[3], bx);
150         smlal(&accum1, ax, bx = bm[5]);
151         smlal(&accum3, ax = a[4], bx);
152         smlal(&accum1, ax, bx = bm[4]);
153         smlal(&accum3, ax = a[5], bx);
154         smlal(&accum1, ax, bx = bm[3]);
155         smlal(&accum3, ax = a[6], bx);
156         smlal(&accum1, ax, bx = bm[2]);
157         smlal(&accum3, ax = a[7], bx);
158         smlal(&accum1, ax, bx = bm[1]);
159         
160         /* 1 terms */
161         smlal(&accum2, ax = a[0], bx);
162         smlal(&accum0, ax, bx = bm[0]);
163         smlal(&accum2, ax = a[1], bx);
164         
165         accum2 += accum0 >> 28;
166         accum3 += accum1 >> 28;
167         
168         c[0] = ((uint32_t)(accum0)) & mask;
169         c[1] = ((uint32_t)(accum2)) & mask;
170         c[8] = ((uint32_t)(accum1)) & mask;
171         c[9] = ((uint32_t)(accum3)) & mask;
172         
173         accumC0 = accum2 >> 28;
174         accumC1 = accum3 >> 28;
175     }
176     {
177         /* t^3 terms */
178         smull(&accum1, ax = aa[3], bx = b[15]);
179         smull(&accum3, ax = aa[4], bx);
180         smlal(&accum1, ax, bx = b[14]);
181         smlal(&accum3, ax = aa[5], bx);
182         smlal(&accum1, ax, bx = b[13]);
183         smlal(&accum3, ax = aa[6], bx);
184         smlal(&accum1, ax, bx = b[12]);
185         smlal(&accum3, ax = aa[7], bx);
186         smlal(&accum1, ax, bx = b[11]);
187         
188         accum0 = accum1;
189         accum2 = accum3;
190         
191         /* t^2 terms */
192         smlal(&accum2, ax = aa[0], bx);
193         smlal(&accum0, ax, bx = b[10]);
194         smlal(&accum2, ax = aa[1], bx);
195         smlal(&accum0, ax, bx = b[9]);
196         smlal(&accum2, ax = aa[2], bx);
197         smlal(&accum0, ax, bx = b[8]);
198         smlal(&accum2, ax = aa[3], bx);
199         
200         smlal(&accum0, ax = a[11], bx = b[7]);
201         smlal(&accum2, ax = a[12], bx);
202         smlal(&accum0, ax, bx = b[6]);
203         smlal(&accum2, ax = a[13], bx);
204         smlal(&accum0, ax, bx = b[5]);
205         smlal(&accum2, ax = a[14], bx);
206         smlal(&accum0, ax, bx = b[4]);
207         smlal(&accum2, ax = a[15], bx);
208         smlal(&accum0, ax, bx = b[3]);
209         
210         /* t terms */
211         accum1 += accum0;
212         accum3 += accum2;
213         smlal(&accum3, ax = a[8], bx);
214         smlal(&accum1, ax, bx = b[2]);
215         smlal(&accum3, ax = a[9], bx);
216         smlal(&accum1, ax, bx = b[1]);
217         smlal(&accum3, ax = a[10], bx);
218         smlal(&accum1, ax, bx = b[0]);
219         smlal(&accum3, ax = a[11], bx);
220         
221         smlal(&accum1, ax = a[3], bx = bm[7]);
222         smlal(&accum3, ax = a[4], bx);
223         smlal(&accum1, ax, bx = bm[6]);
224         smlal(&accum3, ax = a[5], bx);
225         smlal(&accum1, ax, bx = bm[5]);
226         smlal(&accum3, ax = a[6], bx);
227         smlal(&accum1, ax, bx = bm[4]);
228         smlal(&accum3, ax = a[7], bx);
229         smlal(&accum1, ax, bx = bm[3]);
230         
231         /* 1 terms */
232         smlal(&accum2, ax = a[0], bx);
233         smlal(&accum0, ax, bx = bm[2]);
234         smlal(&accum2, ax = a[1], bx);
235         smlal(&accum0, ax, bx = bm[1]);
236         smlal(&accum2, ax = a[2], bx);
237         smlal(&accum0, ax, bx = bm[0]);
238         smlal(&accum2, ax = a[3], bx);
239         
240         accum0 += accumC0;
241         accum1 += accumC1;
242         accum2 += accum0 >> 28;
243         accum3 += accum1 >> 28;
244         
245         c[2] = ((uint32_t)(accum0)) & mask;
246         c[3] = ((uint32_t)(accum2)) & mask;
247         c[10] = ((uint32_t)(accum1)) & mask;
248         c[11] = ((uint32_t)(accum3)) & mask;
249         
250         accumC0 = accum2 >> 28;
251         accumC1 = accum3 >> 28;
252     }
253     {
254         
255         /* t^3 terms */
256         smull(&accum1, ax = aa[5], bx = b[15]);
257         smull(&accum3, ax = aa[6], bx);
258         smlal(&accum1, ax, bx = b[14]);
259         smlal(&accum3, ax = aa[7], bx);
260         smlal(&accum1, ax, bx = b[13]);
261         
262         accum0 = accum1;
263         accum2 = accum3;
264         
265         /* t^2 terms */
266         
267         smlal(&accum2, ax = aa[0], bx);
268         smlal(&accum0, ax, bx = b[12]);
269         smlal(&accum2, ax = aa[1], bx);
270         smlal(&accum0, ax, bx = b[11]);
271         smlal(&accum2, ax = aa[2], bx);
272         smlal(&accum0, ax, bx = b[10]);
273         smlal(&accum2, ax = aa[3], bx);
274         smlal(&accum0, ax, bx = b[9]);
275         smlal(&accum2, ax = aa[4], bx);
276         smlal(&accum0, ax, bx = b[8]);
277         smlal(&accum2, ax = aa[5], bx);
278         
279         
280         smlal(&accum0, ax = a[13], bx = b[7]);
281         smlal(&accum2, ax = a[14], bx);
282         smlal(&accum0, ax, bx = b[6]);
283         smlal(&accum2, ax = a[15], bx);
284         smlal(&accum0, ax, bx = b[5]);
285         
286         /* t terms */
287         accum1 += accum0;
288         accum3 += accum2;
289         
290         smlal(&accum3, ax = a[8], bx);
291         smlal(&accum1, ax, bx = b[4]);
292         smlal(&accum3, ax = a[9], bx);
293         smlal(&accum1, ax, bx = b[3]);
294         smlal(&accum3, ax = a[10], bx);
295         smlal(&accum1, ax, bx = b[2]);
296         smlal(&accum3, ax = a[11], bx);
297         smlal(&accum1, ax, bx = b[1]);
298         smlal(&accum3, ax = a[12], bx);
299         smlal(&accum1, ax, bx = b[0]);
300         smlal(&accum3, ax = a[13], bx);
301         
302         
303         smlal(&accum1, ax = a[5], bx = bm[7]);
304         smlal(&accum3, ax = a[6], bx);
305         smlal(&accum1, ax, bx = bm[6]);
306         smlal(&accum3, ax = a[7], bx);
307         smlal(&accum1, ax, bx = bm[5]);
308         
309         /* 1 terms */
310         
311         smlal(&accum2, ax = a[0], bx);
312         smlal(&accum0, ax, bx = bm[4]);
313         smlal(&accum2, ax = a[1], bx);
314         smlal(&accum0, ax, bx = bm[3]);
315         smlal(&accum2, ax = a[2], bx);
316         smlal(&accum0, ax, bx = bm[2]);
317         smlal(&accum2, ax = a[3], bx);
318         smlal(&accum0, ax, bx = bm[1]);
319         smlal(&accum2, ax = a[4], bx);
320         smlal(&accum0, ax, bx = bm[0]);
321         smlal(&accum2, ax = a[5], bx);
322         
323         accum0 += accumC0;
324         accum1 += accumC1;
325         accum2 += accum0 >> 28;
326         accum3 += accum1 >> 28;
327         
328         c[4] = ((uint32_t)(accum0)) & mask;
329         c[5] = ((uint32_t)(accum2)) & mask;
330         c[12] = ((uint32_t)(accum1)) & mask;
331         c[13] = ((uint32_t)(accum3)) & mask;
332         
333         accumC0 = accum2 >> 28;
334         accumC1 = accum3 >> 28;
335     }
336     {
337         
338         /* t^3 terms */
339         smull(&accum1, ax = aa[7], bx = b[15]);
340         accum0 = accum1;
341         
342         /* t^2 terms */
343         
344         smull(&accum2, ax = aa[0], bx);
345         smlal(&accum0, ax, bx = b[14]);
346         smlal(&accum2, ax = aa[1], bx);
347         smlal(&accum0, ax, bx = b[13]);
348         smlal(&accum2, ax = aa[2], bx);
349         smlal(&accum0, ax, bx = b[12]);
350         smlal(&accum2, ax = aa[3], bx);
351         smlal(&accum0, ax, bx = b[11]);
352         smlal(&accum2, ax = aa[4], bx);
353         smlal(&accum0, ax, bx = b[10]);
354         smlal(&accum2, ax = aa[5], bx);
355         smlal(&accum0, ax, bx = b[9]);
356         smlal(&accum2, ax = aa[6], bx);
357         smlal(&accum0, ax, bx = b[8]);
358         smlal(&accum2, ax = aa[7], bx);
359         
360         
361         smlal(&accum0, ax = a[15], bx = b[7]);
362         
363         /* t terms */
364         accum1 += accum0;
365         accum3 = accum2;
366         
367         smlal(&accum3, ax = a[8], bx);
368         smlal(&accum1, ax, bx = b[6]);
369         smlal(&accum3, ax = a[9], bx);
370         smlal(&accum1, ax, bx = b[5]);
371         smlal(&accum3, ax = a[10], bx);
372         smlal(&accum1, ax, bx = b[4]);
373         smlal(&accum3, ax = a[11], bx);
374         smlal(&accum1, ax, bx = b[3]);
375         smlal(&accum3, ax = a[12], bx);
376         smlal(&accum1, ax, bx = b[2]);
377         smlal(&accum3, ax = a[13], bx);
378         smlal(&accum1, ax, bx = b[1]);
379         smlal(&accum3, ax = a[14], bx);
380         smlal(&accum1, ax, bx = b[0]);
381         smlal(&accum3, ax = a[15], bx);
382         
383         
384         smlal(&accum1, ax = a[7], bx = bm[7]);
385         
386         /* 1 terms */
387         
388         smlal(&accum2, ax = a[0], bx);
389         smlal(&accum0, ax, bx = bm[6]);
390         smlal(&accum2, ax = a[1], bx);
391         smlal(&accum0, ax, bx = bm[5]);
392         smlal(&accum2, ax = a[2], bx);
393         smlal(&accum0, ax, bx = bm[4]);
394         smlal(&accum2, ax = a[3], bx);
395         smlal(&accum0, ax, bx = bm[3]);
396         smlal(&accum2, ax = a[4], bx);
397         smlal(&accum0, ax, bx = bm[2]);
398         smlal(&accum2, ax = a[5], bx);
399         smlal(&accum0, ax, bx = bm[1]);
400         smlal(&accum2, ax = a[6], bx);
401         smlal(&accum0, ax, bx = bm[0]);
402         smlal(&accum2, ax = a[7], bx);
403         
404         accum0 += accumC0;
405         accum1 += accumC1;
406         accum2 += accum0 >> 28;
407         accum3 += accum1 >> 28;
408         
409         c[6] = ((uint32_t)(accum0)) & mask;
410         c[7] = ((uint32_t)(accum2)) & mask;
411         c[14] = ((uint32_t)(accum1)) & mask;
412         c[15] = ((uint32_t)(accum3)) & mask;
413         
414         accum0 = accum2 >> 28;
415         accum1 = accum3 >> 28;
416     }
417
418     accum0 += accum1;
419     accum0 += c[8];
420     accum1 += c[0];
421     c[8] = ((uint32_t)(accum0)) & mask;
422     c[0] = ((uint32_t)(accum1)) & mask;
423     
424     accum0 >>= 28;
425     accum1 >>= 28;
426     c[9] += ((uint32_t)(accum0));
427     c[1] += ((uint32_t)(accum1));
428 }
429
430 void gf_sqr (gf_s *__restrict__ cs, const gf as) {
431     const uint32_t *a = as->limb;
432     uint32_t *c = cs->limb;
433
434     uint64_t accum0 = 0, accum1 = 0, accum2, accum3, accumC0, accumC1, tmp;
435     uint32_t mask = (1<<28) - 1;  
436
437     uint32_t bm[8];
438     
439     int i;
440     for (i=0; i<8; i++) {
441         bm[i] = a[i] - a[i+8];
442     }
443
444     uint32_t ax,bx;
445     {
446         /* t^3 terms */
447         smull2(&accum1, ax = a[9], bx = a[15]);
448         smull2(&accum3, ax = a[10], bx);
449         smlal2(&accum1, ax, bx = a[14]);
450         smlal2(&accum3, ax = a[11], bx);
451         smlal2(&accum1, ax, bx = a[13]);
452         smlal2(&accum3, ax = a[12], bx);
453         smlal(&accum1, ax, ax);
454         
455         accum0 = accum1;
456         accum2 = accum3;
457         
458         /* t^2 terms */
459         smlal2(&accum2, ax = a[8], a[9]);
460         smlal(&accum0, ax, ax);
461         
462         smlal2(&accum0, ax = a[1], bx = a[7]);
463         smlal2(&accum2, ax = a[2], bx);
464         smlal2(&accum0, ax, bx = a[6]);
465         smlal2(&accum2, ax = a[3], bx);
466         smlal2(&accum0, ax, bx = a[5]);
467         smlal2(&accum2, ax = a[4], bx);
468         smlal(&accum0, ax, ax);
469         
470         /* t terms */
471         accum1 += accum0;
472         accum3 += accum2;
473         smlal2(&accum3, ax = a[0], bx = a[1]);
474         smlal(&accum1, ax, ax);
475         
476         accum1 = -accum1;
477         accum3 = -accum3;
478         accum2 = -accum2;
479         accum0 = -accum0;
480         
481         smlal2(&accum1, ax = bm[1], bx = bm[7]);
482         smlal2(&accum3, ax = bm[2], bx);
483         smlal2(&accum1, ax, bx = bm[6]);
484         smlal2(&accum3, ax = bm[3], bx);
485         smlal2(&accum1, ax, bx = bm[5]);
486         smlal2(&accum3, ax = bm[4], bx);
487         smlal(&accum1, ax, ax);
488         
489         /* 1 terms */
490         smlal2(&accum2, ax = bm[0], bx = bm[1]);
491         smlal(&accum0, ax, ax);
492         
493         tmp = -accum3; accum3 = tmp-accum2; accum2 = tmp;
494         tmp = -accum1; accum1 = tmp-accum0; accum0 = tmp;
495         
496         accum2 += accum0 >> 28;
497         accum3 += accum1 >> 28;
498         
499         c[0] = ((uint32_t)(accum0)) & mask;
500         c[1] = ((uint32_t)(accum2)) & mask;
501         c[8] = ((uint32_t)(accum1)) & mask;
502         c[9] = ((uint32_t)(accum3)) & mask;
503         
504         accumC0 = accum2 >> 28;
505         accumC1 = accum3 >> 28;
506     }
507     {
508         /* t^3 terms */
509         smull2(&accum1, ax = a[11], bx = a[15]);
510         smull2(&accum3, ax = a[12], bx);
511         smlal2(&accum1, ax, bx = a[14]);
512         smlal2(&accum3, ax = a[13], bx);
513         smlal(&accum1, ax, ax);
514         
515         accum0 = accum1;
516         accum2 = accum3;
517         
518         /* t^2 terms */
519         smlal2(&accum2, ax = a[8], bx = a[11]);
520         smlal2(&accum0, ax, bx = a[10]);
521         smlal2(&accum2, ax = a[9], bx);
522         smlal(&accum0, ax, ax);
523         
524         smlal2(&accum0, ax = a[3], bx = a[7]);
525         smlal2(&accum2, ax = a[4], bx);
526         smlal2(&accum0, ax, bx = a[6]);
527         smlal2(&accum2, ax = a[5], bx);
528         smlal(&accum0, ax, ax);
529         
530         /* t terms */
531         accum1 += accum0;
532         accum3 += accum2;
533         smlal2(&accum3, ax = a[0], bx = a[3]);
534         smlal2(&accum1, ax, bx = a[2]);
535         smlal2(&accum3, ax = a[1], bx);
536         smlal(&accum1, ax, ax);
537         
538         accum1 = -accum1;
539         accum3 = -accum3;
540         accum2 = -accum2;
541         accum0 = -accum0;
542         
543         smlal2(&accum1, ax = bm[3], bx = bm[7]);
544         smlal2(&accum3, ax = bm[4], bx);
545         smlal2(&accum1, ax, bx = bm[6]);
546         smlal2(&accum3, ax = bm[5], bx);
547         smlal(&accum1, ax, ax);
548         
549         /* 1 terms */
550         smlal2(&accum2, ax = bm[0], bx = bm[3]);
551         smlal2(&accum0, ax, bx = bm[2]);
552         smlal2(&accum2, ax = bm[1], bx);
553         smlal(&accum0, ax, ax);
554         
555         
556         tmp = -accum3; accum3 = tmp-accum2; accum2 = tmp;
557         tmp = -accum1; accum1 = tmp-accum0; accum0 = tmp;
558         
559         accum0 += accumC0;
560         accum1 += accumC1;
561         accum2 += accum0 >> 28;
562         accum3 += accum1 >> 28;
563         
564         c[2] = ((uint32_t)(accum0)) & mask;
565         c[3] = ((uint32_t)(accum2)) & mask;
566         c[10] = ((uint32_t)(accum1)) & mask;
567         c[11] = ((uint32_t)(accum3)) & mask;
568         
569         accumC0 = accum2 >> 28;
570         accumC1 = accum3 >> 28;
571     }
572     {
573         
574         /* t^3 terms */
575         smull2(&accum1, ax = a[13], bx = a[15]);
576         smull2(&accum3, ax = a[14], bx);
577         smlal(&accum1, ax, ax);
578         
579         accum0 = accum1;
580         accum2 = accum3;
581         
582         /* t^2 terms */
583         
584         smlal2(&accum2, ax = a[8], bx = a[13]);
585         smlal2(&accum0, ax, bx = a[12]);
586         smlal2(&accum2, ax = a[9], bx);
587         smlal2(&accum0, ax, bx = a[11]);
588         smlal2(&accum2, ax = a[10], bx);
589         smlal(&accum0, ax, ax);
590         
591         
592         smlal2(&accum0, ax = a[5], bx = a[7]);
593         smlal2(&accum2, ax = a[6], bx);
594         smlal(&accum0, ax, ax);
595         
596         /* t terms */
597         accum1 += accum0;
598         accum3 += accum2;
599         
600         smlal2(&accum3, ax = a[0], bx = a[5]);
601         smlal2(&accum1, ax, bx = a[4]);
602         smlal2(&accum3, ax = a[1], bx);
603         smlal2(&accum1, ax, bx = a[3]);
604         smlal2(&accum3, ax = a[2], bx);
605         smlal(&accum1, ax, ax);
606         
607         accum1 = -accum1;
608         accum3 = -accum3;
609         accum2 = -accum2;
610         accum0 = -accum0;
611         
612         smlal2(&accum1, ax = bm[5], bx = bm[7]);
613         smlal2(&accum3, ax = bm[6], bx);
614         smlal(&accum1, ax, ax);
615         
616         /* 1 terms */
617         
618         smlal2(&accum2, ax = bm[0], bx = bm[5]);
619         smlal2(&accum0, ax, bx = bm[4]);
620         smlal2(&accum2, ax = bm[1], bx);
621         smlal2(&accum0, ax, bx = bm[3]);
622         smlal2(&accum2, ax = bm[2], bx);
623         smlal(&accum0, ax, ax);
624         
625         
626         tmp = -accum3; accum3 = tmp-accum2; accum2 = tmp;
627         tmp = -accum1; accum1 = tmp-accum0; accum0 = tmp;
628         
629         accum0 += accumC0;
630         accum1 += accumC1;
631         accum2 += accum0 >> 28;
632         accum3 += accum1 >> 28;
633         
634         c[4] = ((uint32_t)(accum0)) & mask;
635         c[5] = ((uint32_t)(accum2)) & mask;
636         c[12] = ((uint32_t)(accum1)) & mask;
637         c[13] = ((uint32_t)(accum3)) & mask;
638         
639         accumC0 = accum2 >> 28;
640         accumC1 = accum3 >> 28;
641     }
642     {
643         
644         /* t^3 terms */
645         smull(&accum1, ax = a[15], bx = a[15]);
646         accum0 = accum1;
647         
648         /* t^2 terms */
649         
650         smull2(&accum2, ax = a[8], bx);
651         smlal2(&accum0, ax, bx = a[14]);
652         smlal2(&accum2, ax = a[9], bx);
653         smlal2(&accum0, ax, bx = a[13]);
654         smlal2(&accum2, ax = a[10], bx);
655         smlal2(&accum0, ax, bx = a[12]);
656         smlal2(&accum2, ax = a[11], bx);
657         smlal(&accum0, ax, ax);
658         
659         
660         smlal(&accum0, ax = a[7], bx = a[7]);
661         
662         /* t terms */
663         accum1 += accum0;
664         accum3 = accum2;
665         
666         smlal2(&accum3, ax = a[0], bx);
667         smlal2(&accum1, ax, bx = a[6]);
668         smlal2(&accum3, ax = a[1], bx);
669         smlal2(&accum1, ax, bx = a[5]);
670         smlal2(&accum3, ax = a[2], bx);
671         smlal2(&accum1, ax, bx = a[4]);
672         smlal2(&accum3, ax = a[3], bx);
673         smlal(&accum1, ax, ax);
674         
675         accum1 = -accum1;
676         accum3 = -accum3;
677         accum2 = -accum2;
678         accum0 = -accum0;
679         
680         bx = bm[7];
681         smlal(&accum1, bx, bx);
682         
683         /* 1 terms */
684         
685         smlal2(&accum2, ax = bm[0], bx);
686         smlal2(&accum0, ax, bx = bm[6]);
687         smlal2(&accum2, ax = bm[1], bx);
688         smlal2(&accum0, ax, bx = bm[5]);
689         smlal2(&accum2, ax = bm[2], bx);
690         smlal2(&accum0, ax, bx = bm[4]);
691         smlal2(&accum2, ax = bm[3], bx);
692         smlal(&accum0, ax, ax);
693         
694         tmp = -accum3; accum3 = tmp-accum2; accum2 = tmp;
695         tmp = -accum1; accum1 = tmp-accum0; accum0 = tmp;
696         
697         
698         accum0 += accumC0;
699         accum1 += accumC1;
700         accum2 += accum0 >> 28;
701         accum3 += accum1 >> 28;
702         
703         c[6] = ((uint32_t)(accum0)) & mask;
704         c[7] = ((uint32_t)(accum2)) & mask;
705         c[14] = ((uint32_t)(accum1)) & mask;
706         c[15] = ((uint32_t)(accum3)) & mask;
707         
708         accum0 = accum2 >> 28;
709         accum1 = accum3 >> 28;
710     }
711
712     accum0 += accum1;
713     accum0 += c[8];
714     accum1 += c[0];
715     c[8] = ((uint32_t)(accum0)) & mask;
716     c[0] = ((uint32_t)(accum1)) & mask;
717     
718     accum0 >>= 28;
719     accum1 >>= 28;
720     c[9] += ((uint32_t)(accum0));
721     c[1] += ((uint32_t)(accum1));
722 }
723
724 void gf_mulw_unsigned (
725     gf_s *__restrict__ cs,
726     const gf as,
727     uint32_t b
728 ) {
729     uint32_t mask = (1ull<<28)-1;  
730     assert(b <= mask);
731     
732     const uint32_t *a = as->limb;
733     uint32_t *c = cs->limb;
734
735     uint64_t accum0, accum8;
736
737     int i;
738
739     uint32_t c0, c8, n0, n8;
740     c0 = a[0]; c8 = a[8];
741     accum0 = widemul(b, c0);
742     accum8 = widemul(b, c8);
743
744     c[0] = accum0 & mask; accum0 >>= 28;
745     c[8] = accum8 & mask; accum8 >>= 28;
746     
747     i=1;
748     {
749         n0 = a[i]; n8 = a[i+8];
750         smlal(&accum0, b, n0);
751         smlal(&accum8, b, n8);
752         
753         c[i] = accum0 & mask; accum0 >>= 28;
754         c[i+8] = accum8 & mask; accum8 >>= 28;
755         i++;
756     }
757     {
758         c0 = a[i]; c8 = a[i+8];
759         smlal(&accum0, b, c0);
760         smlal(&accum8, b, c8);
761
762         c[i] = accum0 & mask; accum0 >>= 28;
763         c[i+8] = accum8 & mask; accum8 >>= 28;
764         i++;
765     }
766     {
767         n0 = a[i]; n8 = a[i+8];
768         smlal(&accum0, b, n0);
769         smlal(&accum8, b, n8);
770
771         c[i] = accum0 & mask; accum0 >>= 28;
772         c[i+8] = accum8 & mask; accum8 >>= 28;
773         i++;
774     }
775     {
776         c0 = a[i]; c8 = a[i+8];
777         smlal(&accum0, b, c0);
778         smlal(&accum8, b, c8);
779
780         c[i] = accum0 & mask; accum0 >>= 28;
781         c[i+8] = accum8 & mask; accum8 >>= 28;
782         i++;
783     }
784     {
785         n0 = a[i]; n8 = a[i+8];
786         smlal(&accum0, b, n0);
787         smlal(&accum8, b, n8);
788
789         c[i] = accum0 & mask; accum0 >>= 28;
790         c[i+8] = accum8 & mask; accum8 >>= 28;
791         i++;
792     }
793     {
794         c0 = a[i]; c8 = a[i+8];
795         smlal(&accum0, b, c0);
796         smlal(&accum8, b, c8);
797         
798         c[i] = accum0 & mask; accum0 >>= 28;
799         c[i+8] = accum8 & mask; accum8 >>= 28;
800         i++;
801     }
802     {
803         n0 = a[i]; n8 = a[i+8];
804         smlal(&accum0, b, n0);
805         smlal(&accum8, b, n8);
806
807         c[i] = accum0 & mask; accum0 >>= 28;
808         c[i+8] = accum8 & mask; accum8 >>= 28;
809         i++;
810     }
811
812     accum0 += accum8 + c[8];
813     c[8] = accum0 & mask;
814     c[9] += accum0 >> 28;
815
816     accum8 += c[0];
817     c[0] = accum8 & mask;
818     c[1] += accum8 >> 28;
819 }