poly1305/asm/poly1305-x86_64.pl: minor AVX512 optimization.
[openssl.git] / crypto / poly1305 / poly1305.c
1 /*
2  * Copyright 2015-2016 The OpenSSL Project Authors. All Rights Reserved.
3  *
4  * Licensed under the OpenSSL license (the "License").  You may not use
5  * this file except in compliance with the License.  You can obtain a copy
6  * in the file LICENSE in the source distribution or at
7  * https://www.openssl.org/source/license.html
8  */
9
10 #include <stdlib.h>
11 #include <string.h>
12 #include <openssl/crypto.h>
13
14 #include "internal/poly1305.h"
15 #include "poly1305_local.h"
16
17 size_t Poly1305_ctx_size ()
18 {
19     return sizeof(struct poly1305_context);
20 }
21
22 /* pick 32-bit unsigned integer in little endian order */
23 static unsigned int U8TOU32(const unsigned char *p)
24 {
25     return (((unsigned int)(p[0] & 0xff)) |
26             ((unsigned int)(p[1] & 0xff) << 8) |
27             ((unsigned int)(p[2] & 0xff) << 16) |
28             ((unsigned int)(p[3] & 0xff) << 24));
29 }
30
31 /*
32  * Implementations can be classified by amount of significant bits in
33  * words making up the multi-precision value, or in other words radix
34  * or base of numerical representation, e.g. base 2^64, base 2^32,
35  * base 2^26. Complementary characteristic is how wide is the result of
36  * multiplication of pair of digits, e.g. it would take 128 bits to
37  * accommodate multiplication result in base 2^64 case. These are used
38  * interchangeably. To describe implementation that is. But interface
39  * is designed to isolate this so that low-level primitives implemented
40  * in assembly can be self-contained/self-coherent.
41  */
42 #ifndef POLY1305_ASM
43 /*
44  * Even though there is __int128 reference implementation targeting
45  * 64-bit platforms provided below, it's not obvious that it's optimal
46  * choice for every one of them. Depending on instruction set overall
47  * amount of instructions can be comparable to one in __int64
48  * implementation. Amount of multiplication instructions would be lower,
49  * but not necessarily overall. And in out-of-order execution context,
50  * it is the latter that can be crucial...
51  *
52  * On related note. Poly1305 author, D. J. Bernstein, discusses and
53  * provides floating-point implementations of the algorithm in question.
54  * It made a lot of sense by the time of introduction, because most
55  * then-modern processors didn't have pipelined integer multiplier.
56  * [Not to mention that some had non-constant timing for integer
57  * multiplications.] Floating-point instructions on the other hand could
58  * be issued every cycle, which allowed to achieve better performance.
59  * Nowadays, with SIMD and/or out-or-order execution, shared or
60  * even emulated FPU, it's more complicated, and floating-point
61  * implementation is not necessarily optimal choice in every situation,
62  * rather contrary...
63  *
64  *                                              <appro@openssl.org>
65  */
66
67 typedef unsigned int u32;
68
69 /*
70  * poly1305_blocks processes a multiple of POLY1305_BLOCK_SIZE blocks
71  * of |inp| no longer than |len|. Behaviour for |len| not divisible by
72  * block size is unspecified in general case, even though in reference
73  * implementation the trailing chunk is simply ignored. Per algorithm
74  * specification, every input block, complete or last partial, is to be
75  * padded with a bit past most significant byte. The latter kind is then
76  * padded with zeros till block size. This last partial block padding
77  * is caller(*)'s responsibility, and because of this the last partial
78  * block is always processed with separate call with |len| set to
79  * POLY1305_BLOCK_SIZE and |padbit| to 0. In all other cases |padbit|
80  * should be set to 1 to perform implicit padding with 128th bit.
81  * poly1305_blocks does not actually check for this constraint though,
82  * it's caller(*)'s responsibility to comply.
83  *
84  * (*)  In the context "caller" is not application code, but higher
85  *      level Poly1305_* from this very module, so that quirks are
86  *      handled locally.
87  */
88 static void
89 poly1305_blocks(void *ctx, const unsigned char *inp, size_t len, u32 padbit);
90
91 /*
92  * Type-agnostic "rip-off" from constant_time_locl.h
93  */
94 # define CONSTANT_TIME_CARRY(a,b) ( \
95          (a ^ ((a ^ b) | ((a - b) ^ b))) >> (sizeof(a) * 8 - 1) \
96          )
97
98 # if !defined(PEDANTIC) && \
99      (defined(__SIZEOF_INT128__) && __SIZEOF_INT128__==16) && \
100      (defined(__SIZEOF_LONG__) && __SIZEOF_LONG__==8)
101
102 typedef unsigned long u64;
103 typedef unsigned __int128 u128;
104
105 typedef struct {
106     u64 h[3];
107     u64 r[2];
108 } poly1305_internal;
109
110 /* pick 32-bit unsigned integer in little endian order */
111 static u64 U8TOU64(const unsigned char *p)
112 {
113     return (((u64)(p[0] & 0xff)) |
114             ((u64)(p[1] & 0xff) << 8) |
115             ((u64)(p[2] & 0xff) << 16) |
116             ((u64)(p[3] & 0xff) << 24) |
117             ((u64)(p[4] & 0xff) << 32) |
118             ((u64)(p[5] & 0xff) << 40) |
119             ((u64)(p[6] & 0xff) << 48) |
120             ((u64)(p[7] & 0xff) << 56));
121 }
122
123 /* store a 32-bit unsigned integer in little endian */
124 static void U64TO8(unsigned char *p, u64 v)
125 {
126     p[0] = (unsigned char)((v) & 0xff);
127     p[1] = (unsigned char)((v >> 8) & 0xff);
128     p[2] = (unsigned char)((v >> 16) & 0xff);
129     p[3] = (unsigned char)((v >> 24) & 0xff);
130     p[4] = (unsigned char)((v >> 32) & 0xff);
131     p[5] = (unsigned char)((v >> 40) & 0xff);
132     p[6] = (unsigned char)((v >> 48) & 0xff);
133     p[7] = (unsigned char)((v >> 56) & 0xff);
134 }
135
136 static void poly1305_init(void *ctx, const unsigned char key[16])
137 {
138     poly1305_internal *st = (poly1305_internal *) ctx;
139
140     /* h = 0 */
141     st->h[0] = 0;
142     st->h[1] = 0;
143     st->h[2] = 0;
144
145     /* r &= 0xffffffc0ffffffc0ffffffc0fffffff */
146     st->r[0] = U8TOU64(&key[0]) & 0x0ffffffc0fffffff;
147     st->r[1] = U8TOU64(&key[8]) & 0x0ffffffc0ffffffc;
148 }
149
150 static void
151 poly1305_blocks(void *ctx, const unsigned char *inp, size_t len, u32 padbit)
152 {
153     poly1305_internal *st = (poly1305_internal *)ctx;
154     u64 r0, r1;
155     u64 s1;
156     u64 h0, h1, h2, c;
157     u128 d0, d1;
158
159     r0 = st->r[0];
160     r1 = st->r[1];
161
162     s1 = r1 + (r1 >> 2);
163
164     h0 = st->h[0];
165     h1 = st->h[1];
166     h2 = st->h[2];
167
168     while (len >= POLY1305_BLOCK_SIZE) {
169         /* h += m[i] */
170         h0 = (u64)(d0 = (u128)h0 + U8TOU64(inp + 0));
171         h1 = (u64)(d1 = (u128)h1 + (d0 >> 64) + U8TOU64(inp + 8));
172         /*
173          * padbit can be zero only when original len was
174          * POLY1306_BLOCK_SIZE, but we don't check
175          */
176         h2 += (u64)(d1 >> 64) + padbit;
177
178         /* h *= r "%" p, where "%" stands for "partial remainder" */
179         d0 = ((u128)h0 * r0) +
180              ((u128)h1 * s1);
181         d1 = ((u128)h0 * r1) +
182              ((u128)h1 * r0) +
183              (h2 * s1);
184         h2 = (h2 * r0);
185
186         /* last reduction step: */
187         /* a) h2:h0 = h2<<128 + d1<<64 + d0 */
188         h0 = (u64)d0;
189         h1 = (u64)(d1 += d0 >> 64);
190         h2 += (u64)(d1 >> 64);
191         /* b) (h2:h0 += (h2:h0>>130) * 5) %= 2^130 */
192         c = (h2 >> 2) + (h2 & ~3UL);
193         h2 &= 3;
194         h0 += c;
195         h1 += (c = CONSTANT_TIME_CARRY(h0,c));
196         h2 += CONSTANT_TIME_CARRY(h1,c);
197         /*
198          * Occasional overflows to 3rd bit of h2 are taken care of
199          * "naturally". If after this point we end up at the top of
200          * this loop, then the overflow bit will be accounted for
201          * in next iteration. If we end up in poly1305_emit, then
202          * comparison to modulus below will still count as "carry
203          * into 131st bit", so that properly reduced value will be
204          * picked in conditional move.
205          */
206
207         inp += POLY1305_BLOCK_SIZE;
208         len -= POLY1305_BLOCK_SIZE;
209     }
210
211     st->h[0] = h0;
212     st->h[1] = h1;
213     st->h[2] = h2;
214 }
215
216 static void poly1305_emit(void *ctx, unsigned char mac[16],
217                           const u32 nonce[4])
218 {
219     poly1305_internal *st = (poly1305_internal *) ctx;
220     u64 h0, h1, h2;
221     u64 g0, g1, g2;
222     u128 t;
223     u64 mask;
224
225     h0 = st->h[0];
226     h1 = st->h[1];
227     h2 = st->h[2];
228
229     /* compare to modulus by computing h + -p */
230     g0 = (u64)(t = (u128)h0 + 5);
231     g1 = (u64)(t = (u128)h1 + (t >> 64));
232     g2 = h2 + (u64)(t >> 64);
233
234     /* if there was carry into 131st bit, h1:h0 = g1:g0 */
235     mask = 0 - (g2 >> 2);
236     g0 &= mask;
237     g1 &= mask;
238     mask = ~mask;
239     h0 = (h0 & mask) | g0;
240     h1 = (h1 & mask) | g1;
241
242     /* mac = (h + nonce) % (2^128) */
243     h0 = (u64)(t = (u128)h0 + nonce[0] + ((u64)nonce[1]<<32));
244     h1 = (u64)(t = (u128)h1 + nonce[2] + ((u64)nonce[3]<<32) + (t >> 64));
245
246     U64TO8(mac + 0, h0);
247     U64TO8(mac + 8, h1);
248 }
249
250 # else
251
252 #  if defined(_WIN32) && !defined(__MINGW32__)
253 typedef unsigned __int64 u64;
254 #  elif defined(__arch64__)
255 typedef unsigned long u64;
256 #  else
257 typedef unsigned long long u64;
258 #  endif
259
260 typedef struct {
261     u32 h[5];
262     u32 r[4];
263 } poly1305_internal;
264
265 /* store a 32-bit unsigned integer in little endian */
266 static void U32TO8(unsigned char *p, unsigned int v)
267 {
268     p[0] = (unsigned char)((v) & 0xff);
269     p[1] = (unsigned char)((v >> 8) & 0xff);
270     p[2] = (unsigned char)((v >> 16) & 0xff);
271     p[3] = (unsigned char)((v >> 24) & 0xff);
272 }
273
274 static void poly1305_init(void *ctx, const unsigned char key[16])
275 {
276     poly1305_internal *st = (poly1305_internal *) ctx;
277
278     /* h = 0 */
279     st->h[0] = 0;
280     st->h[1] = 0;
281     st->h[2] = 0;
282     st->h[3] = 0;
283     st->h[4] = 0;
284
285     /* r &= 0xffffffc0ffffffc0ffffffc0fffffff */
286     st->r[0] = U8TOU32(&key[0]) & 0x0fffffff;
287     st->r[1] = U8TOU32(&key[4]) & 0x0ffffffc;
288     st->r[2] = U8TOU32(&key[8]) & 0x0ffffffc;
289     st->r[3] = U8TOU32(&key[12]) & 0x0ffffffc;
290 }
291
292 static void
293 poly1305_blocks(void *ctx, const unsigned char *inp, size_t len, u32 padbit)
294 {
295     poly1305_internal *st = (poly1305_internal *)ctx;
296     u32 r0, r1, r2, r3;
297     u32 s1, s2, s3;
298     u32 h0, h1, h2, h3, h4, c;
299     u64 d0, d1, d2, d3;
300
301     r0 = st->r[0];
302     r1 = st->r[1];
303     r2 = st->r[2];
304     r3 = st->r[3];
305
306     s1 = r1 + (r1 >> 2);
307     s2 = r2 + (r2 >> 2);
308     s3 = r3 + (r3 >> 2);
309
310     h0 = st->h[0];
311     h1 = st->h[1];
312     h2 = st->h[2];
313     h3 = st->h[3];
314     h4 = st->h[4];
315
316     while (len >= POLY1305_BLOCK_SIZE) {
317         /* h += m[i] */
318         h0 = (u32)(d0 = (u64)h0 + U8TOU32(inp + 0));
319         h1 = (u32)(d1 = (u64)h1 + (d0 >> 32) + U8TOU32(inp + 4));
320         h2 = (u32)(d2 = (u64)h2 + (d1 >> 32) + U8TOU32(inp + 8));
321         h3 = (u32)(d3 = (u64)h3 + (d2 >> 32) + U8TOU32(inp + 12));
322         h4 += (u32)(d3 >> 32) + padbit;
323
324         /* h *= r "%" p, where "%" stands for "partial remainder" */
325         d0 = ((u64)h0 * r0) +
326              ((u64)h1 * s3) +
327              ((u64)h2 * s2) +
328              ((u64)h3 * s1);
329         d1 = ((u64)h0 * r1) +
330              ((u64)h1 * r0) +
331              ((u64)h2 * s3) +
332              ((u64)h3 * s2) +
333              (h4 * s1);
334         d2 = ((u64)h0 * r2) +
335              ((u64)h1 * r1) +
336              ((u64)h2 * r0) +
337              ((u64)h3 * s3) +
338              (h4 * s2);
339         d3 = ((u64)h0 * r3) +
340              ((u64)h1 * r2) +
341              ((u64)h2 * r1) +
342              ((u64)h3 * r0) +
343              (h4 * s3);
344         h4 = (h4 * r0);
345
346         /* last reduction step: */
347         /* a) h4:h0 = h4<<128 + d3<<96 + d2<<64 + d1<<32 + d0 */
348         h0 = (u32)d0;
349         h1 = (u32)(d1 += d0 >> 32);
350         h2 = (u32)(d2 += d1 >> 32);
351         h3 = (u32)(d3 += d2 >> 32);
352         h4 += (u32)(d3 >> 32);
353         /* b) (h4:h0 += (h4:h0>>130) * 5) %= 2^130 */
354         c = (h4 >> 2) + (h4 & ~3U);
355         h4 &= 3;
356         h0 += c;
357         h1 += (c = CONSTANT_TIME_CARRY(h0,c));
358         h2 += (c = CONSTANT_TIME_CARRY(h1,c));
359         h3 += (c = CONSTANT_TIME_CARRY(h2,c));
360         h4 += CONSTANT_TIME_CARRY(h3,c);
361         /*
362          * Occasional overflows to 3rd bit of h4 are taken care of
363          * "naturally". If after this point we end up at the top of
364          * this loop, then the overflow bit will be accounted for
365          * in next iteration. If we end up in poly1305_emit, then
366          * comparison to modulus below will still count as "carry
367          * into 131st bit", so that properly reduced value will be
368          * picked in conditional move.
369          */
370
371         inp += POLY1305_BLOCK_SIZE;
372         len -= POLY1305_BLOCK_SIZE;
373     }
374
375     st->h[0] = h0;
376     st->h[1] = h1;
377     st->h[2] = h2;
378     st->h[3] = h3;
379     st->h[4] = h4;
380 }
381
382 static void poly1305_emit(void *ctx, unsigned char mac[16],
383                           const u32 nonce[4])
384 {
385     poly1305_internal *st = (poly1305_internal *) ctx;
386     u32 h0, h1, h2, h3, h4;
387     u32 g0, g1, g2, g3, g4;
388     u64 t;
389     u32 mask;
390
391     h0 = st->h[0];
392     h1 = st->h[1];
393     h2 = st->h[2];
394     h3 = st->h[3];
395     h4 = st->h[4];
396
397     /* compare to modulus by computing h + -p */
398     g0 = (u32)(t = (u64)h0 + 5);
399     g1 = (u32)(t = (u64)h1 + (t >> 32));
400     g2 = (u32)(t = (u64)h2 + (t >> 32));
401     g3 = (u32)(t = (u64)h3 + (t >> 32));
402     g4 = h4 + (u32)(t >> 32);
403
404     /* if there was carry into 131st bit, h3:h0 = g3:g0 */
405     mask = 0 - (g4 >> 2);
406     g0 &= mask;
407     g1 &= mask;
408     g2 &= mask;
409     g3 &= mask;
410     mask = ~mask;
411     h0 = (h0 & mask) | g0;
412     h1 = (h1 & mask) | g1;
413     h2 = (h2 & mask) | g2;
414     h3 = (h3 & mask) | g3;
415
416     /* mac = (h + nonce) % (2^128) */
417     h0 = (u32)(t = (u64)h0 + nonce[0]);
418     h1 = (u32)(t = (u64)h1 + (t >> 32) + nonce[1]);
419     h2 = (u32)(t = (u64)h2 + (t >> 32) + nonce[2]);
420     h3 = (u32)(t = (u64)h3 + (t >> 32) + nonce[3]);
421
422     U32TO8(mac + 0, h0);
423     U32TO8(mac + 4, h1);
424     U32TO8(mac + 8, h2);
425     U32TO8(mac + 12, h3);
426 }
427 # endif
428 #else
429 int poly1305_init(void *ctx, const unsigned char key[16], void *func);
430 void poly1305_blocks(void *ctx, const unsigned char *inp, size_t len,
431                      unsigned int padbit);
432 void poly1305_emit(void *ctx, unsigned char mac[16],
433                    const unsigned int nonce[4]);
434 #endif
435
436 void Poly1305_Init(POLY1305 *ctx, const unsigned char key[32])
437 {
438     ctx->nonce[0] = U8TOU32(&key[16]);
439     ctx->nonce[1] = U8TOU32(&key[20]);
440     ctx->nonce[2] = U8TOU32(&key[24]);
441     ctx->nonce[3] = U8TOU32(&key[28]);
442
443 #ifndef POLY1305_ASM
444     poly1305_init(ctx->opaque, key);
445 #else
446     /*
447      * Unlike reference poly1305_init assembly counterpart is expected
448      * to return a value: non-zero if it initializes ctx->func, and zero
449      * otherwise. Latter is to simplify assembly in cases when there no
450      * multiple code paths to switch between.
451      */
452     if (!poly1305_init(ctx->opaque, key, &ctx->func)) {
453         ctx->func.blocks = poly1305_blocks;
454         ctx->func.emit = poly1305_emit;
455     }
456 #endif
457
458     ctx->num = 0;
459
460 }
461
462 #ifdef POLY1305_ASM
463 /*
464  * This "eclipses" poly1305_blocks and poly1305_emit, but it's
465  * conscious choice imposed by -Wshadow compiler warnings.
466  */
467 # define poly1305_blocks (*poly1305_blocks_p)
468 # define poly1305_emit   (*poly1305_emit_p)
469 #endif
470
471 void Poly1305_Update(POLY1305 *ctx, const unsigned char *inp, size_t len)
472 {
473 #ifdef POLY1305_ASM
474     /*
475      * As documented, poly1305_blocks is never called with input
476      * longer than single block and padbit argument set to 0. This
477      * property is fluently used in assembly modules to optimize
478      * padbit handling on loop boundary.
479      */
480     poly1305_blocks_f poly1305_blocks_p = ctx->func.blocks;
481 #endif
482     size_t rem, num;
483
484     if ((num = ctx->num)) {
485         rem = POLY1305_BLOCK_SIZE - num;
486         if (len >= rem) {
487             memcpy(ctx->data + num, inp, rem);
488             poly1305_blocks(ctx->opaque, ctx->data, POLY1305_BLOCK_SIZE, 1);
489             inp += rem;
490             len -= rem;
491         } else {
492             /* Still not enough data to process a block. */
493             memcpy(ctx->data + num, inp, len);
494             ctx->num = num + len;
495             return;
496         }
497     }
498
499     rem = len % POLY1305_BLOCK_SIZE;
500     len -= rem;
501
502     if (len >= POLY1305_BLOCK_SIZE) {
503         poly1305_blocks(ctx->opaque, inp, len, 1);
504         inp += len;
505     }
506
507     if (rem)
508         memcpy(ctx->data, inp, rem);
509
510     ctx->num = rem;
511 }
512
513 void Poly1305_Final(POLY1305 *ctx, unsigned char mac[16])
514 {
515 #ifdef POLY1305_ASM
516     poly1305_blocks_f poly1305_blocks_p = ctx->func.blocks;
517     poly1305_emit_f poly1305_emit_p = ctx->func.emit;
518 #endif
519     size_t num;
520
521     if ((num = ctx->num)) {
522         ctx->data[num++] = 1;   /* pad bit */
523         while (num < POLY1305_BLOCK_SIZE)
524             ctx->data[num++] = 0;
525         poly1305_blocks(ctx->opaque, ctx->data, POLY1305_BLOCK_SIZE, 0);
526     }
527
528     poly1305_emit(ctx->opaque, mac, ctx->nonce);
529
530     /* zero out the state */
531     OPENSSL_cleanse(ctx, sizeof(*ctx));
532 }