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