minor
[tinc] / src / ed25519 / sc.c
1 #include "fixedint.h"
2 #include "sc.h"
3
4 static uint64_t load_3(const unsigned char *in) {
5         uint64_t result;
6
7         result = in[0];
8         result |= shlu64(in[1], 8);
9         result |= shlu64(in[2], 16);
10
11         return result;
12 }
13
14 static uint64_t load_4(const unsigned char *in) {
15         uint64_t result;
16
17         result = in[0];
18         result |= shlu64(in[1], 8);
19         result |= shlu64(in[2], 16);
20         result |= shlu64(in[3], 24);
21
22         return result;
23 }
24
25 /*
26 Input:
27   s[0]+256*s[1]+...+256^63*s[63] = s
28
29 Output:
30   s[0]+256*s[1]+...+256^31*s[31] = s mod l
31   where l = 2^252 + 27742317777372353535851937790883648493.
32   Overwrites s in place.
33 */
34
35 void sc_reduce(unsigned char *s) {
36         int64_t s0 = 2097151 & load_3(s);
37         int64_t s1 = 2097151 & (load_4(s + 2) >> 5);
38         int64_t s2 = 2097151 & (load_3(s + 5) >> 2);
39         int64_t s3 = 2097151 & (load_4(s + 7) >> 7);
40         int64_t s4 = 2097151 & (load_4(s + 10) >> 4);
41         int64_t s5 = 2097151 & (load_3(s + 13) >> 1);
42         int64_t s6 = 2097151 & (load_4(s + 15) >> 6);
43         int64_t s7 = 2097151 & (load_3(s + 18) >> 3);
44         int64_t s8 = 2097151 & load_3(s + 21);
45         int64_t s9 = 2097151 & (load_4(s + 23) >> 5);
46         int64_t s10 = 2097151 & (load_3(s + 26) >> 2);
47         int64_t s11 = 2097151 & (load_4(s + 28) >> 7);
48         int64_t s12 = 2097151 & (load_4(s + 31) >> 4);
49         int64_t s13 = 2097151 & (load_3(s + 34) >> 1);
50         int64_t s14 = 2097151 & (load_4(s + 36) >> 6);
51         int64_t s15 = 2097151 & (load_3(s + 39) >> 3);
52         int64_t s16 = 2097151 & load_3(s + 42);
53         int64_t s17 = 2097151 & (load_4(s + 44) >> 5);
54         int64_t s18 = 2097151 & (load_3(s + 47) >> 2);
55         int64_t s19 = 2097151 & (load_4(s + 49) >> 7);
56         int64_t s20 = 2097151 & (load_4(s + 52) >> 4);
57         int64_t s21 = 2097151 & (load_3(s + 55) >> 1);
58         int64_t s22 = 2097151 & (load_4(s + 57) >> 6);
59         int64_t s23 = (load_4(s + 60) >> 3);
60         int64_t carry0;
61         int64_t carry1;
62         int64_t carry2;
63         int64_t carry3;
64         int64_t carry4;
65         int64_t carry5;
66         int64_t carry6;
67         int64_t carry7;
68         int64_t carry8;
69         int64_t carry9;
70         int64_t carry10;
71         int64_t carry11;
72         int64_t carry12;
73         int64_t carry13;
74         int64_t carry14;
75         int64_t carry15;
76         int64_t carry16;
77
78         s11 += s23 * 666643;
79         s12 += s23 * 470296;
80         s13 += s23 * 654183;
81         s14 -= s23 * 997805;
82         s15 += s23 * 136657;
83         s16 -= s23 * 683901;
84         s10 += s22 * 666643;
85         s11 += s22 * 470296;
86         s12 += s22 * 654183;
87         s13 -= s22 * 997805;
88         s14 += s22 * 136657;
89         s15 -= s22 * 683901;
90         s9 += s21 * 666643;
91         s10 += s21 * 470296;
92         s11 += s21 * 654183;
93         s12 -= s21 * 997805;
94         s13 += s21 * 136657;
95         s14 -= s21 * 683901;
96         s8 += s20 * 666643;
97         s9 += s20 * 470296;
98         s10 += s20 * 654183;
99         s11 -= s20 * 997805;
100         s12 += s20 * 136657;
101         s13 -= s20 * 683901;
102         s7 += s19 * 666643;
103         s8 += s19 * 470296;
104         s9 += s19 * 654183;
105         s10 -= s19 * 997805;
106         s11 += s19 * 136657;
107         s12 -= s19 * 683901;
108         s6 += s18 * 666643;
109         s7 += s18 * 470296;
110         s8 += s18 * 654183;
111         s9 -= s18 * 997805;
112         s10 += s18 * 136657;
113         s11 -= s18 * 683901;
114         carry6 = (s6 + (1 << 20)) >> 21;
115         s7 += carry6;
116         s6 -= shl64(carry6, 21);
117         carry8 = (s8 + (1 << 20)) >> 21;
118         s9 += carry8;
119         s8 -= shl64(carry8, 21);
120         carry10 = (s10 + (1 << 20)) >> 21;
121         s11 += carry10;
122         s10 -= shl64(carry10, 21);
123         carry12 = (s12 + (1 << 20)) >> 21;
124         s13 += carry12;
125         s12 -= shl64(carry12, 21);
126         carry14 = (s14 + (1 << 20)) >> 21;
127         s15 += carry14;
128         s14 -= shl64(carry14, 21);
129         carry16 = (s16 + (1 << 20)) >> 21;
130         s17 += carry16;
131         s16 -= shl64(carry16, 21);
132         carry7 = (s7 + (1 << 20)) >> 21;
133         s8 += carry7;
134         s7 -= shl64(carry7, 21);
135         carry9 = (s9 + (1 << 20)) >> 21;
136         s10 += carry9;
137         s9 -= shl64(carry9, 21);
138         carry11 = (s11 + (1 << 20)) >> 21;
139         s12 += carry11;
140         s11 -= shl64(carry11, 21);
141         carry13 = (s13 + (1 << 20)) >> 21;
142         s14 += carry13;
143         s13 -= shl64(carry13, 21);
144         carry15 = (s15 + (1 << 20)) >> 21;
145         s16 += carry15;
146         s15 -= shl64(carry15, 21);
147         s5 += s17 * 666643;
148         s6 += s17 * 470296;
149         s7 += s17 * 654183;
150         s8 -= s17 * 997805;
151         s9 += s17 * 136657;
152         s10 -= s17 * 683901;
153         s4 += s16 * 666643;
154         s5 += s16 * 470296;
155         s6 += s16 * 654183;
156         s7 -= s16 * 997805;
157         s8 += s16 * 136657;
158         s9 -= s16 * 683901;
159         s3 += s15 * 666643;
160         s4 += s15 * 470296;
161         s5 += s15 * 654183;
162         s6 -= s15 * 997805;
163         s7 += s15 * 136657;
164         s8 -= s15 * 683901;
165         s2 += s14 * 666643;
166         s3 += s14 * 470296;
167         s4 += s14 * 654183;
168         s5 -= s14 * 997805;
169         s6 += s14 * 136657;
170         s7 -= s14 * 683901;
171         s1 += s13 * 666643;
172         s2 += s13 * 470296;
173         s3 += s13 * 654183;
174         s4 -= s13 * 997805;
175         s5 += s13 * 136657;
176         s6 -= s13 * 683901;
177         s0 += s12 * 666643;
178         s1 += s12 * 470296;
179         s2 += s12 * 654183;
180         s3 -= s12 * 997805;
181         s4 += s12 * 136657;
182         s5 -= s12 * 683901;
183         s12 = 0;
184         carry0 = (s0 + (1 << 20)) >> 21;
185         s1 += carry0;
186         s0 -= shl64(carry0, 21);
187         carry2 = (s2 + (1 << 20)) >> 21;
188         s3 += carry2;
189         s2 -= shl64(carry2, 21);
190         carry4 = (s4 + (1 << 20)) >> 21;
191         s5 += carry4;
192         s4 -= shl64(carry4, 21);
193         carry6 = (s6 + (1 << 20)) >> 21;
194         s7 += carry6;
195         s6 -= shl64(carry6, 21);
196         carry8 = (s8 + (1 << 20)) >> 21;
197         s9 += carry8;
198         s8 -= shl64(carry8, 21);
199         carry10 = (s10 + (1 << 20)) >> 21;
200         s11 += carry10;
201         s10 -= shl64(carry10, 21);
202         carry1 = (s1 + (1 << 20)) >> 21;
203         s2 += carry1;
204         s1 -= shl64(carry1, 21);
205         carry3 = (s3 + (1 << 20)) >> 21;
206         s4 += carry3;
207         s3 -= shl64(carry3, 21);
208         carry5 = (s5 + (1 << 20)) >> 21;
209         s6 += carry5;
210         s5 -= shl64(carry5, 21);
211         carry7 = (s7 + (1 << 20)) >> 21;
212         s8 += carry7;
213         s7 -= shl64(carry7, 21);
214         carry9 = (s9 + (1 << 20)) >> 21;
215         s10 += carry9;
216         s9 -= shl64(carry9, 21);
217         carry11 = (s11 + (1 << 20)) >> 21;
218         s12 += carry11;
219         s11 -= shl64(carry11, 21);
220         s0 += s12 * 666643;
221         s1 += s12 * 470296;
222         s2 += s12 * 654183;
223         s3 -= s12 * 997805;
224         s4 += s12 * 136657;
225         s5 -= s12 * 683901;
226         s12 = 0;
227         carry0 = s0 >> 21;
228         s1 += carry0;
229         s0 -= shl64(carry0, 21);
230         carry1 = s1 >> 21;
231         s2 += carry1;
232         s1 -= shl64(carry1, 21);
233         carry2 = s2 >> 21;
234         s3 += carry2;
235         s2 -= shl64(carry2, 21);
236         carry3 = s3 >> 21;
237         s4 += carry3;
238         s3 -= shl64(carry3, 21);
239         carry4 = s4 >> 21;
240         s5 += carry4;
241         s4 -= shl64(carry4, 21);
242         carry5 = s5 >> 21;
243         s6 += carry5;
244         s5 -= shl64(carry5, 21);
245         carry6 = s6 >> 21;
246         s7 += carry6;
247         s6 -= shl64(carry6, 21);
248         carry7 = s7 >> 21;
249         s8 += carry7;
250         s7 -= shl64(carry7, 21);
251         carry8 = s8 >> 21;
252         s9 += carry8;
253         s8 -= shl64(carry8, 21);
254         carry9 = s9 >> 21;
255         s10 += carry9;
256         s9 -= shl64(carry9, 21);
257         carry10 = s10 >> 21;
258         s11 += carry10;
259         s10 -= shl64(carry10, 21);
260         carry11 = s11 >> 21;
261         s12 += carry11;
262         s11 -= shl64(carry11, 21);
263         s0 += s12 * 666643;
264         s1 += s12 * 470296;
265         s2 += s12 * 654183;
266         s3 -= s12 * 997805;
267         s4 += s12 * 136657;
268         s5 -= s12 * 683901;
269         carry0 = s0 >> 21;
270         s1 += carry0;
271         s0 -= shl64(carry0, 21);
272         carry1 = s1 >> 21;
273         s2 += carry1;
274         s1 -= shl64(carry1, 21);
275         carry2 = s2 >> 21;
276         s3 += carry2;
277         s2 -= shl64(carry2, 21);
278         carry3 = s3 >> 21;
279         s4 += carry3;
280         s3 -= shl64(carry3, 21);
281         carry4 = s4 >> 21;
282         s5 += carry4;
283         s4 -= shl64(carry4, 21);
284         carry5 = s5 >> 21;
285         s6 += carry5;
286         s5 -= shl64(carry5, 21);
287         carry6 = s6 >> 21;
288         s7 += carry6;
289         s6 -= shl64(carry6, 21);
290         carry7 = s7 >> 21;
291         s8 += carry7;
292         s7 -= shl64(carry7, 21);
293         carry8 = s8 >> 21;
294         s9 += carry8;
295         s8 -= shl64(carry8, 21);
296         carry9 = s9 >> 21;
297         s10 += carry9;
298         s9 -= shl64(carry9, 21);
299         carry10 = s10 >> 21;
300         s11 += carry10;
301         s10 -= shl64(carry10, 21);
302
303         s[0] = (unsigned char)(s0 >> 0);
304         s[1] = (unsigned char)(s0 >> 8);
305         s[2] = (unsigned char)((s0 >> 16) | shl64(s1, 5));
306         s[3] = (unsigned char)(s1 >> 3);
307         s[4] = (unsigned char)(s1 >> 11);
308         s[5] = (unsigned char)((s1 >> 19) | shl64(s2, 2));
309         s[6] = (unsigned char)(s2 >> 6);
310         s[7] = (unsigned char)((s2 >> 14) | shl64(s3, 7));
311         s[8] = (unsigned char)(s3 >> 1);
312         s[9] = (unsigned char)(s3 >> 9);
313         s[10] = (unsigned char)((s3 >> 17) | shl64(s4, 4));
314         s[11] = (unsigned char)(s4 >> 4);
315         s[12] = (unsigned char)(s4 >> 12);
316         s[13] = (unsigned char)((s4 >> 20) | shl64(s5, 1));
317         s[14] = (unsigned char)(s5 >> 7);
318         s[15] = (unsigned char)((s5 >> 15) | shl64(s6, 6));
319         s[16] = (unsigned char)(s6 >> 2);
320         s[17] = (unsigned char)(s6 >> 10);
321         s[18] = (unsigned char)((s6 >> 18) | shl64(s7, 3));
322         s[19] = (unsigned char)(s7 >> 5);
323         s[20] = (unsigned char)(s7 >> 13);
324         s[21] = (unsigned char)(s8 >> 0);
325         s[22] = (unsigned char)(s8 >> 8);
326         s[23] = (unsigned char)((s8 >> 16) | shl64(s9, 5));
327         s[24] = (unsigned char)(s9 >> 3);
328         s[25] = (unsigned char)(s9 >> 11);
329         s[26] = (unsigned char)((s9 >> 19) | shl64(s10, 2));
330         s[27] = (unsigned char)(s10 >> 6);
331         s[28] = (unsigned char)((s10 >> 14) | shl64(s11, 7));
332         s[29] = (unsigned char)(s11 >> 1);
333         s[30] = (unsigned char)(s11 >> 9);
334         s[31] = (unsigned char)(s11 >> 17);
335 }
336
337
338
339 /*
340 Input:
341   a[0]+256*a[1]+...+256^31*a[31] = a
342   b[0]+256*b[1]+...+256^31*b[31] = b
343   c[0]+256*c[1]+...+256^31*c[31] = c
344
345 Output:
346   s[0]+256*s[1]+...+256^31*s[31] = (ab+c) mod l
347   where l = 2^252 + 27742317777372353535851937790883648493.
348 */
349
350 void sc_muladd(unsigned char *s, const unsigned char *a, const unsigned char *b, const unsigned char *c) {
351         int64_t a0 = 2097151 & load_3(a);
352         int64_t a1 = 2097151 & (load_4(a + 2) >> 5);
353         int64_t a2 = 2097151 & (load_3(a + 5) >> 2);
354         int64_t a3 = 2097151 & (load_4(a + 7) >> 7);
355         int64_t a4 = 2097151 & (load_4(a + 10) >> 4);
356         int64_t a5 = 2097151 & (load_3(a + 13) >> 1);
357         int64_t a6 = 2097151 & (load_4(a + 15) >> 6);
358         int64_t a7 = 2097151 & (load_3(a + 18) >> 3);
359         int64_t a8 = 2097151 & load_3(a + 21);
360         int64_t a9 = 2097151 & (load_4(a + 23) >> 5);
361         int64_t a10 = 2097151 & (load_3(a + 26) >> 2);
362         int64_t a11 = (load_4(a + 28) >> 7);
363         int64_t b0 = 2097151 & load_3(b);
364         int64_t b1 = 2097151 & (load_4(b + 2) >> 5);
365         int64_t b2 = 2097151 & (load_3(b + 5) >> 2);
366         int64_t b3 = 2097151 & (load_4(b + 7) >> 7);
367         int64_t b4 = 2097151 & (load_4(b + 10) >> 4);
368         int64_t b5 = 2097151 & (load_3(b + 13) >> 1);
369         int64_t b6 = 2097151 & (load_4(b + 15) >> 6);
370         int64_t b7 = 2097151 & (load_3(b + 18) >> 3);
371         int64_t b8 = 2097151 & load_3(b + 21);
372         int64_t b9 = 2097151 & (load_4(b + 23) >> 5);
373         int64_t b10 = 2097151 & (load_3(b + 26) >> 2);
374         int64_t b11 = (load_4(b + 28) >> 7);
375         int64_t c0 = 2097151 & load_3(c);
376         int64_t c1 = 2097151 & (load_4(c + 2) >> 5);
377         int64_t c2 = 2097151 & (load_3(c + 5) >> 2);
378         int64_t c3 = 2097151 & (load_4(c + 7) >> 7);
379         int64_t c4 = 2097151 & (load_4(c + 10) >> 4);
380         int64_t c5 = 2097151 & (load_3(c + 13) >> 1);
381         int64_t c6 = 2097151 & (load_4(c + 15) >> 6);
382         int64_t c7 = 2097151 & (load_3(c + 18) >> 3);
383         int64_t c8 = 2097151 & load_3(c + 21);
384         int64_t c9 = 2097151 & (load_4(c + 23) >> 5);
385         int64_t c10 = 2097151 & (load_3(c + 26) >> 2);
386         int64_t c11 = (load_4(c + 28) >> 7);
387         int64_t s0;
388         int64_t s1;
389         int64_t s2;
390         int64_t s3;
391         int64_t s4;
392         int64_t s5;
393         int64_t s6;
394         int64_t s7;
395         int64_t s8;
396         int64_t s9;
397         int64_t s10;
398         int64_t s11;
399         int64_t s12;
400         int64_t s13;
401         int64_t s14;
402         int64_t s15;
403         int64_t s16;
404         int64_t s17;
405         int64_t s18;
406         int64_t s19;
407         int64_t s20;
408         int64_t s21;
409         int64_t s22;
410         int64_t s23;
411         int64_t carry0;
412         int64_t carry1;
413         int64_t carry2;
414         int64_t carry3;
415         int64_t carry4;
416         int64_t carry5;
417         int64_t carry6;
418         int64_t carry7;
419         int64_t carry8;
420         int64_t carry9;
421         int64_t carry10;
422         int64_t carry11;
423         int64_t carry12;
424         int64_t carry13;
425         int64_t carry14;
426         int64_t carry15;
427         int64_t carry16;
428         int64_t carry17;
429         int64_t carry18;
430         int64_t carry19;
431         int64_t carry20;
432         int64_t carry21;
433         int64_t carry22;
434
435         s0 = c0 + a0 * b0;
436         s1 = c1 + a0 * b1 + a1 * b0;
437         s2 = c2 + a0 * b2 + a1 * b1 + a2 * b0;
438         s3 = c3 + a0 * b3 + a1 * b2 + a2 * b1 + a3 * b0;
439         s4 = c4 + a0 * b4 + a1 * b3 + a2 * b2 + a3 * b1 + a4 * b0;
440         s5 = c5 + a0 * b5 + a1 * b4 + a2 * b3 + a3 * b2 + a4 * b1 + a5 * b0;
441         s6 = c6 + a0 * b6 + a1 * b5 + a2 * b4 + a3 * b3 + a4 * b2 + a5 * b1 + a6 * b0;
442         s7 = c7 + a0 * b7 + a1 * b6 + a2 * b5 + a3 * b4 + a4 * b3 + a5 * b2 + a6 * b1 + a7 * b0;
443         s8 = c8 + a0 * b8 + a1 * b7 + a2 * b6 + a3 * b5 + a4 * b4 + a5 * b3 + a6 * b2 + a7 * b1 + a8 * b0;
444         s9 = c9 + a0 * b9 + a1 * b8 + a2 * b7 + a3 * b6 + a4 * b5 + a5 * b4 + a6 * b3 + a7 * b2 + a8 * b1 + a9 * b0;
445         s10 = c10 + a0 * b10 + a1 * b9 + a2 * b8 + a3 * b7 + a4 * b6 + a5 * b5 + a6 * b4 + a7 * b3 + a8 * b2 + a9 * b1 + a10 * b0;
446         s11 = c11 + a0 * b11 + a1 * b10 + a2 * b9 + a3 * b8 + a4 * b7 + a5 * b6 + a6 * b5 + a7 * b4 + a8 * b3 + a9 * b2 + a10 * b1 + a11 * b0;
447         s12 = a1 * b11 + a2 * b10 + a3 * b9 + a4 * b8 + a5 * b7 + a6 * b6 + a7 * b5 + a8 * b4 + a9 * b3 + a10 * b2 + a11 * b1;
448         s13 = a2 * b11 + a3 * b10 + a4 * b9 + a5 * b8 + a6 * b7 + a7 * b6 + a8 * b5 + a9 * b4 + a10 * b3 + a11 * b2;
449         s14 = a3 * b11 + a4 * b10 + a5 * b9 + a6 * b8 + a7 * b7 + a8 * b6 + a9 * b5 + a10 * b4 + a11 * b3;
450         s15 = a4 * b11 + a5 * b10 + a6 * b9 + a7 * b8 + a8 * b7 + a9 * b6 + a10 * b5 + a11 * b4;
451         s16 = a5 * b11 + a6 * b10 + a7 * b9 + a8 * b8 + a9 * b7 + a10 * b6 + a11 * b5;
452         s17 = a6 * b11 + a7 * b10 + a8 * b9 + a9 * b8 + a10 * b7 + a11 * b6;
453         s18 = a7 * b11 + a8 * b10 + a9 * b9 + a10 * b8 + a11 * b7;
454         s19 = a8 * b11 + a9 * b10 + a10 * b9 + a11 * b8;
455         s20 = a9 * b11 + a10 * b10 + a11 * b9;
456         s21 = a10 * b11 + a11 * b10;
457         s22 = a11 * b11;
458         s23 = 0;
459         carry0 = (s0 + (1 << 20)) >> 21;
460         s1 += carry0;
461         s0 -= shl64(carry0, 21);
462         carry2 = (s2 + (1 << 20)) >> 21;
463         s3 += carry2;
464         s2 -= shl64(carry2, 21);
465         carry4 = (s4 + (1 << 20)) >> 21;
466         s5 += carry4;
467         s4 -= shl64(carry4, 21);
468         carry6 = (s6 + (1 << 20)) >> 21;
469         s7 += carry6;
470         s6 -= shl64(carry6, 21);
471         carry8 = (s8 + (1 << 20)) >> 21;
472         s9 += carry8;
473         s8 -= shl64(carry8, 21);
474         carry10 = (s10 + (1 << 20)) >> 21;
475         s11 += carry10;
476         s10 -= shl64(carry10, 21);
477         carry12 = (s12 + (1 << 20)) >> 21;
478         s13 += carry12;
479         s12 -= shl64(carry12, 21);
480         carry14 = (s14 + (1 << 20)) >> 21;
481         s15 += carry14;
482         s14 -= shl64(carry14, 21);
483         carry16 = (s16 + (1 << 20)) >> 21;
484         s17 += carry16;
485         s16 -= shl64(carry16, 21);
486         carry18 = (s18 + (1 << 20)) >> 21;
487         s19 += carry18;
488         s18 -= shl64(carry18, 21);
489         carry20 = (s20 + (1 << 20)) >> 21;
490         s21 += carry20;
491         s20 -= shl64(carry20, 21);
492         carry22 = (s22 + (1 << 20)) >> 21;
493         s23 += carry22;
494         s22 -= shl64(carry22, 21);
495         carry1 = (s1 + (1 << 20)) >> 21;
496         s2 += carry1;
497         s1 -= shl64(carry1, 21);
498         carry3 = (s3 + (1 << 20)) >> 21;
499         s4 += carry3;
500         s3 -= shl64(carry3, 21);
501         carry5 = (s5 + (1 << 20)) >> 21;
502         s6 += carry5;
503         s5 -= shl64(carry5, 21);
504         carry7 = (s7 + (1 << 20)) >> 21;
505         s8 += carry7;
506         s7 -= shl64(carry7, 21);
507         carry9 = (s9 + (1 << 20)) >> 21;
508         s10 += carry9;
509         s9 -= shl64(carry9, 21);
510         carry11 = (s11 + (1 << 20)) >> 21;
511         s12 += carry11;
512         s11 -= shl64(carry11, 21);
513         carry13 = (s13 + (1 << 20)) >> 21;
514         s14 += carry13;
515         s13 -= shl64(carry13, 21);
516         carry15 = (s15 + (1 << 20)) >> 21;
517         s16 += carry15;
518         s15 -= shl64(carry15, 21);
519         carry17 = (s17 + (1 << 20)) >> 21;
520         s18 += carry17;
521         s17 -= shl64(carry17, 21);
522         carry19 = (s19 + (1 << 20)) >> 21;
523         s20 += carry19;
524         s19 -= shl64(carry19, 21);
525         carry21 = (s21 + (1 << 20)) >> 21;
526         s22 += carry21;
527         s21 -= shl64(carry21, 21);
528         s11 += s23 * 666643;
529         s12 += s23 * 470296;
530         s13 += s23 * 654183;
531         s14 -= s23 * 997805;
532         s15 += s23 * 136657;
533         s16 -= s23 * 683901;
534         s10 += s22 * 666643;
535         s11 += s22 * 470296;
536         s12 += s22 * 654183;
537         s13 -= s22 * 997805;
538         s14 += s22 * 136657;
539         s15 -= s22 * 683901;
540         s9 += s21 * 666643;
541         s10 += s21 * 470296;
542         s11 += s21 * 654183;
543         s12 -= s21 * 997805;
544         s13 += s21 * 136657;
545         s14 -= s21 * 683901;
546         s8 += s20 * 666643;
547         s9 += s20 * 470296;
548         s10 += s20 * 654183;
549         s11 -= s20 * 997805;
550         s12 += s20 * 136657;
551         s13 -= s20 * 683901;
552         s7 += s19 * 666643;
553         s8 += s19 * 470296;
554         s9 += s19 * 654183;
555         s10 -= s19 * 997805;
556         s11 += s19 * 136657;
557         s12 -= s19 * 683901;
558         s6 += s18 * 666643;
559         s7 += s18 * 470296;
560         s8 += s18 * 654183;
561         s9 -= s18 * 997805;
562         s10 += s18 * 136657;
563         s11 -= s18 * 683901;
564         carry6 = (s6 + (1 << 20)) >> 21;
565         s7 += carry6;
566         s6 -= shl64(carry6, 21);
567         carry8 = (s8 + (1 << 20)) >> 21;
568         s9 += carry8;
569         s8 -= shl64(carry8, 21);
570         carry10 = (s10 + (1 << 20)) >> 21;
571         s11 += carry10;
572         s10 -= shl64(carry10, 21);
573         carry12 = (s12 + (1 << 20)) >> 21;
574         s13 += carry12;
575         s12 -= shl64(carry12, 21);
576         carry14 = (s14 + (1 << 20)) >> 21;
577         s15 += carry14;
578         s14 -= shl64(carry14, 21);
579         carry16 = (s16 + (1 << 20)) >> 21;
580         s17 += carry16;
581         s16 -= shl64(carry16, 21);
582         carry7 = (s7 + (1 << 20)) >> 21;
583         s8 += carry7;
584         s7 -= shl64(carry7, 21);
585         carry9 = (s9 + (1 << 20)) >> 21;
586         s10 += carry9;
587         s9 -= shl64(carry9, 21);
588         carry11 = (s11 + (1 << 20)) >> 21;
589         s12 += carry11;
590         s11 -= shl64(carry11, 21);
591         carry13 = (s13 + (1 << 20)) >> 21;
592         s14 += carry13;
593         s13 -= shl64(carry13, 21);
594         carry15 = (s15 + (1 << 20)) >> 21;
595         s16 += carry15;
596         s15 -= shl64(carry15, 21);
597         s5 += s17 * 666643;
598         s6 += s17 * 470296;
599         s7 += s17 * 654183;
600         s8 -= s17 * 997805;
601         s9 += s17 * 136657;
602         s10 -= s17 * 683901;
603         s4 += s16 * 666643;
604         s5 += s16 * 470296;
605         s6 += s16 * 654183;
606         s7 -= s16 * 997805;
607         s8 += s16 * 136657;
608         s9 -= s16 * 683901;
609         s3 += s15 * 666643;
610         s4 += s15 * 470296;
611         s5 += s15 * 654183;
612         s6 -= s15 * 997805;
613         s7 += s15 * 136657;
614         s8 -= s15 * 683901;
615         s2 += s14 * 666643;
616         s3 += s14 * 470296;
617         s4 += s14 * 654183;
618         s5 -= s14 * 997805;
619         s6 += s14 * 136657;
620         s7 -= s14 * 683901;
621         s1 += s13 * 666643;
622         s2 += s13 * 470296;
623         s3 += s13 * 654183;
624         s4 -= s13 * 997805;
625         s5 += s13 * 136657;
626         s6 -= s13 * 683901;
627         s0 += s12 * 666643;
628         s1 += s12 * 470296;
629         s2 += s12 * 654183;
630         s3 -= s12 * 997805;
631         s4 += s12 * 136657;
632         s5 -= s12 * 683901;
633         s12 = 0;
634         carry0 = (s0 + (1 << 20)) >> 21;
635         s1 += carry0;
636         s0 -= shl64(carry0, 21);
637         carry2 = (s2 + (1 << 20)) >> 21;
638         s3 += carry2;
639         s2 -= shl64(carry2, 21);
640         carry4 = (s4 + (1 << 20)) >> 21;
641         s5 += carry4;
642         s4 -= shl64(carry4, 21);
643         carry6 = (s6 + (1 << 20)) >> 21;
644         s7 += carry6;
645         s6 -= shl64(carry6, 21);
646         carry8 = (s8 + (1 << 20)) >> 21;
647         s9 += carry8;
648         s8 -= shl64(carry8, 21);
649         carry10 = (s10 + (1 << 20)) >> 21;
650         s11 += carry10;
651         s10 -= shl64(carry10, 21);
652         carry1 = (s1 + (1 << 20)) >> 21;
653         s2 += carry1;
654         s1 -= shl64(carry1, 21);
655         carry3 = (s3 + (1 << 20)) >> 21;
656         s4 += carry3;
657         s3 -= shl64(carry3, 21);
658         carry5 = (s5 + (1 << 20)) >> 21;
659         s6 += carry5;
660         s5 -= shl64(carry5, 21);
661         carry7 = (s7 + (1 << 20)) >> 21;
662         s8 += carry7;
663         s7 -= shl64(carry7, 21);
664         carry9 = (s9 + (1 << 20)) >> 21;
665         s10 += carry9;
666         s9 -= shl64(carry9, 21);
667         carry11 = (s11 + (1 << 20)) >> 21;
668         s12 += carry11;
669         s11 -= shl64(carry11, 21);
670         s0 += s12 * 666643;
671         s1 += s12 * 470296;
672         s2 += s12 * 654183;
673         s3 -= s12 * 997805;
674         s4 += s12 * 136657;
675         s5 -= s12 * 683901;
676         s12 = 0;
677         carry0 = s0 >> 21;
678         s1 += carry0;
679         s0 -= shl64(carry0, 21);
680         carry1 = s1 >> 21;
681         s2 += carry1;
682         s1 -= shl64(carry1, 21);
683         carry2 = s2 >> 21;
684         s3 += carry2;
685         s2 -= shl64(carry2, 21);
686         carry3 = s3 >> 21;
687         s4 += carry3;
688         s3 -= shl64(carry3, 21);
689         carry4 = s4 >> 21;
690         s5 += carry4;
691         s4 -= shl64(carry4, 21);
692         carry5 = s5 >> 21;
693         s6 += carry5;
694         s5 -= shl64(carry5, 21);
695         carry6 = s6 >> 21;
696         s7 += carry6;
697         s6 -= shl64(carry6, 21);
698         carry7 = s7 >> 21;
699         s8 += carry7;
700         s7 -= shl64(carry7, 21);
701         carry8 = s8 >> 21;
702         s9 += carry8;
703         s8 -= shl64(carry8, 21);
704         carry9 = s9 >> 21;
705         s10 += carry9;
706         s9 -= shl64(carry9, 21);
707         carry10 = s10 >> 21;
708         s11 += carry10;
709         s10 -= shl64(carry10, 21);
710         carry11 = s11 >> 21;
711         s12 += carry11;
712         s11 -= shl64(carry11, 21);
713         s0 += s12 * 666643;
714         s1 += s12 * 470296;
715         s2 += s12 * 654183;
716         s3 -= s12 * 997805;
717         s4 += s12 * 136657;
718         s5 -= s12 * 683901;
719         carry0 = s0 >> 21;
720         s1 += carry0;
721         s0 -= shl64(carry0, 21);
722         carry1 = s1 >> 21;
723         s2 += carry1;
724         s1 -= shl64(carry1, 21);
725         carry2 = s2 >> 21;
726         s3 += carry2;
727         s2 -= shl64(carry2, 21);
728         carry3 = s3 >> 21;
729         s4 += carry3;
730         s3 -= shl64(carry3, 21);
731         carry4 = s4 >> 21;
732         s5 += carry4;
733         s4 -= shl64(carry4, 21);
734         carry5 = s5 >> 21;
735         s6 += carry5;
736         s5 -= shl64(carry5, 21);
737         carry6 = s6 >> 21;
738         s7 += carry6;
739         s6 -= shl64(carry6, 21);
740         carry7 = s7 >> 21;
741         s8 += carry7;
742         s7 -= shl64(carry7, 21);
743         carry8 = s8 >> 21;
744         s9 += carry8;
745         s8 -= shl64(carry8, 21);
746         carry9 = s9 >> 21;
747         s10 += carry9;
748         s9 -= shl64(carry9, 21);
749         carry10 = s10 >> 21;
750         s11 += carry10;
751         s10 -= shl64(carry10, 21);
752
753         s[0] = (unsigned char)(s0 >> 0);
754         s[1] = (unsigned char)(s0 >> 8);
755         s[2] = (unsigned char)((s0 >> 16) | shl64(s1, 5));
756         s[3] = (unsigned char)(s1 >> 3);
757         s[4] = (unsigned char)(s1 >> 11);
758         s[5] = (unsigned char)((s1 >> 19) | shl64(s2, 2));
759         s[6] = (unsigned char)(s2 >> 6);
760         s[7] = (unsigned char)((s2 >> 14) | shl64(s3, 7));
761         s[8] = (unsigned char)(s3 >> 1);
762         s[9] = (unsigned char)(s3 >> 9);
763         s[10] = (unsigned char)((s3 >> 17) | shl64(s4, 4));
764         s[11] = (unsigned char)(s4 >> 4);
765         s[12] = (unsigned char)(s4 >> 12);
766         s[13] = (unsigned char)((s4 >> 20) | shl64(s5, 1));
767         s[14] = (unsigned char)(s5 >> 7);
768         s[15] = (unsigned char)((s5 >> 15) | shl64(s6, 6));
769         s[16] = (unsigned char)(s6 >> 2);
770         s[17] = (unsigned char)(s6 >> 10);
771         s[18] = (unsigned char)((s6 >> 18) | shl64(s7, 3));
772         s[19] = (unsigned char)(s7 >> 5);
773         s[20] = (unsigned char)(s7 >> 13);
774         s[21] = (unsigned char)(s8 >> 0);
775         s[22] = (unsigned char)(s8 >> 8);
776         s[23] = (unsigned char)((s8 >> 16) | shl64(s9, 5));
777         s[24] = (unsigned char)(s9 >> 3);
778         s[25] = (unsigned char)(s9 >> 11);
779         s[26] = (unsigned char)((s9 >> 19) | shl64(s10, 2));
780         s[27] = (unsigned char)(s10 >> 6);
781         s[28] = (unsigned char)((s10 >> 14) | shl64(s11, 7));
782         s[29] = (unsigned char)(s11 >> 1);
783         s[30] = (unsigned char)(s11 >> 9);
784         s[31] = (unsigned char)(s11 >> 17);
785 }