From ee39a683a2b623a1da0747ec20f20b63470a2db6 Mon Sep 17 00:00:00 2001 From: Mark Wooding Date: Tue, 18 Apr 2017 00:39:24 +0100 Subject: [PATCH] math/f25519.c: Implementation for arithmetic in GF(2^255 - 19). There's both a fast implementation for platforms with 64-bit arithmetic, and a slow baseline for minimal C89 platforms. The code works better on two's complement systems with arithmetic right shifts, but it works portably. * Arithmetic shifts are implemented with hairy masking and exact division, but GCC notices and optimizes accordingly. * Two's complement is used in the conditional-swap machinery, but there's a fallback using multiplication if the `configure' script can't detect it. --- configure.ac | 15 + math/Makefile.am | 18 + math/f25519.c | 1117 +++++++++++++++++++++++++++++++++++++++++++++++++ math/f25519.h | 209 +++++++++ math/qfarith.h | 229 ++++++++++ math/t/f25519 | 594 ++++++++++++++++++++++++++ utils/curve25519.sage | 43 ++ utils/qfarith-test | 109 +++++ 8 files changed, 2334 insertions(+) create mode 100644 math/f25519.c create mode 100644 math/f25519.h create mode 100644 math/qfarith.h create mode 100644 math/t/f25519 create mode 100644 utils/curve25519.sage create mode 100755 utils/qfarith-test diff --git a/configure.ac b/configure.ac index 6171efd2..bd2c3b5c 100644 --- a/configure.ac +++ b/configure.ac @@ -307,6 +307,21 @@ fi catacomb_LIMIT([SIZET], [=0], [~(size_t)0]) AC_SUBST([limits]) +dnl Figure out other aspects of the implementation's arithmetic. +AC_CACHE_CHECK([whether negative numbers use two's complement], + [catacomb_cv_neg_twoc], +[AC_TRY_COMPILE( +[#include ], +[int check[2*!!(-INT_MAX == INT_MIN + 1) - 1];], +[catacomb_cv_neg_twoc=yes], +[catacomb_cv_neg_twoc=no])]) +case $catacomb_cv_neg_twoc in + yes) + AC_DEFINE([NEG_TWOC], [1], + [Define if signed numbers are represented in two's complement.]) + ;; +esac + dnl Functions used for noise-gathering. AC_CHECK_FUNCS([setgroups]) AC_CHECK_HEADERS([linux/random.h]) diff --git a/math/Makefile.am b/math/Makefile.am index e8b11195..2a60ee0f 100644 --- a/math/Makefile.am +++ b/math/Makefile.am @@ -413,4 +413,22 @@ ectab.c: $(mpgen) typeinfo.py ectab.in $(MPGEN) ectab $(srcdir)/ectab.in >ectab.c.new && \ mv ectab.c.new ectab.c +###-------------------------------------------------------------------------- +### Other strange things. + +pkginclude_HEADERS += qfarith.h + +pkginclude_HEADERS += f25519.h +libmath_la_SOURCES += f25519.c +TESTS += f25519.t$(EXEEXT) +TESTS += f25519-p10.t$(EXEEXT) +EXTRA_DIST += t/f25519 + +check_PROGRAMS += f25519-p10.t +f25519_p10_t_SOURCES = f25519.c +f25519_p10_t_CPPFLAGS = $(AM_CPPFLAGS) -DTEST_RIG -DSRCDIR="\"$(srcdir)\"" +f25519_p10_t_CPPFLAGS += -DF25519_IMPL=10 +f25519_p10_t_LDADD = $(TEST_LIBS) $(top_builddir)/libcatacomb.la +f25519_p10_t_LDADD += $(mLib_LIBS) $(CATACOMB_LIBS) $(LIBS) + ###----- That's all, folks -------------------------------------------------- diff --git a/math/f25519.c b/math/f25519.c new file mode 100644 index 00000000..220a73c8 --- /dev/null +++ b/math/f25519.c @@ -0,0 +1,1117 @@ +/* -*-c-*- + * + * Arithmetic modulo 2^255 - 19 + * + * (c) 2017 Straylight/Edgeware + */ + +/*----- Licensing notice --------------------------------------------------* + * + * This file is part of Catacomb. + * + * Catacomb is free software; you can redistribute it and/or modify + * it under the terms of the GNU Library General Public License as + * published by the Free Software Foundation; either version 2 of the + * License, or (at your option) any later version. + * + * Catacomb is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU Library General Public License for more details. + * + * You should have received a copy of the GNU Library General Public + * License along with Catacomb; if not, write to the Free + * Software Foundation, Inc., 59 Temple Place - Suite 330, Boston, + * MA 02111-1307, USA. + */ + +/*----- Header files ------------------------------------------------------*/ + +#include "config.h" + +#include "f25519.h" + +/*----- Basic setup -------------------------------------------------------*/ + +#if F25519_IMPL == 26 +/* Elements x of GF(2^255 - 19) are represented by ten signed integers x_i: x + * = SUM_{0<=i<10} x_i 2^ceil(51i/2), mostly following Bernstein's original + * paper. + */ + +typedef int32 piece; typedef int64 dblpiece; +typedef uint32 upiece; typedef uint64 udblpiece; +#define P p26 +#define PIECEWD(i) ((i)%2 ? 25 : 26) +#define NPIECE 10 + +#define M26 0x03ffffffu +#define M25 0x01ffffffu +#define B26 0x04000000u +#define B25 0x02000000u +#define B24 0x01000000u + +#define PIECES(v) v##0, v##1, v##2, v##3, v##4, v##5, v##6, v##7, v##8, v##9 +#define FETCH(v, w) do { \ + v##0 = (w)->P[0]; v##1 = (w)->P[1]; \ + v##2 = (w)->P[2]; v##3 = (w)->P[3]; \ + v##4 = (w)->P[4]; v##5 = (w)->P[5]; \ + v##6 = (w)->P[6]; v##7 = (w)->P[7]; \ + v##8 = (w)->P[8]; v##9 = (w)->P[9]; \ +} while (0) +#define STASH(w, v) do { \ + (w)->P[0] = v##0; (w)->P[1] = v##1; \ + (w)->P[2] = v##2; (w)->P[3] = v##3; \ + (w)->P[4] = v##4; (w)->P[5] = v##5; \ + (w)->P[6] = v##6; (w)->P[7] = v##7; \ + (w)->P[8] = v##8; (w)->P[9] = v##9; \ +} while (0) + +#elif F25519_IMPL == 10 +/* Elements x of GF(2^255 - 19) are represented by 26 signed integers x_i: x + * = SUM_{0<=i<26} x_i 2^ceil(255i/26); i.e., most pieces are 10 bits wide, + * except for pieces 5, 10, 15, 20, and 25 which have 9 bits. + */ + +typedef int16 piece; typedef int32 dblpiece; +typedef uint16 upiece; typedef uint32 udblpiece; +#define P p10 +#define PIECEWD(i) \ + ((i) == 5 || (i) == 10 || (i) == 15 || (i) == 20 || (i) == 25 ? 9 : 10) +#define NPIECE 26 + +#define B10 0x0400 +#define B9 0x200 +#define B8 0x100 +#define M10 0x3ff +#define M9 0x1ff + +#endif + +/*----- Debugging machinery -----------------------------------------------*/ + +#if defined(F25519_DEBUG) || defined(TEST_RIG) + +#include + +#include "mp.h" +#include "mptext.h" + +static mp *get_2p255m91(void) +{ + mpw w19 = 19; + mp *p = MP_NEW, m19; + + p = mp_setbit(p, MP_ZERO, 255); + mp_build(&m19, &w19, &w19 + 1); + p = mp_sub(p, p, &m19); + return (p); +} + +DEF_FDUMP(fdump, piece, PIECEWD, NPIECE, 32, get_2p255m91()) + +#endif + +/*----- Loading and storing -----------------------------------------------*/ + +/* --- @f25519_load@ --- * + * + * Arguments: @f25519 *z@ = where to store the result + * @const octet xv[32]@ = source to read + * + * Returns: --- + * + * Use: Reads an element of %$\gf{2^{255} - 19}$% in external + * representation from @xv@ and stores it in @z@. + * + * External representation is little-endian base-256. Some + * elements have multiple encodings, which are not produced by + * correct software; use of noncanonical encodings is not an + * error, and toleration of them is considered a performance + * feature. + */ + +void f25519_load(f25519 *z, const octet xv[32]) +{ +#if F25519_IMPL == 26 + + uint32 xw0 = LOAD32_L(xv + 0), xw1 = LOAD32_L(xv + 4), + xw2 = LOAD32_L(xv + 8), xw3 = LOAD32_L(xv + 12), + xw4 = LOAD32_L(xv + 16), xw5 = LOAD32_L(xv + 20), + xw6 = LOAD32_L(xv + 24), xw7 = LOAD32_L(xv + 28); + piece PIECES(x), b, c; + + /* First, split the 32-bit words into the irregularly-sized pieces we need + * for the field representation. These pieces are all positive. We'll do + * the sign correction afterwards. + * + * It may be that the top bit of the input is set. Avoid trouble by + * folding that back round into the bottom piece of the representation. + * + * Here, we briefly have 0 <= x_0 < 2^26 + 19, but will resolve this later. + * Otherwise, we have 0 <= x_{2i} < 2^26, and 0 <= x_{2i+1} < 2^25. + */ + x0 = ((xw0 << 0)&0x03ffffff) + 19*((xw7 >> 31)&0x00000001); + x1 = ((xw1 << 6)&0x01ffffc0) | ((xw0 >> 26)&0x0000003f); + x2 = ((xw2 << 13)&0x03ffe000) | ((xw1 >> 19)&0x00001fff); + x3 = ((xw3 << 19)&0x01f80000) | ((xw2 >> 13)&0x0007ffff); + x4 = ((xw3 >> 6)&0x03ffffff); + x5 = (xw4 << 0)&0x01ffffff; + x6 = ((xw5 << 7)&0x03ffff80) | ((xw4 >> 25)&0x0000007f); + x7 = ((xw6 << 13)&0x01ffe000) | ((xw5 >> 19)&0x00001fff); + x8 = ((xw7 << 20)&0x03f00000) | ((xw6 >> 12)&0x000fffff); + x9 = ((xw7 >> 6)&0x01ffffff); + + /* Next, we convert these pieces into a roughly balanced signed + * representation. For each piece, check to see if its top bit is set. If + * it is, then lend a bit to the next piece over. For x_9, this needs to + * be carried around, which is a little fiddly. In particular, we delay + * the carry until after all of the pieces have been balanced. If we don't + * do this, then we have to do a more expensive test for nonzeroness to + * decide whether to lend a bit leftwards rather than just testing a single + * bit. + * + * This fixes up the anomalous size of x_0: the loan of a bit becomes an + * actual carry if x_0 >= 2^26. By the end, then, we have: + * + * { 2^25 if i even + * |x_i| <= { + * { 2^24 if i odd + * + * Note that we don't try for a canonical representation here: both upper + * and lower bounds are achievable. + * + * All of the x_i at this point are positive, so we don't need to do + * anything wierd when masking them. + */ + b = x9&B24; c = 19&((b >> 19) - (b >> 24)); x9 -= b << 1; + b = x8&B25; x9 += b >> 25; x8 -= b << 1; + b = x7&B24; x8 += b >> 24; x7 -= b << 1; + b = x6&B25; x7 += b >> 25; x6 -= b << 1; + b = x5&B24; x6 += b >> 24; x5 -= b << 1; + b = x4&B25; x5 += b >> 25; x4 -= b << 1; + b = x3&B24; x4 += b >> 24; x3 -= b << 1; + b = x2&B25; x3 += b >> 25; x2 -= b << 1; + b = x1&B24; x2 += b >> 24; x1 -= b << 1; + b = x0&B25; x1 += (b >> 25) + (x0 >> 26); x0 = (x0&M26) - (b << 1); + x0 += c; + + /* And with that, we're done. */ + STASH(z, x); + +#elif F25519_IMPL == 10 + + piece x[NPIECE]; + unsigned i, j, n, wd; + uint32 a; + int b, c; + + /* First, just get the content out of the buffer. */ + for (i = j = a = n = 0, wd = 10; j < NPIECE; i++) { + a |= (uint32)xv[i] << n; n += 8; + if (n >= wd) { + x[j++] = a&MASK(wd); + a >>= wd; n -= wd; + wd = PIECEWD(j); + } + } + + /* There's a little bit left over from the top byte. Carry it into the low + * piece. + */ + x[0] += 19*(int)(a&MASK(n)); + + /* Next, convert the pieces into a roughly balanced signed representation. + * If a piece's top bit is set, lend a bit to the next piece over. For + * x_25, this needs to be carried around, which is a bit fiddly. + */ + b = x[NPIECE - 1]&B8; + c = 19&((b >> 3) - (b >> 8)); + x[NPIECE - 1] -= b << 1; + for (i = NPIECE - 2; i > 0; i--) { + wd = PIECEWD(i) - 1; + b = x[i]&BIT(wd); + x[i + 1] += b >> wd; + x[i] -= b << 1; + } + b = x[0]&B9; + x[1] += (b >> 9) + (x[0] >> 10); + x[0] = (x[0]&M10) - (b << 1) + c; + + /* And we're done. */ + for (i = 0; i < NPIECE; i++) z->P[i] = x[i]; + +#endif +} + +/* --- @f25519_store@ --- * + * + * Arguments: @octet zv[32]@ = where to write the result + * @const f25519 *x@ = the field element to write + * + * Returns: --- + * + * Use: Stores a field element in the given octet vector in external + * representation. A canonical encoding is always stored, so, + * in particular, the top bit of @xv[31]@ is always left clear. + */ + +void f25519_store(octet zv[32], const f25519 *x) +{ +#if F25519_IMPL == 26 + + piece PIECES(x), PIECES(y), c, d; + uint32 zw0, zw1, zw2, zw3, zw4, zw5, zw6, zw7; + mask32 m; + + FETCH(x, x); + + /* First, propagate the carries throughout the pieces. By the end of this, + * we'll have all of the pieces canonically sized and positive, and maybe + * there'll be (signed) carry out. The carry c is in { -1, 0, +1 }, and + * the remaining value will be in the half-open interval [0, 2^255). The + * whole represented value is then x + 2^255 c. + * + * It's worth paying careful attention to the bounds. We assume that we + * start out with |x_i| <= 2^30. We start by cutting off and reducing the + * carry c_9 from the topmost piece, x_9. This leaves 0 <= x_9 < 2^25; and + * we'll have |c_9| <= 2^5. We multiply this by 19 and we'll add it onto + * x_0 and propagate the carries: but what bounds can we calculate on x + * before this? + * + * Let o_i = floor(51 i/2). We have X_i = SUM_{0<=j= 2^255 - 19, then we should subtract 2^255 - 19 from the whole + * value; if c = -1 then we should add 2^255 - 19; and otherwise we should + * do nothing. + * + * But conditional behaviour is bad, m'kay. So here's what we do instead. + * + * The first job is to sort out what we wanted to do. If c = -1 then we + * want to (a) invert the constant addend and (b) feed in a carry-in; + * otherwise, we don't. + */ + m = SIGN(c); + d = m&1; + + /* Now do the addition/subtraction. Remember that all of the x_i are + * nonnegative, so shifting and masking are safe and easy. + */ + d += x0 + (19 ^ (M26&m)); y0 = d&M26; d >>= 26; + d += x1 + (M25&m); y1 = d&M25; d >>= 25; + d += x2 + (M26&m); y2 = d&M26; d >>= 26; + d += x3 + (M25&m); y3 = d&M25; d >>= 25; + d += x4 + (M26&m); y4 = d&M26; d >>= 26; + d += x5 + (M25&m); y5 = d&M25; d >>= 25; + d += x6 + (M26&m); y6 = d&M26; d >>= 26; + d += x7 + (M25&m); y7 = d&M25; d >>= 25; + d += x8 + (M26&m); y8 = d&M26; d >>= 26; + d += x9 + (M25&m); y9 = d&M25; d >>= 25; + + /* The final carry-out is in d; since we only did addition, and the x_i are + * nonnegative, then d is in { 0, 1 }. We want to keep y, rather than x, + * if (a) c /= 0 (in which case we know that the old value was + * unsatisfactory), or (b) if d = 1 (in which case, if c = 0, we know that + * the subtraction didn't cause a borrow, so we must be in the case where + * 2^255 - 19 <= x < 2^255). + */ + m = NONZEROP(c) | ~NONZEROP(d - 1); + x0 = (y0&m) | (x0&~m); x1 = (y1&m) | (x1&~m); + x2 = (y2&m) | (x2&~m); x3 = (y3&m) | (x3&~m); + x4 = (y4&m) | (x4&~m); x5 = (y5&m) | (x5&~m); + x6 = (y6&m) | (x6&~m); x7 = (y7&m) | (x7&~m); + x8 = (y8&m) | (x8&~m); x9 = (y9&m) | (x9&~m); + + /* Extract 32-bit words from the value. */ + zw0 = ((x0 >> 0)&0x03ffffff) | (((uint32)x1 << 26)&0xfc000000); + zw1 = ((x1 >> 6)&0x0007ffff) | (((uint32)x2 << 19)&0xfff80000); + zw2 = ((x2 >> 13)&0x00001fff) | (((uint32)x3 << 13)&0xffffe000); + zw3 = ((x3 >> 19)&0x0000003f) | (((uint32)x4 << 6)&0xffffffc0); + zw4 = ((x5 >> 0)&0x01ffffff) | (((uint32)x6 << 25)&0xfe000000); + zw5 = ((x6 >> 7)&0x0007ffff) | (((uint32)x7 << 19)&0xfff80000); + zw6 = ((x7 >> 13)&0x00000fff) | (((uint32)x8 << 12)&0xfffff000); + zw7 = ((x8 >> 20)&0x0000003f) | (((uint32)x9 << 6)&0x7fffffc0); + + /* Store the result as an octet string. */ + STORE32_L(zv + 0, zw0); STORE32_L(zv + 4, zw1); + STORE32_L(zv + 8, zw2); STORE32_L(zv + 12, zw3); + STORE32_L(zv + 16, zw4); STORE32_L(zv + 20, zw5); + STORE32_L(zv + 24, zw6); STORE32_L(zv + 28, zw7); + +#elif F25519_IMPL == 10 + + piece y[NPIECE], yy[NPIECE], c, d; + unsigned i, j, n, wd; + uint32 m, a; + + /* Before we do anything, copy the input so we can hack on it. */ + for (i = 0; i < NPIECE; i++) y[i] = x->P[i]; + + /* First, propagate the carries throughout the pieces. + * + * It's worth paying careful attention to the bounds. We assume that we + * start out with |y_i| <= 2^14. We start by cutting off and reducing the + * carry c_25 from the topmost piece, y_25. This leaves 0 <= y_25 < 2^9; + * and we'll have |c_25| <= 2^5. We multiply this by 19 and we'll ad it + * onto y_0 and propagte the carries: but what bounds can we calculate on + * y before this? + * + * Let o_i = floor(255 i/26). We have Y_i = SUM_{0<=j>= 10; + for (i = 1; i < NPIECE; i++) { + wd = PIECEWD(i); + d += y[i] + (MASK(wd)&m); + yy[i] = d&MASK(wd); + d >>= wd; + } + + /* Choose which value to keep. */ + m = NONZEROP(c) | ~NONZEROP(d - 1); + for (i = 0; i < NPIECE; i++) y[i] = (yy[i]&m) | (y[i]&~m); + + /* Store the result as an octet string. */ + for (i = j = a = n = 0; i < NPIECE; i++) { + a |= (upiece)y[i] << n; n += PIECEWD(i); + while (n >= 8) { + zv[j++] = a&0xff; + a >>= 8; n -= 8; + } + } + zv[j++] = a; + +#endif +} + +/* --- @f25519_set@ --- * + * + * Arguments: @f25519 *z@ = where to write the result + * @int a@ = a small-ish constant + * + * Returns: --- + * + * Use: Sets @z@ to equal @a@. + */ + +void f25519_set(f25519 *x, int a) +{ + unsigned i; + + x->P[0] = a; + for (i = 1; i < NPIECE; i++) x->P[i] = 0; +} + +/*----- Basic arithmetic --------------------------------------------------*/ + +/* --- @f25519_add@ --- * + * + * Arguments: @f25519 *z@ = where to put the result (may alias @x@ or @y@) + * @const f25519 *x, *y@ = two operands + * + * Returns: --- + * + * Use: Set @z@ to the sum %$x + y$%. + */ + +void f25519_add(f25519 *z, const f25519 *x, const f25519 *y) +{ +#if F25519_IMPL == 26 + z->P[0] = x->P[0] + y->P[0]; z->P[1] = x->P[1] + y->P[1]; + z->P[2] = x->P[2] + y->P[2]; z->P[3] = x->P[3] + y->P[3]; + z->P[4] = x->P[4] + y->P[4]; z->P[5] = x->P[5] + y->P[5]; + z->P[6] = x->P[6] + y->P[6]; z->P[7] = x->P[7] + y->P[7]; + z->P[8] = x->P[8] + y->P[8]; z->P[9] = x->P[9] + y->P[9]; +#elif F25519_IMPL == 10 + unsigned i; + for (i = 0; i < NPIECE; i++) z->P[i] = x->P[i] + y->P[i]; +#endif +} + +/* --- @f25519_sub@ --- * + * + * Arguments: @f25519 *z@ = where to put the result (may alias @x@ or @y@) + * @const f25519 *x, *y@ = two operands + * + * Returns: --- + * + * Use: Set @z@ to the difference %$x - y$%. + */ + +void f25519_sub(f25519 *z, const f25519 *x, const f25519 *y) +{ +#if F25519_IMPL == 26 + z->P[0] = x->P[0] - y->P[0]; z->P[1] = x->P[1] - y->P[1]; + z->P[2] = x->P[2] - y->P[2]; z->P[3] = x->P[3] - y->P[3]; + z->P[4] = x->P[4] - y->P[4]; z->P[5] = x->P[5] - y->P[5]; + z->P[6] = x->P[6] - y->P[6]; z->P[7] = x->P[7] - y->P[7]; + z->P[8] = x->P[8] - y->P[8]; z->P[9] = x->P[9] - y->P[9]; +#elif F25519_IMPL == 10 + unsigned i; + for (i = 0; i < NPIECE; i++) z->P[i] = x->P[i] - y->P[i]; +#endif +} + +/*----- Constant-time utilities -------------------------------------------*/ + +/* --- @f25519_condswap@ --- * + * + * Arguments: @f25519 *x, *y@ = two operands + * @uint32 m@ = a mask + * + * Returns: --- + * + * Use: If @m@ is zero, do nothing; if @m@ is all-bits-set, then + * exchange @x@ and @y@. If @m@ has some other value, then + * scramble @x@ and @y@ in an unhelpful way. + */ + +void f25519_condswap(f25519 *x, f25519 *y, uint32 m) +{ + mask32 mm = FIX_MASK32(m); + +#if F25519_IMPL == 26 + CONDSWAP(x->P[0], y->P[0], mm); + CONDSWAP(x->P[1], y->P[1], mm); + CONDSWAP(x->P[2], y->P[2], mm); + CONDSWAP(x->P[3], y->P[3], mm); + CONDSWAP(x->P[4], y->P[4], mm); + CONDSWAP(x->P[5], y->P[5], mm); + CONDSWAP(x->P[6], y->P[6], mm); + CONDSWAP(x->P[7], y->P[7], mm); + CONDSWAP(x->P[8], y->P[8], mm); + CONDSWAP(x->P[9], y->P[9], mm); +#elif F25519_IMPL == 10 + unsigned i; + for (i = 0; i < NPIECE; i++) CONDSWAP(x->P[i], y->P[i], mm); +#endif +} + +/*----- Multiplication ----------------------------------------------------*/ + +#if F25519_IMPL == 26 + +/* Let B = 2^63 - 1 be the largest value such that +B and -B can be + * represented in a double-precision piece. On entry, it must be the case + * that |X_i| <= M <= B - 2^25 for some M. If this is the case, then, on + * exit, we will have |Z_i| <= 2^25 + 19 M/2^25. + */ +#define CARRYSTEP(z, x, m, b, f, xx, n) do { \ + (z) = (dblpiece)((udblpiece)(x)&(m)) - (b) + \ + (f)*ASR(dblpiece, (xx), (n)); \ +} while (0) +#define CARRY_REDUCE(z, x) do { \ + dblpiece PIECES(_t); \ + \ + /* Bias the input pieces. This keeps the carries and so on centred \ + * around zero rather than biased positive. \ + */ \ + _t0 = (x##0) + B25; _t1 = (x##1) + B24; \ + _t2 = (x##2) + B25; _t3 = (x##3) + B24; \ + _t4 = (x##4) + B25; _t5 = (x##5) + B24; \ + _t6 = (x##6) + B25; _t7 = (x##7) + B24; \ + _t8 = (x##8) + B25; _t9 = (x##9) + B24; \ + \ + /* Calculate the reduced pieces. Careful with the bithacking. */ \ + CARRYSTEP(z##0, _t0, M26, B25, 19, _t9, 25); \ + CARRYSTEP(z##1, _t1, M25, B24, 1, _t0, 26); \ + CARRYSTEP(z##2, _t2, M26, B25, 1, _t1, 25); \ + CARRYSTEP(z##3, _t3, M25, B24, 1, _t2, 26); \ + CARRYSTEP(z##4, _t4, M26, B25, 1, _t3, 25); \ + CARRYSTEP(z##5, _t5, M25, B24, 1, _t4, 26); \ + CARRYSTEP(z##6, _t6, M26, B25, 1, _t5, 25); \ + CARRYSTEP(z##7, _t7, M25, B24, 1, _t6, 26); \ + CARRYSTEP(z##8, _t8, M26, B25, 1, _t7, 25); \ + CARRYSTEP(z##9, _t9, M25, B24, 1, _t8, 26); \ +} while (0) + +#elif F25519_IMPL == 10 + +/* Perform carry propagation on X. */ +static void carry_reduce(dblpiece x[NPIECE]) +{ + /* Initial bounds: we assume |x_i| < 2^31 - 2^27. */ + + unsigned i, j; + dblpiece c; + + /* The result is nearly canonical, because we do sequential carry + * propagation, because smaller processors are more likely to prefer the + * smaller working set than the instruction-level parallelism. + * + * Start at x_23; truncate it to 10 bits, and propagate the carry to x_24. + * Truncate x_24 to 10 bits, and add the carry onto x_25. Truncate x_25 to + * 9 bits, and add 19 times the carry onto x_0. And so on. + * + * Let c_i be the portion of x_i to be carried onto x_{i+1}. I claim that + * |c_i| <= 2^22. Then the carry /into/ any x_i has magnitude at most + * 19*2^22 < 2^27 (allowing for the reduction as we carry from x_25 to + * x_0), and x_i after carry is bounded above by 2^31. Hence, the carry + * out is at most 2^22, as claimed. + * + * Once we reach x_23 for the second time, we start with |x_23| <= 2^9. + * The carry into x_23 is at most 2^27 as calculated above; so the carry + * out into x_24 has magnitude at most 2^17. In turn, |x_24| <= 2^9 before + * the carry, so is now no more than 2^18 in magnitude, and the carry out + * into x_25 is at most 2^8. This leaves |x_25| < 2^9 after carry + * propagation. + * + * Be careful with the bit hacking because the quantities involved are + * signed. + */ + + /*For each piece, we bias it so that floor division (as done by an + * arithmetic right shift) and modulus (as done by bitwise-AND) does the + * right thing. + */ +#define CARRY(i, wd, b, m) do { \ + x[i] += (b); \ + c = ASR(dblpiece, x[i], (wd)); \ + x[i] = (dblpiece)((udblpiece)x[i]&(m)) - (b); \ +} while (0) + + { CARRY(23, 10, B9, M10); } + { x[24] += c; CARRY(24, 10, B9, M10); } + { x[25] += c; CARRY(25, 9, B8, M9); } + { x[0] += 19*c; CARRY( 0, 10, B9, M10); } + for (i = 1; i < 21; ) { + for (j = i + 4; i < j; ) { x[i] += c; CARRY( i, 10, B9, M10); i++; } + { x[i] += c; CARRY( i, 9, B8, M9); i++; } + } + while (i < 25) { x[i] += c; CARRY( i, 10, B9, M10); i++; } + x[25] += c; + +#undef CARRY +} + +#endif + +/* --- @f25519_mulconst@ --- * + * + * Arguments: @f25519 *z@ = where to put the result (may alias @x@) + * @const f25519 *x@ = an operand + * @long a@ = a small-ish constant; %$|a| < 2^{20}$%. + * + * Returns: --- + * + * Use: Set @z@ to the product %$a x$%. + */ + +void f25519_mulconst(f25519 *z, const f25519 *x, long a) +{ +#if F25519_IMPL == 26 + + piece PIECES(x); + dblpiece PIECES(z), aa = a; + + FETCH(x, x); + + /* Suppose that |x_i| <= 2^27, and |a| <= 2^23. Then we'll have + * |z_i| <= 2^50. + */ + z0 = aa*x0; z1 = aa*x1; z2 = aa*x2; z3 = aa*x3; z4 = aa*x4; + z5 = aa*x5; z6 = aa*x6; z7 = aa*x7; z8 = aa*x8; z9 = aa*x9; + + /* Following `CARRY_REDUCE', we'll have |z_i| <= 2^26. */ + CARRY_REDUCE(z, z); + STASH(z, z); + +#elif F25519_IMPL == 10 + + dblpiece y[NPIECE]; + unsigned i; + + for (i = 0; i < NPIECE; i++) y[i] = a*x->P[i]; + carry_reduce(y); + for (i = 0; i < NPIECE; i++) z->P[i] = y[i]; + +#endif +} + +/* --- @f25519_mul@ --- * + * + * Arguments: @f25519 *z@ = where to put the result (may alias @x@ or @y@) + * @const f25519 *x, *y@ = two operands + * + * Returns: --- + * + * Use: Set @z@ to the product %$x y$%. + */ + +void f25519_mul(f25519 *z, const f25519 *x, const f25519 *y) +{ +#if F25519_IMPL == 26 + + piece PIECES(x), PIECES(y); + dblpiece PIECES(z); + unsigned i; + + FETCH(x, x); FETCH(y, y); + + /* Suppose that |x_i|, |y_i| <= 2^27. Then we'll have + * + * |z_0| <= 267*2^54 + * |z_1| <= 154*2^54 + * |z_2| <= 213*2^54 + * |z_3| <= 118*2^54 + * |z_4| <= 159*2^54 + * |z_5| <= 82*2^54 + * |z_6| <= 105*2^54 + * |z_7| <= 46*2^54 + * |z_8| <= 51*2^54 + * |z_9| <= 10*2^54 + * + * all of which are less than 2^63 - 2^25. + */ + +#define M(a, b) ((dblpiece)(a)*(b)) + z0 = M(x0, y0) + + 19*(M(x2, y8) + M(x4, y6) + M(x6, y4) + M(x8, y2)) + + 38*(M(x1, y9) + M(x3, y7) + M(x5, y5) + M(x7, y3) + M(x9, y1)); + z1 = M(x0, y1) + M(x1, y0) + + 19*(M(x2, y9) + M(x3, y8) + M(x4, y7) + M(x5, y6) + + M(x6, y5) + M(x7, y4) + M(x8, y3) + M(x9, y2)); + z2 = M(x0, y2) + M(x2, y0) + + 2* M(x1, y1) + + 19*(M(x4, y8) + M(x6, y6) + M(x8, y4)) + + 38*(M(x3, y9) + M(x5, y7) + M(x7, y5) + M(x9, y3)); + z3 = M(x0, y3) + M(x1, y2) + M(x2, y1) + M(x3, y0) + + 19*(M(x4, y9) + M(x5, y8) + M(x6, y7) + + M(x7, y6) + M(x8, y5) + M(x9, y4)); + z4 = M(x0, y4) + M(x2, y2) + M(x4, y0) + + 2*(M(x1, y3) + M(x3, y1)) + + 19*(M(x6, y8) + M(x8, y6)) + + 38*(M(x5, y9) + M(x7, y7) + M(x9, y5)); + z5 = M(x0, y5) + M(x1, y4) + M(x2, y3) + + M(x3, y2) + M(x4, y1) + M(x5, y0) + + 19*(M(x6, y9) + M(x7, y8) + M(x8, y7) + M(x9, y6)); + z6 = M(x0, y6) + M(x2, y4) + M(x4, y2) + M(x6, y0) + + 2*(M(x1, y5) + M(x3, y3) + M(x5, y1)) + + 19* M(x8, y8) + + 38*(M(x7, y9) + M(x9, y7)); + z7 = M(x0, y7) + M(x1, y6) + M(x2, y5) + M(x3, y4) + + M(x4, y3) + M(x5, y2) + M(x6, y1) + M(x7, y0) + + 19*(M(x8, y9) + M(x9, y8)); + z8 = M(x0, y8) + M(x2, y6) + M(x4, y4) + M(x6, y2) + M(x8, y0) + + 2*(M(x1, y7) + M(x3, y5) + M(x5, y3) + M(x7, y1)) + + 38* M(x9, y9); + z9 = M(x0, y9) + M(x1, y8) + M(x2, y7) + M(x3, y6) + M(x4, y5) + + M(x5, y4) + M(x6, y3) + M(x7, y2) + M(x8, y1) + M(x9, y0); +#undef M + + /* From above, we have |z_i| <= 2^63 - 2^25. A pass of `CARRY_REDUCE' will + * leave |z_i| <= 2^38 + 2^25; and a second pass will leave |z_i| <= 2^25 + + * 2^13, which is comfortable for an addition prior to the next + * multiplication. + */ + for (i = 0; i < 2; i++) CARRY_REDUCE(z, z); + STASH(z, z); + +#elif F25519_IMPL == 10 + + dblpiece u[NPIECE], t, tt, p; + unsigned i, j, k; + + /* This is unpleasant. Honestly, this table seems to be the best way of + * doing it. + */ + static const unsigned short off[NPIECE] = { + 0, 10, 20, 30, 40, 50, 59, 69, 79, 89, 99, 108, 118, + 128, 138, 148, 157, 167, 177, 187, 197, 206, 216, 226, 236, 246 + }; + + /* First pass: things we must multiply by 19 or 38. */ + for (i = 0; i < NPIECE - 1; i++) { + t = tt = 0; + for (j = i + 1; j < NPIECE; j++) { + k = NPIECE + i - j; p = (dblpiece)x->P[j]*y->P[k]; + if (off[i] < off[j] + off[k] - 255) tt += p; + else t += p; + } + u[i] = 19*(t + 2*tt); + } + u[NPIECE - 1] = 0; + + /* Second pass: things we must multiply by 1 or 2. */ + for (i = 0; i < NPIECE; i++) { + t = tt = 0; + for (j = 0; j <= i; j++) { + k = i - j; p = (dblpiece)x->P[j]*y->P[k]; + if (off[i] < off[j] + off[k]) tt += p; + else t += p; + } + u[i] += t + 2*tt; + } + + /* And we're done. */ + carry_reduce(u); + for (i = 0; i < NPIECE; i++) z->P[i] = u[i]; + +#endif +} + +/* --- @f25519_sqr@ --- * + * + * Arguments: @f25519 *z@ = where to put the result (may alias @x@ or @y@) + * @const f25519 *x@ = an operand + * + * Returns: --- + * + * Use: Set @z@ to the square %$x^2$%. + */ + +void f25519_sqr(f25519 *z, const f25519 *x) +{ +#if F25519_IMPL == 26 + + piece PIECES(x); + dblpiece PIECES(z); + unsigned i; + + FETCH(x, x); + + /* See `f25519_mul' for bounds. */ + +#define M(a, b) ((dblpiece)(a)*(b)) + z0 = M(x0, x0) + + 38*(M(x2, x8) + M(x4, x6) + M(x5, x5)) + + 76*(M(x1, x9) + M(x3, x7)); + z1 = 2* M(x0, x1) + + 38*(M(x2, x9) + M(x3, x8) + M(x4, x7) + M(x5, x6)); + z2 = 2*(M(x0, x2) + M(x1, x1)) + + 19* M(x6, x6) + + 38* M(x4, x8) + + 76*(M(x3, x9) + M(x5, x7)); + z3 = 2*(M(x0, x3) + M(x1, x2)) + + 38*(M(x4, x9) + M(x5, x8) + M(x6, x7)); + z4 = M(x2, x2) + + 2* M(x0, x4) + + 4* M(x1, x3) + + 38*(M(x6, x8) + M(x7, x7)) + + 76* M(x5, x9); + z5 = 2*(M(x0, x5) + M(x1, x4) + M(x2, x3)) + + 38*(M(x6, x9) + M(x7, x8)); + z6 = 2*(M(x0, x6) + M(x2, x4) + M(x3, x3)) + + 4* M(x1, x5) + + 19* M(x8, x8) + + 76* M(x7, x9); + z7 = 2*(M(x0, x7) + M(x1, x6) + M(x2, x5) + M(x3, x4)) + + 38* M(x8, x9); + z8 = M(x4, x4) + + 2*(M(x0, x8) + M(x2, x6)) + + 4*(M(x1, x7) + M(x3, x5)) + + 38* M(x9, x9); + z9 = 2*(M(x0, x9) + M(x1, x8) + M(x2, x7) + M(x3, x6) + M(x4, x5)); +#undef M + + /* See `f25519_mul' for details. */ + for (i = 0; i < 2; i++) CARRY_REDUCE(z, z); + STASH(z, z); + +#elif F25519_IMPL == 10 + f25519_mul(z, x, x); +#endif +} + +/*----- More complicated things -------------------------------------------*/ + +/* --- @f25519_inv@ --- * + * + * Arguments: @f25519 *z@ = where to put the result (may alias @x@) + * @const f25519 *x@ = an operand + * + * Returns: --- + * + * Use: Stores in @z@ the multiplicative inverse %$x^{-1}$%. If + * %$x = 0$% then @z@ is set to zero. This is considered a + * feature. + */ + +void f25519_inv(f25519 *z, const f25519 *x) +{ + f25519 t, u, t2, t11, t2p10m1, t2p50m1; + unsigned i; + +#define SQRN(z, x, n) do { \ + f25519_sqr((z), (x)); \ + for (i = 1; i < (n); i++) f25519_sqr((z), (z)); \ +} while (0) + + /* Calculate x^-1 = x^(p - 2) = x^(2^255 - 21), which also handles x = 0 as + * intended. The addition chain here is from Bernstein's implementation; I + * couldn't find a better one. + */ /* step | value */ + f25519_sqr(&t2, x); /* 1 | 2 */ + SQRN(&u, &t2, 2); /* 3 | 8 */ + f25519_mul(&t, &u, x); /* 4 | 9 */ + f25519_mul(&t11, &t, &t2); /* 5 | 11 = 2^5 - 21 */ + f25519_sqr(&u, &t11); /* 6 | 22 */ + f25519_mul(&t, &t, &u); /* 7 | 31 = 2^5 - 1 */ + SQRN(&u, &t, 5); /* 12 | 2^10 - 2^5 */ + f25519_mul(&t2p10m1, &t, &u); /* 13 | 2^10 - 1 */ + SQRN(&u, &t2p10m1, 10); /* 23 | 2^20 - 2^10 */ + f25519_mul(&t, &t2p10m1, &u); /* 24 | 2^20 - 1 */ + SQRN(&u, &t, 20); /* 44 | 2^40 - 2^20 */ + f25519_mul(&t, &t, &u); /* 45 | 2^40 - 1 */ + SQRN(&u, &t, 10); /* 55 | 2^50 - 2^10 */ + f25519_mul(&t2p50m1, &t2p10m1, &u); /* 56 | 2^50 - 1 */ + SQRN(&u, &t2p50m1, 50); /* 106 | 2^100 - 2^50 */ + f25519_mul(&t, &t2p50m1, &u); /* 107 | 2^100 - 1 */ + SQRN(&u, &t, 100); /* 207 | 2^200 - 2^100 */ + f25519_mul(&t, &t, &u); /* 208 | 2^200 - 1 */ + SQRN(&u, &t, 50); /* 258 | 2^250 - 2^50 */ + f25519_mul(&t, &t2p50m1, &u); /* 259 | 2^250 - 1 */ + SQRN(&u, &t, 5); /* 264 | 2^255 - 2^5 */ + f25519_mul(z, &u, &t11); /* 265 | 2^255 - 21 */ + +#undef SQRN +} + +/*----- Test rig ----------------------------------------------------------*/ + +#ifdef TEST_RIG + +#include +#include + +static void fixdstr(dstr *d) +{ + if (d->len > 32) + die(1, "invalid length for f25519"); + else if (d->len < 32) { + dstr_ensure(d, 32); + memset(d->buf + d->len, 0, 32 - d->len); + d->len = 32; + } +} + +static void cvt_f25519(const char *buf, dstr *d) +{ + dstr dd = DSTR_INIT; + + type_hex.cvt(buf, &dd); fixdstr(&dd); + dstr_ensure(d, sizeof(f25519)); d->len = sizeof(f25519); + f25519_load((f25519 *)d->buf, (const octet *)dd.buf); + dstr_destroy(&dd); +} + +static void dump_f25519(dstr *d, FILE *fp) + { fdump(stderr, "???", (const piece *)d->buf); } + +static void cvt_f25519_ref(const char *buf, dstr *d) + { type_hex.cvt(buf, d); fixdstr(d); } + +static void dump_f25519_ref(dstr *d, FILE *fp) +{ + f25519 x; + + f25519_load(&x, (const octet *)d->buf); + fdump(stderr, "???", x.P); +} + +static int eq(const f25519 *x, dstr *d) + { octet b[32]; f25519_store(b, x); return (memcmp(b, d->buf, 32) == 0); } + +static const test_type + type_f25519 = { cvt_f25519, dump_f25519 }, + type_f25519_ref = { cvt_f25519_ref, dump_f25519_ref }; + +#define TEST_UNOP(op) \ + static int vrf_##op(dstr dv[]) \ + { \ + f25519 *x = (f25519 *)dv[0].buf; \ + f25519 z, zz; \ + int ok = 1; \ + \ + f25519_##op(&z, x); \ + if (!eq(&z, &dv[1])) { \ + ok = 0; \ + fprintf(stderr, "failed!\n"); \ + fdump(stderr, "x", x->P); \ + fdump(stderr, "calc", z.P); \ + f25519_load(&zz, (const octet *)dv[1].buf); \ + fdump(stderr, "z", zz.P); \ + } \ + \ + return (ok); \ + } + +TEST_UNOP(sqr) +TEST_UNOP(inv) + +#define TEST_BINOP(op) \ + static int vrf_##op(dstr dv[]) \ + { \ + f25519 *x = (f25519 *)dv[0].buf, *y = (f25519 *)dv[1].buf; \ + f25519 z, zz; \ + int ok = 1; \ + \ + f25519_##op(&z, x, y); \ + if (!eq(&z, &dv[2])) { \ + ok = 0; \ + fprintf(stderr, "failed!\n"); \ + fdump(stderr, "x", x->P); \ + fdump(stderr, "y", y->P); \ + fdump(stderr, "calc", z.P); \ + f25519_load(&zz, (const octet *)dv[2].buf); \ + fdump(stderr, "z", zz.P); \ + } \ + \ + return (ok); \ + } + +TEST_BINOP(add) +TEST_BINOP(sub) +TEST_BINOP(mul) + +static int vrf_mulc(dstr dv[]) +{ + f25519 *x = (f25519 *)dv[0].buf; + long a = *(const long *)dv[1].buf; + f25519 z, zz; + int ok = 1; + + f25519_mulconst(&z, x, a); + if (!eq(&z, &dv[2])) { + ok = 0; + fprintf(stderr, "failed!\n"); + fdump(stderr, "x", x->P); + fprintf(stderr, "a = %ld\n", a); + fdump(stderr, "calc", z.P); + f25519_load(&zz, (const octet *)dv[2].buf); + fdump(stderr, "z", zz.P); + } + + return (ok); +} + +static int vrf_condswap(dstr dv[]) +{ + f25519 *x = (f25519 *)dv[0].buf, *y = (f25519 *)dv[1].buf; + f25519 xx = *x, yy = *y; + uint32 m = *(uint32 *)dv[2].buf; + int ok = 1; + + f25519_condswap(&xx, &yy, m); + if (!eq(&xx, &dv[3]) || !eq(&yy, &dv[4])) { + ok = 0; + fprintf(stderr, "failed!\n"); + fdump(stderr, "x", x->P); + fdump(stderr, "y", y->P); + fprintf(stderr, "m = 0x%08lx\n", (unsigned long)m); + fdump(stderr, "calc xx", xx.P); + fdump(stderr, "calc yy", yy.P); + f25519_load(&xx, (const octet *)dv[3].buf); + f25519_load(&yy, (const octet *)dv[4].buf); + fdump(stderr, "want xx", xx.P); + fdump(stderr, "want yy", yy.P); + } + + return (ok); +} + +static int vrf_sub_mulc_add_sub_mul(dstr dv[]) +{ + f25519 *u = (f25519 *)dv[0].buf, *v = (f25519 *)dv[1].buf, + *w = (f25519 *)dv[3].buf, *x = (f25519 *)dv[4].buf, + *y = (f25519 *)dv[5].buf; + long a = *(const long *)dv[2].buf; + f25519 umv, aumv, wpaumv, xmy, z, zz; + int ok = 1; + + f25519_sub(&umv, u, v); + f25519_mulconst(&aumv, &umv, a); + f25519_add(&wpaumv, w, &aumv); + f25519_sub(&xmy, x, y); + f25519_mul(&z, &wpaumv, &xmy); + + if (!eq(&z, &dv[6])) { + ok = 0; + fprintf(stderr, "failed!\n"); + fdump(stderr, "u", u->P); + fdump(stderr, "v", v->P); + fdump(stderr, "u - v", umv.P); + fprintf(stderr, "a = %ld\n", a); + fdump(stderr, "a (u - v)", aumv.P); + fdump(stderr, "w + a (u - v)", wpaumv.P); + fdump(stderr, "x", x->P); + fdump(stderr, "y", y->P); + fdump(stderr, "x - y", xmy.P); + fdump(stderr, "(x - y) (w + a (u - v))", z.P); + f25519_load(&zz, (const octet *)dv[6].buf); fdump(stderr, "z", zz.P); + } + + return (ok); +} + +static test_chunk tests[] = { + { "add", vrf_add, { &type_f25519, &type_f25519, &type_f25519_ref } }, + { "sub", vrf_sub, { &type_f25519, &type_f25519, &type_f25519_ref } }, + { "mul", vrf_mul, { &type_f25519, &type_f25519, &type_f25519_ref } }, + { "mulconst", vrf_mulc, { &type_f25519, &type_long, &type_f25519_ref } }, + { "condswap", vrf_condswap, + { &type_f25519, &type_f25519, &type_uint32, + &type_f25519_ref, &type_f25519_ref } }, + { "sqr", vrf_sqr, { &type_f25519, &type_f25519_ref } }, + { "inv", vrf_inv, { &type_f25519, &type_f25519_ref } }, + { "sub-mulc-add-sub-mul", vrf_sub_mulc_add_sub_mul, + { &type_f25519, &type_f25519, &type_long, &type_f25519, + &type_f25519, &type_f25519, &type_f25519_ref } }, + { 0, 0, { 0 } } +}; + +int main(int argc, char *argv[]) +{ + test_run(argc, argv, tests, SRCDIR "/t/f25519"); + return (0); +} + +#endif + +/*----- That's all, folks -------------------------------------------------*/ diff --git a/math/f25519.h b/math/f25519.h new file mode 100644 index 00000000..a70ca2df --- /dev/null +++ b/math/f25519.h @@ -0,0 +1,209 @@ +/* -*-c-*- + * + * Arithmetic modulo 2^255 - 19 + * + * (c) 2017 Straylight/Edgeware + */ + +/*----- Licensing notice --------------------------------------------------* + * + * This file is part of Catacomb. + * + * Catacomb is free software; you can redistribute it and/or modify + * it under the terms of the GNU Library General Public License as + * published by the Free Software Foundation; either version 2 of the + * License, or (at your option) any later version. + * + * Catacomb is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU Library General Public License for more details. + * + * You should have received a copy of the GNU Library General Public + * License along with Catacomb; if not, write to the Free + * Software Foundation, Inc., 59 Temple Place - Suite 330, Boston, + * MA 02111-1307, USA. + */ + +#ifndef CATACOMB_F25519_H +#define CATACOMB_F25519_H + +#ifdef __cplusplus + extern "C" { +#endif + +/*----- Header files ------------------------------------------------------*/ + +#include + +#ifndef CATACOMB_QFARITH_H +# include "qfarith.h" +#endif + +/*----- Data structures ---------------------------------------------------*/ + +typedef union { + int32 p26[10]; + int16 p10[26]; +} f25519; + +#if !defined(F25519_IMPL) && defined(HAVE_INT64) +# define F25519_IMPL 26 +#endif + +#ifndef F25519_IMPL +# define F25519_IMPL 10 +#endif + +/*----- Functions provided ------------------------------------------------*/ + +/* --- @f25519_set@ --- * + * + * Arguments: @f25519 *z@ = where to write the result + * @int a@ = a small-ish constant + * + * Returns: --- + * + * Use: Sets @z@ to equal @a@. + */ + +extern void f25519_set(f25519 */*x*/, int /*a*/); + +/* --- @f25519_load@ --- * + * + * Arguments: @f25519 *z@ = where to store the result + * @const octet xv[32]@ = source to read + * + * Returns: --- + * + * Use: Reads an element of %$\gf{2^{255} - 19}$% in external + * representation from @xv@ and stores it in @z@. + * + * External representation is little-endian base-256. Elements + * have multiple encodings, which are not produced by correct + * software; use of noncanonical encodings is not an error, and + * toleration of them is considered a performance feature. + * + * Some specifications, e.g., RFC7748, require the topmost bit + * (i.e., bit 7 of @wv[31]@) to be ignored. Callers + * implementing such specifications should clear the bit + * explicitly. (It's much easier for a caller who wants the bit + * to be ignored to clear it than for a caller who wants the bit + * to be significant to apply the necessary change by hand.) + */ + +extern void f25519_load(f25519 */*z*/, const octet /*xv*/[32]); + +/* --- @f25519_store@ --- * + * + * Arguments: @octet zv[32]@ = where to write the result + * @const f25519 *x@ = the field element to write + * + * Returns: --- + * + * Use: Stores a field element in the given octet vector in external + * representation. A canonical encoding is always stored, so, + * in particular, the top bit of @xv[31]@ is always left clear. + */ + +extern void f25519_store(octet /*zv*/[32], const f25519 */*x*/); + +/* --- @f25519_condswap@ --- * + * + * Arguments: @f25519 *x, *y@ = two operands + * @uint32 m@ = a mask + * + * Returns: --- + * + * Use: If @m@ is zero, do nothing; if @m@ is all-bits-set, then + * exchange @x@ and @y@. If @m@ has some other value, then + * scramble @x@ and @y@ in an unhelpful way. + */ + +extern void f25519_condswap(f25519 */*x*/, f25519 */*y*/, uint32 /*m*/); + +/* --- @f25519_add@ --- * + * + * Arguments: @f25519 *z@ = where to put the result (may alias @x@ or @y@) + * @const f25519 *x, *y@ = two operands + * + * Returns: --- + * + * Use: Set @z@ to the sum %$x + y$%. + */ + +extern void f25519_add(f25519 */*z*/, + const f25519 */*x*/, const f25519 */*y*/); + +/* --- @f25519_sub@ --- * + * + * Arguments: @f25519 *z@ = where to put the result (may alias @x@ or @y@) + * @const f25519 *x, *y@ = two operands + * + * Returns: --- + * + * Use: Set @z@ to the difference %$x - y$%. + */ + +extern void f25519_sub(f25519 */*z*/, + const f25519 */*x*/, const f25519 */*y*/); + +/* --- @f25519_mulconst@ --- * + * + * Arguments: @f25519 *z@ = where to put the result (may alias @x@) + * @const f25519 *x@ = an operand + * @long a@ = a small-ish constant; %$|a| < 2^{20}$%. + * + * Returns: --- + * + * Use: Set @z@ to the product %$a x$%. + */ + +extern void f25519_mulconst(f25519 */*z*/, const f25519 */*x*/, long /*a*/); + +/* --- @f25519_mul@ --- * + * + * Arguments: @f25519 *z@ = where to put the result (may alias @x@ or @y@) + * @const f25519 *x, *y@ = two operands + * + * Returns: --- + * + * Use: Set @z@ to the product %$x y$%. + */ + +extern void f25519_mul(f25519 */*z*/, + const f25519 */*x*/, const f25519 */*y*/); + +/* --- @f25519_sqr@ --- * + * + * Arguments: @f25519 *z@ = where to put the result (may alias @x@ or @y@) + * @const f25519 *x@ = an operand + * + * Returns: --- + * + * Use: Set @z@ to the square %$x^2$%. + */ + +extern void f25519_sqr(f25519 */*z*/, const f25519 */*x*/); + +/* --- @f25519_inv@ --- * + * + * Arguments: @f25519 *z@ = where to put the result (may alias @x@) + * @const f25519 *x@ = an operand + * + * Returns: --- + * + * Use: Stores in @z@ the multiplicative inverse %$x^{-1}$%. If + * %$x = 0$% then @z@ is set to zero. This is considered a + * feature. + */ + +extern void f25519_inv(f25519 */*z*/, const f25519 */*x*/); + +/*----- That's all, folks -------------------------------------------------*/ + +#ifdef __cplusplus + } +#endif + +#endif diff --git a/math/qfarith.h b/math/qfarith.h new file mode 100644 index 00000000..abcd293c --- /dev/null +++ b/math/qfarith.h @@ -0,0 +1,229 @@ +/* -*-c-*- + * + * Utilities for quick field arithmetic + * + * (c) 2017 Straylight/Edgeware + */ + +/*----- Licensing notice --------------------------------------------------* + * + * This file is part of Catacomb. + * + * Catacomb is free software; you can redistribute it and/or modify + * it under the terms of the GNU Library General Public License as + * published by the Free Software Foundation; either version 2 of the + * License, or (at your option) any later version. + * + * Catacomb is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU Library General Public License for more details. + * + * You should have received a copy of the GNU Library General Public + * License along with Catacomb; if not, write to the Free + * Software Foundation, Inc., 59 Temple Place - Suite 330, Boston, + * MA 02111-1307, USA. + */ + +#ifndef CATACOMB_QFARITH_H +#define CATACOMB_QFARITH_H + +#ifdef __cplusplus + extern "C" { +#endif + +/*----- Header files ------------------------------------------------------*/ + +#include + +#include + +/*----- Signed integer types ----------------------------------------------*/ + +/* See if we can find a suitable 64-bit or wider type. Don't bother if we + * don't have a corresponding unsigned type, because we really need both. + */ +#ifdef HAVE_UINT64 +# if INT_MAX >> 31 == 0xffffffff +# define HAVE_INT64 + typedef int int64; +# elif LONG_MAX >> 31 == 0xffffffff +# define HAVE_INT64 + typedef long int64; +# elif defined(LLONG_MAX) +# define HAVE_INT64 + MLIB_BITS_EXTENSION typedef long long int64; +# endif +#endif + +/* Choose suitable 32- and 16-bit types. */ +#if INT_MAX >= 0x7fffffff + typedef int int32; +#else + typedef long int32; +#endif + +typedef short int16; + +/*----- General bit-hacking utilities -------------------------------------*/ + +/* Individual bits, and masks for low bits. */ +#define BIT(n) (1ul << (n)) +#define MASK(n) (BIT(n) - 1) + +/* Arithmetic right shift. If X is a value of type TY, and N is a + * nonnegative integer, then return the value of X shifted right by N bits; + * alternatively, this is floor(X/2^N). + * + * GCC manages to compile this into a simple shift operation if one is + * available, but it's correct for all C targets. + */ +#define ASR(ty, x, n) (((x) - (ty)((u##ty)(x)&MASK(n)))/(ty)BIT(n)) + +/*----- Constant-time utilities -------------------------------------------*/ + +/* The following have better implementations on a two's complement target. */ +#ifdef NEG_TWOC + + /* If we have two's complement arithmetic then masks are signed; this + * avoids a switch to unsigned representation, with the consequent problem + * of overflow when we convert back. + */ + typedef int32 mask32; + + /* Convert an unsigned mask M into a `mask32'. This is a hairy-looking + * no-op on many targets, but, given that we have two's complement + * integers, it is free of arithmetic overflow. + */ +# define FIX_MASK32(m) \ + ((mask32)((m)&0x7fffffffu) + (-(mask32)0x7fffffff - 1)*(((m) >> 31)&1u)) + + /* If Z is zero and M has its low 32 bits set, then copy (at least) the low + * 32 bits of X to Z; if M is zero, do nothing. Otherwise, scramble Z + * unhelpfully. + */ +# define CONDPICK(z, x, m) do { (z) |= (x)&(m); } while (0) + + /* If M has its low 32 bits set, then return (at least) the low 32 bits of + * X in Z; if M is zero, then return (at least) the low 32 bits of Y in Z. + * Otherwise, return an unhelful result. + */ +# define PICK2(x, y, m) (((x)&(m)) | ((y)&~(m))) + + /* If M has its low 32 bits set then swap (at least) the low 32 bits of X + * and Y; if M is zero, then do nothing. Otherwise, scramble X and Y + * unhelpfully. + */ +# define CONDSWAP(x, y, m) \ + do { mask32 t_ = ((x) ^ (y))&(m); (x) ^= t_; (y) ^= t_; } while (0) +#else + + /* We don't have two's complement arithmetic. We can't use bithacking at + * all: if we try to hack on the bits of signed numbers we'll come unstuck + * when we hit the other representation of zero; and if we switch to + * unsigned arithmetic then we'll have overflow when we try to convert a + * negative number back. So fall back to simple arithmetic. + */ + typedef uint32 mask32; + + /* Convert an unsigned mask M into a `mask32'. Our masks are unsigned, so + * this does nothing. + */ +# define FIX_MASK32(m) ((mask32)(m)) + + /* If Z is zero and M has its low 32 bits set, then copy (at least) the low + * 32 bits of X to Z; if M is zero, do nothing. Otherwise, scramble Z + * unhelpfully. + */ +# define CONDPICK(z, x, m) \ + do { (z) += (x)*(int)((unsigned)(m)&1u); } while (0) + + /* If M has its low 32 bits set, then return (at least) the low 32 bits of + * X in Z; if M is zero, then return (at least) the low 32 bits of Y in Z. + * Otherwise, return an unhelful result. + */ +# define PICK2(x, y, m) \ + ((x)*(int)((unsigned)(m)&1u) + (y)*(int)(1 - ((unsigned)(m)&1u))) + + /* If M has its low 32 bits set then swap (at least) the low 32 bits of X + * and Y; if M is zero, then do nothing. Otherwise, scramble X and Y + * unhelpfully. + */ +# define CONDSWAP(x, y, m) do { \ + int32 x_ = PICK2((y), (x), (m)), y_ = PICK2((x), (y), (m)); \ + (x) = x_; (y) = y_; \ + } while (0) +#endif + +/* Return zero if bit 31 of X is clear, or a mask with (at least) the low 32 + * bits set if bit 31 of X is set. + */ +#define SIGN(x) (-(mask32)(((uint32)(x) >> 31)&1)) + +/* Return zero if X is zero, or a mask with (at least) the low 32 bits set if + * X is nonzero. + */ +#define NONZEROP(x) SIGN((U32(x) >> 1) - U32(x)) + +/*----- Debugging utilities -----------------------------------------------*/ + +/* Define a debugging function DUMPFN, which will dump an integer represented + * modulo M. The integer is represented as a vector of NPIECE pieces of type + * PIECETY. The pieces are assembled at possibly irregular offsets: piece i + * logically has width PIECEWD(i), but may overhang the next piece. The + * pieces may be signed. GETMOD is an expression which calculates and + * returns the value of M, as an `mp *'. + * + * The generated function writes the value of such an integer X to the stream + * FP, labelled with the string WHAT. + * + * The definition assumes that , , and + * have been included. + */ +#define DEF_FDUMP(dumpfn, piecety, piecewd, npiece, noctet, getmod) \ + static void dumpfn(FILE *fp, const char *what, const piecety *x) \ + { \ + mpw w; \ + mp m, *y = MP_ZERO, *t = MP_NEW, *p; \ + octet b[noctet]; \ + unsigned i, o; \ + \ + p = getmod; \ + mp_build(&m, &w, &w + 1); \ + for (i = o = 0; i < npiece; i++) { \ + if (x[i] >= 0) { w = x[i]; m.f &= ~MP_NEG; } \ + else { w = -x[i]; m.f |= MP_NEG; } \ + t = mp_lsl(t, &m, o); \ + y = mp_add(y, y, t); \ + o += piecewd(i); \ + } \ + \ + fprintf(fp, "%s = <", what); \ + for (i = 0; i < npiece; i++) { \ + if (i) fputs(", ", fp); \ + fprintf(fp, "%ld", (long)x[i]); \ + } \ + fputs(">\n\t= ", fp); \ + mp_writefile(y, fp, 10); \ + fputs("\n\t== ", fp); \ + mp_div(0, &y, y, p); \ + mp_writefile(y, fp, 10); \ + fputs("\n\t= 0x", fp); \ + mp_writefile(y, fp, 16); \ + fputs(" (mod 2^255 - 19)\n\t= [", fp); \ + mp_storel(y, b, sizeof(b)); \ + for (i = 0; i < 32; i++) { \ + if (i && !(i%4)) fputc(':', fp); \ + fprintf(fp, "%02x", b[i]); \ + } \ + fputs("]\n", fp); \ + mp_drop(y); mp_drop(p); mp_drop(t); \ + } + +/*----- That's all, folks -------------------------------------------------*/ + +#ifdef __cplusplus + } +#endif + +#endif diff --git a/math/t/f25519 b/math/t/f25519 new file mode 100644 index 00000000..a2264bc6 --- /dev/null +++ b/math/t/f25519 @@ -0,0 +1,594 @@ +### Test cases for arithmetic mod 2^255 - 19. + +add { + ## Some load and store tests. + 0102030405060708090a0b0c0d0e0f101112131415161718191a1b1c1d1e1f20 + 0000000000000000000000000000000000000000000000000000000000000000 + 0102030405060708090a0b0c0d0e0f101112131415161718191a1b1c1d1e1f20; + ecffffffffffffffffffffffffffffffffffffffffffffffffffffffffffff7f + 0000000000000000000000000000000000000000000000000000000000000000 + ecffffffffffffffffffffffffffffffffffffffffffffffffffffffffffff7f; + edffffffffffffffffffffffffffffffffffffffffffffffffffffffffffff7f + 00000000000000000000000000000000000000000000000000000000000000 + 00; + + ## Random tests. + bb8cd45b6185ee8ca34c4d7560425fd5c2fc22871113c6427198e8f9439a8523 + 98f16cc5e44b7228b7f3b6ed91930d68bba5071157d19259c48034d2817e30d5 + 667e412146d160b55a400463f2d56c3d7ea22a9868e4589c35191dccc518b678; + 8ef938a7b07a668c234c2c5bb65115b8e2e22a3703cda51c70e7a56a075d4f34 + 1f3cd4bf4f768058c4de1aa5aeb90e227338e308c7a0a199050ca1c38d1e034c + c0350d6700f1e6e4e72a4700650b24da551b0e40ca6d47b675f3462e957b5200; + 8b2c4678fd4f8b7254eea925e4da8d472db7f97d349d1b2134ced884b747ea7d + 58dca7954f843b7e4dc56b048903e7728df66ef6824f29932b2f568c995ddf2c + f608ee0d4dd4c6f0a1b3152a6dde74babaad6874b7ec44b45ffd2e1151a5c92a; + 9cd64f0cf439bccd224f7a789396b9f888998d1ba9759c93a5d6a27f45db0e67 + 12f724a127052612081e27e3915db9a540de9b64c566edc313914ad94424464e + c1cd74ad1b3fe2df2a6da15b25f4729ec97729806edc8957b967ed588aff5435; + 82e59fa9703b2557303665fe24e3d46259da974d78155f05c6b8891f92fc808d + 98a073ec62baf1354b681a3c78d62df1698ed5077dbe2bd68ba4a6e0233e3447 + 2d861396d3f5168d7b9e7f3a9db90254c3686d55f5d38adb515d3000b63ab554; + 686ceca221d0f4c6d98287e97c506085064e7b172672b88ee555dd0308bfd309 + 60adc3675be5e6c2921a74b5b8b24b0a26341a19567dae713a720786bc72df16 + c819b00a7db5db896c9dfb9e3503ac8f2c8295307cef660020c8e489c431b320; + a50a5b82d2451ed9a5844585712202a6dc0e0230501d57b8c00f4bd4771aed9e + 98d71b8b140aea02ad47d9b19e47d2a03257adc1d5c47611adc667652b78c87c + 63e2760de74f08dc52cc1e37106ad4460f66aff125e2cdc96dd6b239a392b51b; + 1637e3eb4c5d173b8cba29bdb13b930f48ae0599e354819f8ae0cf70f0f7b61b + ffeef1767df70de7a4a75864c03935d7934c6cf5228b9c3000813ed0e55447e2 + 2826d562ca542522316282217275c8e6dbfa718e06e01dd08a610e41d64cfe7d; + 9e8bf571c45bc5f37c9126e70498698978903b94d9c6712baef91894ab924fda + 480f68f8ee644bfb64e3fa514eaf16aea6368da6813802c4c3e6ce8f6540dfa3 + 0c9b5d6ab3c010efe1742139534780371fc7c83a5bff73ef71e0e72311d32e7e; + df9bdbd478bd1ec1ff078cdc76f2a5fcdc52f0198ba412e0b5a10e664b175dad + 080bc47bd2fb2ba3be6a07a02162f76020479f2cc6469c6766df23acf9517896 + 0da79f504bb94a64be72937c98549d5dfd998f4651ebae471c8132124569d543; + b6dd223f5d8a3dc5dd6f7a0679173ff656303017a17f17f4086fd4ecd3a23020 + 560d827aad23c26fddbfd16c3514a9b518a4586d6e3c5739860743662b14f4b4 + 1feba4b90aaeff34bb2f4c73ae2be8ab6fd488840fbc6e2d8f761753ffb62455; + 1d59691a670c40143ca7f7e693389ab64a911b17a2fe4368fba9097b40864162 + cb42d3ba12a8b90ef1622ea09241d3157c515043040914356d6d99b95089fa07 + e89b3cd579b4f9222d0a2687267a6dccc6e26b5aa607589d6817a334910f3c6a; + 82b48877e9c8c18edb52db3e1417d894e572a8edd84335d93c642158f2a926e0 + 93c6a39c1ef97f318c3b13bf8e77bea653b2eb130bdf43b00a0e0e40a362e819 + 287b2c1408c241c0678eeefda28e963b39259401e422798947722f98950c0f7a; + a6137481c5b9b8bcb3c87f3868c5f812da71aeb21c5ce18a508773120ea4f659 + 812eacdc5a8eb2e8fde0c23a19c2e33466d754ab3109d28ff48fb41dab8066e2 + 4d42205e20486ba5b1a942738187dc474049035e4e65b31a45172830b9245d3c; + 598fb4096a7054e317e66bc9633d56dfcaeabc31a7bbc9da13074b242edeca09 + 929b1bddd55c99b06b99a3ffbac5adc1d28823cb8354adcfea60b7a87fe32039 + eb2ad0e63fcded93837f0fc91e0304a19d73e0fc2a1077aafe6702cdadc1eb42; + c42183d915267e655409326b1b77a042d8918feefd75dd3d18551f398ef426e7 + e3729a6110686ffa57d1c69284a60fb36173148a89b1810d04614f03c8f90084 + cd941d3b268eed5facdaf8fd9f1db0f53905a47887275f4b1cb66e3c56ee276b; + 67b8b7ea7e7f9f52d686f906589c919b293bff6d3127e7ac581dc1630af7bcc1 + 99d18d2190ef0da148b843c014fea3efea59f9adbb28fa38407191fcc24e029c + 268a450c0f6fadf31e3f3dc76c9a358b1495f81bed4fe1e5988e5260cd45bf5d; + bc3422b79e01234fefd3637a19bfb2dd0a83bf2071b27bf5389918a50e579294 + b665dffe064aceffaabba94697cf27054d016427ebaf68858c5250cc5c1b770f + 859a01b6a54bf14e9a8f0dc1b08edae2578423485c62e47ac5eb68716b720924; + d419cbfd5514d9cffa1765b9f273fd7492000a3d110586ed0b8137760a2e808e + bcb52f423917c074464b2078d1c3f883fb793b02411b8b13c08266d17658a499 + b6cffa3f8f2b994441638531c437f6f88d7a453f52201101cc039e4781862428; + 2f7eb8051d16a52b797e2618881bd0857eb205bf97e581e09680f86308f7e082 + 9ceda63b29c3bfedbb6f3f0be2f3f2bd285bb9d21620c3edde69629b00b1c3cd + f16b5f4146d9641935ee65236a0fc343a70dbf91ae0545ce75ea5aff08a8a450; +} + +sub { + dfb53b4ff67c742728c564d50ae20e3e3617804478a95e511d37ff4dcf04fa3f + 9bc2006bb07dcf7bf0b461e7599e3c3cbfcdf6b6f41fc1db2219767ce9c43bcd + 1ef33ae445ffa4ab371003eeb043d2017749898d83899d75fa1d89d1e53fbe72; + 432ec65076c94af41e054a43f64fe4eecb23f50a303fd38e0172e4e846e67394 + 3642c84c5f5fcd5a0f3f2e37114beb8e7b614c39d354f1b32fe7e9cc918da5c8 + faebfd03176a7d990fc61b0ce504f95f50c2a8d15ceae1dad18afa1bb558ce4b; + 1725a986c03d9d53346ab2064545ba90ab38191e19ded9ecc2f44d422ab04150 + 11a20caef0239e4ffb5ab2f7d3576633a0d4a08c78a02f35765fafff5398c48b + f3829cd8cf19ff03390f000f71ed535d0b647891a03daab74c959e42d6177d44; + ec5d5cdfc8183f9aa3d47634d16427c3a399e2166bfa4ddfecbc911e5f33ceba + 649f8e29c483765c4232c369fd0e50a0026eac306eb4bd6924e7b4d7d3f73e26 + 9bbecdb50495c83d61a2b3cad355d722a12b36e6fc459075c8d5dc468b3b8f14; + cc9927c331884612d8fab9a2262ad0db38f55e0cc587e45b972ea8b1f504450d + 30c5c9f8d5491a350051bbc6bb78e568d98061b43cb1797cd3d7454e8a7aa1bf + 76d45dca5b3e2cddd7a9fedb6ab1ea725f74fd5788d66adfc35662636b8aa34d; + 662bc96da214e50f13e07a15b8e0b78f512c744a66975b2170ba28079ca7d03b + 84b072295dca0abef806ec5c1e594179dc5e9eb172c178af1902a42ea7a3728e + cf7a5644454ada511ad98eb89987761675cdd598f3d5e27156b884d8f4035e2d; + ed6f947157492eccc5fb1d0542593d58a47e2d87355dc0378e65376163fcf457 + cb27e06c5ebf42a61caadd4fb222799772aaef31487bf61ac84c15f58af265d3 + 0f48b404f989eb25a95140b58f36c4c031d43d55ede1c91cc618226cd8098f04; + 33c2fc57d8ffc29e6d670cfc4e25b695320055cf443ade53f49f6c23c452f644 + 3814860ccbb4391894c107961170d3af1594988a7972b97d0fb2169631b216fa + d5ad764b0d4b8986d9a504663db5e2e51c6cbc44cbc724d6e4ed558d92a0df4a; + d0fd0d5e73be5c534aa37645aecc49b202fadb30733f16d0de719a70bec2c751 + 910c7b4d317397cff70c91e5bf8ec42ba1700b0866adbbfe99b7fbb2c46c7093 + 2cf19210424bc5835296e55fee3d85866189d0280d925ad144ba9ebdf955573e; + e0f8130d6b0125e824b254ed57edc59455a8a4704a3f14a4ae8dbf8721595468 + 4e4f69af6f85742068ad70225e8850c85aef3d3fb3ba3b52c2be979db4eb979d + 7fa9aa5dfb7bb0c7bc04e4caf96475ccfab866319784d851ecce27ea6c6dbc4a; + 204d5f69ebab4597d2aa9295ecc2959480a2ab0b225674cc8a147cf7c1e59717 + 3abdd29c94819ec7507673bd5cc4fb267fcb17c1509abd5e462b46b012f815a8 + c08f8ccc562aa7cf81341fd88ffe996d01d7934ad1bbb66d44e93547afed816f; + e2d75c74060aee615b6cfc815c8618f826b07c4987e10f20d27ff74eab3ac7d0 + 7a06b32d284b5faedb91dac8e94cb7f76c678bad81ed78ac4e75eb9edd8d11e9 + 55d1a946debe8eb37fda21b972396100ba48f19b05f49673830a0cb0cdacb567; + 8d50b0573b0bb945bf128e1a9a37c151e991581673cf6295db5a6223d428ac9d + 47ca477794fdf378f3486355b5348960ee957470b8c74d62199ac48b4478c4ca + 338668e0a60dc5cccbc92ac5e40238f1fafbe3a5ba071533c2c09d978fb0e752; + b66b49190a99456aaac5c33beaec423a3fa79a84df0f00122aed29151185ceef + 13c06f7941743bfa261845d99772fd0d5f1a8258aeda53b161505ed8bdca1bc8 + a3abd99fc8240a7083ad7e62527a452ce08c182c3135ac60c89ccb3c53bab227; + 44ac89bb8d4de1a8756b00401236ea4e316a2fe637f6e5f052576da70c29c4ea + 05dd39ac32768d899b80ed10279aa10d6316d60ab2d061e877dc5f21ce8e1140 + 52cf4f0f5bd7531fdaea122feb9b4841ce5359db85258408db7a0d863e9ab22a; + ef0f6c7f56ff22260c66d583f79c7f581181d3a0f3167c93786bb6c89de4372c + be934afeffaf5e6755a92e6606ffd054627dbe1b18ad6ae543363a01da605512 + 317c2181564fc4beb6bca61df19dae03af031585db6911ae34357cc7c383e219; + 773e0db8f86c5182b5f51d10376f42cde275dac6a646ab1cab065760e045107e + 3f92ddafdfe67d357b4f9b61807bcb8b0f1f6c2dacba27902f48fc15a31db606 + 38ac2f081986d34c3aa682aeb6f37641d3566e99fa8b838c7bbe5a4a3d285a77; + 34860a848d0afb20e43f4a36b3c8941f5821864216da9cb3ce1a11f056565a3d + 9dfed86048acd45ec0009c1ba3b20a1ba645ec141070fd8ae41d25e159bbce9f + 84873123455e26c2233fae1a10168a04b2db992d066a9f28eafceb0efd9a8b1d; + 854ba91672ae1b8f6c29395dcfa4f3514ed97dc2c9e1e419f926825619fac55b + 7a700653ad355caf17f89045ecfb3ee821194b54048e55d8eed41fea2b0a22ab + f8daa2c3c478bfdf5431a817e3a8b4692cc0326ec5538f410a52626cedefa330; + 0c66b3d8c3916940b70d381d232ba8aecfe27e3ae63fa6cbb287216708f224f0 + 477186656245c4083dff93a95c6e5deb707916e52cd74abbd522743c24c53088 + c5f42c73614ca5377a0ea473c6bc4ac35e696855b9685b10dd64ad2ae42cf467; +} + +condswap { + 14e55571df646a69a8280bf8dbf8e9afc15bf5558bb8b8236ebcaa19a96053bd + 8bef5021598a8175565e1f2b522ed1c4306fef4e0e973b50d3a03db1fcf11a43 + 0xffffffff + 8bef5021598a8175565e1f2b522ed1c4306fef4e0e973b50d3a03db1fcf11a43 + 27e55571df646a69a8280bf8dbf8e9afc15bf5558bb8b8236ebcaa19a960533d; + 53d8c54864e18f8fc742b2cb996f1d7595122ca9c90de9f49485cf5a5ef05058 + 9aa95be655b6704a3dd482e8424f82d17443ae26ab41092fad95df9cac678787 + 0x00000000 + 53d8c54864e18f8fc742b2cb996f1d7595122ca9c90de9f49485cf5a5ef05058 + ada95be655b6704a3dd482e8424f82d17443ae26ab41092fad95df9cac678707; + 8c7b4ffbf50edf7edcf8d69b67796de38b56a5d4e3a9fb2201354cca0b1b7c57 + 856e707988c37e59fc6689d024b98bc01b71034484ba194ec7f45106267cb7f6 + 0xffffffff + 986e707988c37e59fc6689d024b98bc01b71034484ba194ec7f45106267cb776 + 8c7b4ffbf50edf7edcf8d69b67796de38b56a5d4e3a9fb2201354cca0b1b7c57; + c8cb24a636f6d8f2c045c8441b21314884a71d55c73d4ef7ad65620ac0d86616 + 02b47be53df7c8521d57589256e69b8bbeb1e838d84c090287c04a25d5c08e1a + 0x00000000 + c8cb24a636f6d8f2c045c8441b21314884a71d55c73d4ef7ad65620ac0d86616 + 02b47be53df7c8521d57589256e69b8bbeb1e838d84c090287c04a25d5c08e1a; + 7d96f47612fb7135edfa71dd4526d8b4c944d74ddf004e1653d3af52485168b4 + a4e83162602d20cf71a279c1071038ecf94f38be11ad1191927b3f480d393e1d + 0x00000000 + 9096f47612fb7135edfa71dd4526d8b4c944d74ddf004e1653d3af5248516834 + a4e83162602d20cf71a279c1071038ecf94f38be11ad1191927b3f480d393e1d; + c3a3d4462091b36bd6fd494e70c3166478be0a65bca7f7361246c52dda402981 + 62d0ae1abbe5cb909328f59e71330269e4f2b1d3c664e7f1f0780e0a4774c262 + 0xffffffff + 62d0ae1abbe5cb909328f59e71330269e4f2b1d3c664e7f1f0780e0a4774c262 + d6a3d4462091b36bd6fd494e70c3166478be0a65bca7f7361246c52dda402901; + 38c9b0d92be50270e29cc3181952e745cbf21c8d145ae0c0a6e59f25394ea59c + caf4896a2d3666200741d817ee1dac661652b9a1741ea966d5ed7cbc4c76217a + 0x00000000 + 4bc9b0d92be50270e29cc3181952e745cbf21c8d145ae0c0a6e59f25394ea51c + caf4896a2d3666200741d817ee1dac661652b9a1741ea966d5ed7cbc4c76217a; + d637a72affa4112bbccd2878208ce34fae442e97175e1b39f43c00b14d312a2b + 76209f8b1f152ba7e1e6294f707831021eb2ac03872ed3774ab8fc5e6a0e864f + 0x00000000 + d637a72affa4112bbccd2878208ce34fae442e97175e1b39f43c00b14d312a2b + 76209f8b1f152ba7e1e6294f707831021eb2ac03872ed3774ab8fc5e6a0e864f; + 12f3bd417fb3860b2a4fd04eb1c51c558a4c5a3dcc1ac457654a45a3cd45b0d4 + 46f802836c7bde099775f3fb4f5f6345364a7e158bae9e16c4c99aee0cfa5e42 + 0xffffffff + 46f802836c7bde099775f3fb4f5f6345364a7e158bae9e16c4c99aee0cfa5e42 + 25f3bd417fb3860b2a4fd04eb1c51c558a4c5a3dcc1ac457654a45a3cd45b054; + bf1809201f83a34cb07aa2e2372516631ba0513e600956e6702903e9084770c4 + 414bbc77d13d7b4e706acc07406e8d624f2485463f9948d81b72268c1e086a14 + 0xffffffff + 414bbc77d13d7b4e706acc07406e8d624f2485463f9948d81b72268c1e086a14 + d21809201f83a34cb07aa2e2372516631ba0513e600956e6702903e908477044; + bb8703f2ddaaaaa6ad19219684b1bde576b38bbe9065178b24c1bc55563fe525 + 1cb3ca7846a3c584b1bc05b25ee91d4779ab9ac64ecef0fcbaea8d311b55d618 + 0xffffffff + 1cb3ca7846a3c584b1bc05b25ee91d4779ab9ac64ecef0fcbaea8d311b55d618 + bb8703f2ddaaaaa6ad19219684b1bde576b38bbe9065178b24c1bc55563fe525; + a07c7110c32e362864a42cd2f371ff420bfd442d291cc15ec079d642b5c85bdf + 497fa7867fcd0617c4cd765aa6f46b89390744b5b57d11ab732153e075fcb607 + 0xffffffff + 497fa7867fcd0617c4cd765aa6f46b89390744b5b57d11ab732153e075fcb607 + b37c7110c32e362864a42cd2f371ff420bfd442d291cc15ec079d642b5c85b5f; + 10fddb7f48015b394512257da026938b67f65dd0a84dace51417d68da003a913 + ebd6706b7cd2ee2d97848ee1adaa990d8ddbfad48e43f51d7ef843fa7597410c + 0xffffffff + ebd6706b7cd2ee2d97848ee1adaa990d8ddbfad48e43f51d7ef843fa7597410c + 10fddb7f48015b394512257da026938b67f65dd0a84dace51417d68da003a913; + 1dd835e2173b936a91da37b3aa11e848e4497964a9a78ea929a19105eb981a24 + 3e6f1ab7d5c460dfb7ad1926c60b64a8dd48ed0115b0655a6d8619666b8c8dde + 0xffffffff + 516f1ab7d5c460dfb7ad1926c60b64a8dd48ed0115b0655a6d8619666b8c8d5e + 1dd835e2173b936a91da37b3aa11e848e4497964a9a78ea929a19105eb981a24; + c4541b1380e796a333d88affdccb4d2bc68bc5a1b3890ef2fc1a2c38dbcac725 + 6f6338ad8dc0562bc02e307743b81013ffff6c135fdf8603f9956b7f94eb55fa + 0xffffffff + 826338ad8dc0562bc02e307743b81013ffff6c135fdf8603f9956b7f94eb557a + c4541b1380e796a333d88affdccb4d2bc68bc5a1b3890ef2fc1a2c38dbcac725; + d71f16494f041d1547e5dd9ed539e6b9c1506cebb66f6424e6c7aaf6d0ace080 + 78c1e6cc0b6f64e161def3b2ede7465f21f6afe61299e9a3c52b1a7a080b9b26 + 0x00000000 + ea1f16494f041d1547e5dd9ed539e6b9c1506cebb66f6424e6c7aaf6d0ace000 + 78c1e6cc0b6f64e161def3b2ede7465f21f6afe61299e9a3c52b1a7a080b9b26; + 4d90838ab9cc0a79f4e2aa6a860bf8cfbce5f834aab428d4a06f8ef4a2da582a + 53124b302e7675ea5784e2e7d9fdecb38d9c8be158312dd81c51757bafb79b42 + 0x00000000 + 4d90838ab9cc0a79f4e2aa6a860bf8cfbce5f834aab428d4a06f8ef4a2da582a + 53124b302e7675ea5784e2e7d9fdecb38d9c8be158312dd81c51757bafb79b42; + c8e764092ea695b0dd8f0ecbca5a4e02096c821b0f5cd57de8e3e50d8b2a40db + db7fb38dfede5259a584dbcc697f259ce380110636e301acb308bf687449110e + 0x00000000 + dbe764092ea695b0dd8f0ecbca5a4e02096c821b0f5cd57de8e3e50d8b2a405b + db7fb38dfede5259a584dbcc697f259ce380110636e301acb308bf687449110e; + 34cf8e55b3d2494362a87a24907fbdab61cb4f452ef5a889a1aa40a22be6d976 + 679c6d6f6dca1b71bb098e3854c63ffb8ddb76f1a1cb246ce956fd71d6477c20 + 0x00000000 + 34cf8e55b3d2494362a87a24907fbdab61cb4f452ef5a889a1aa40a22be6d976 + 679c6d6f6dca1b71bb098e3854c63ffb8ddb76f1a1cb246ce956fd71d6477c20; + 9638ab322ced065068f98597cc61fc2bd846b5849dc39a881209b8efcbc95eec + 32922ff2b62e968417225b765f2787387c8a42fcab1dbec3815298f53813a8cd + 0xffffffff + 45922ff2b62e968417225b765f2787387c8a42fcab1dbec3815298f53813a84d + a938ab322ced065068f98597cc61fc2bd846b5849dc39a881209b8efcbc95e6c; +} + +mulconst { + 8d517a8804d934ae93388e1d59748c3c22866d7be6cec58a357ec43f6dd029db -319312 + 3f2fea85c9195f2246f4f6fd2e8d5de95c5065edbea45a174bc0d1ef3db4aa54; + 1b5d6a6c9cdd0dabd126e5f5ee2513f3d27c56dd3aad89c882fe9fd2d38a93b9 499342 + dd310555966ebe7f409dcce2933f91af8e14c8909f78e4f113606789dd7f504b; + 9e03b44e93743e0cf5768f65bd3c0f255e001de330d5e1a863d6e4dbb11748a9 -284189 + 6e8c7f8d0f1476cf2603f5702fe8e8bba5c01e86b6e008ac72987fdd6721183c; + 9cc65e302736aab9bf02dd451c216ac3dd2066af299d7f6bb361eb3b36777f74 -15738 + 8c067e5dccdafaf56301050be17b8895c07e7f13572bbc5ce1af8b5be043d91b; + 7f36ff91ed165c6ee44aca261dc09e3e124c2ed8167064c20d9b84fa226e2555 423857 + dc7415db1ab5e5a5658d089ced49bf97718bccc55e2f15b83a28157cbe47e559; + e8852c80db69c7e55d476ed00ca7c0735052f7f72f81198a57236feb89a5903c -140934 + 4923ee40bb01314c2fe9f2f87d0743823311b2d9265219eb7a6de5fdfd09a308; + 505ef65d294b77e260757bc867199c4b0922b5022dfd6a1fe7478ea27ec102c5 -307893 + b14d72c991a0bfc22b1b95934bbe589b6a4a1f9c76a89d7012f8cb03379a8844; + e410c260dc8c1838f3e9658ae586425039b373ee2ab88c750a51853b012af0af 471420 + 42ae337a920ee1b83768cbf0750f524cf1f9e98b99f195e869785f90a2356756; + e2d181914b54175190462d11402717f42efa37967b2fe0bc4433fd4c862fd646 492178 + 4c6e1ba6896d2d661fbf99b77f63e3095d3341f04bf258dc3e4f51e39c68f576; + 9cdda0d7da106ef5cfb34774f48308b64953e92869dab0ba804927080509d8a2 -312149 + 4e75e56793563268480efaf6232ad702bc5aef6246bd44e3499ff437810b523b; + 9dac4442858b8dcc324464911e6f23ed046e42612e8105dc3aa2b9eccea15dc8 -88818 + ce52a0897600a05383df0a0b99a451f0d5568146b4323392fc164153715ecb0a; + c9f92ceffa7027ffa71e64b714c36b28a6a18c3608f8170d11d51acbee3c62ba -259755 + 89582517069925572219b2eec2a25029e8e1c78772ca8c50bf23220d02740763; + 8f8bcf28f4920465680c8e15befba4ebd23c1e27b71c224c98461f69caf76eb1 -430350 + e96b93f0aab6ca70f7d2c1121110efc320aff34e1ec693d7bf17bce1a756d622; + 6f4d413a27388fe84a2585e809927bd470ba1b0dc749ae6e2dfa100617a5a7db 11486 + 1cf1f5bd0770d949ee36df858e5842805d121722de30d2ef6bc2b82d021fc14b; + e1d05a072d1ab70d895969271f44b7c9d1654171836cf234724d2f90434c04e6 8057 + e32b3c792cd3b6a770eaa2627ff8de8895873d73d734b961646fbce0643a443d; + d8564ee3ff02c80f0bf3e9cdfce42cd8ac1012ad8ea74f8c0f5754b5775a8de5 -304373 + 884266ea43a6a9b72d113b437dae63c664ef53884c893212fd864412f4e5e457; + bf66769ee0cf873764524e760ee654e8d5090aa774f2e96144d945f8da37af19 -159091 + d7438999b3bdb99a46e2f3dc8492cbaf52b0a19c651cc382c384940861bccb6b; + 72b7d02b1fe23adea87696864c0a9f0466db85e39db397f0eb428a50e277aa2a -521654 + a611bd153336e178a7f3d912adc2fb39b342132458faeedbdb0a35973efae23e; + ff5c39d4b7fff86810a7a3b14a58bc01e6ecf062116adaa78dac0d1bd3f6039b -140759 + eb9e32bf7d027cda198048da50aa79a11bc99f3b37a219a2104c78f4affb584f; + a6f1f7f30e3f9e5467f180f92d2c94165fd9427d0c7d383b9077ca9a7e50dac5 -317637 + 56a32847c50bc9609db9796e1650fbe97524aacea746d3cf285f0f2fced01a49; +} + +mul { + ## Easy multiplication. + 0500000c00003800006000004000000009000006000048000090000040000000 + 070000180000200000e00000800200000400000c000030000090000080020000 + de200094610090df0020fe0180270700650d00dc2a00f84100e0340040570000; + + ## Random tests. + adcc6f10734dae2273304a6aa493ee8f96e05f2402341c0d997dce58ea57ff6f + df3c8d5bb6380fa278d4ff2994c7865bcc146596d3c3126242f0dc3509b57449 + 5ee2615b75a40367ced4ed4bf1578d13e7afa4d415f38f802b7e7b374365b627; + 054bc4d5e9b78b068fd20645eeb1f03aef2e89a1f56cb50e5f1170a86529526a + f976ffa1d1b4a33d0fd866528897cea3eaffdc75e31aa65450a62ff765fe7985 + 0210280745f6c29a04433f07f880bfd2eddd758d82996e79c65de817c30e7309; + e361c6fedbbe081b7cd683d23a2bf0ac889806be7230b56d9398959713a300ec + c8121e53d78c2e9ba2ba7f51b7cb15cec970a63b4a642ac470aa1f2145f07bd0 + b48bb1c7ddd26cd9e5a2cb73bb3885dbe4a47e03886d4dfd59e5729d544d566a; + 463158a2077f931567d45e19deb4454a2ae77045db70a2c078e160ebdb74fd94 + 3f311e3f8d225016448fbce3fbfbe84002736e3b0d5a90f57d705ede8706008f + 44cd97e25c20ae2fe2ef4f2630981a1fe54c79ebd904bc40e39d9fa615ffa76c; + 1c73fca52556e01204d4b957a01d1049f2247c8666d90c1d3703cff16a38c8bf + e9afa0cc85577ec7eb1ac66d87260d5f3659936ff88cde8bd1c3fa1d499e8bcb + 0979a8baf734c86a4f57ae400b7790c4985fba3088cb1b265dfc7a66c794121b; + 6e2bab74e427b1a44f0cb181f9dedd6fed43cfd395d3708c4ae8fd2137215eb6 + ef400834f2dad8763a1a5a58a37d1b36a78873860e4dac808d3166f13724942c + 0d481adbce7bd82ebc5a3e5ed7a492e78a032ebf211a5bd2371c528c77954936; + cb3575b41d9521da5c85c6438af14903a8d6d9c3857aa8e101e5295fb277de7b + 9f76715d2ae3c8ca182406553d3481fc36d67727ca2d51f5db37940aa3208d8d + 0b2fa619265bd424b9e20d38bf4552ad63de8ab552166080d099edea5fbadf7f; + abd66cf3b20dc4c5a866f1a3b965f794a726f96dd00f576d16ead66eea3a1177 + 97f70fc92e0cda0e010f2a9d7608c33db954af18bc18ad55bf66e7d91ae31abd + dcbb2cf34274971ef030827865b093a00761de3da90e46f00e18d0fd934a2b7f; + 3b3c96ff9be8e17e0a0d3979ea5e5623fff5fa18914151abf79fe66639beef06 + 1d29c361677805f1dd2b839e025e56208a10c0b8d9f2a8c490b0a892c0e27080 + 0050873016bc67023d0c479fffa56a5a70512f9864d00b0a3e951f1e85c8b412; + b263d4a71b6d963e0a972819e5e5805b9522f142d3c0ff2738c51cce31e9f0d5 + 64f3fbcff2195c297292df3e1cf900aae457e09864c5a6248ea21ef8efd557d9 + 7e6f6f511cbf2743c0fc1c6431243908606cb257d20f8a8a0885e6ecd000e811; + 0df0497e1e01dbcc4a48d4865127d9cc31a9a0a56e90cc31619e65c1521b2369 + 9277f0c6a5bcd7c7ed8dcf7182795f67a12e856b126d09966689834c43aada47 + 798d270ae34aa2452694b4ff74049a6a5169d71c85bd1a9ac59625e9298cdf1b; + 1dc64bef2189eb5349f5a8983c1854a40141a5390ea4bf8535a15aa71d0289b7 + fb59dc9b08f9163881b10020ce040da296ed8c1a9a935c27ca0047a4d7f6342b + 9fc23104cefd3014a469e6d3a0631338eaf91e24994ecd50bfc0b4547317317e; + 975b108909e87803af98af9425d2a982865b6cc9c144985b32fe5bb842702e24 + a69e9130dc40be50e827dc50dac05bf900946cf67ec35f1bb1705001314e4fcb + 12d042e50d0cff7b43a5ede61ce827a99619bf84995d488fe9eb56a8de769164; + 39c0db5a4b161c368b0f71c7e8ef97dddc444e4a33181e0743519f2a879ae04f + d3cc45ff509ed45998bd5f218408ad51f7457b35e8ad69d5431f174047d7c6e5 + 96b877604f52782ff2a5b1a12d7db642371ca993b648c4971cbf593213503219; + f8c6715e4a327cbbdbb08d74121619abc73dc421d7267029e39b50c1510d2045 + 63870a83cad4c1c64195ae7b92f3b19c830f35f1a3186815565d846c787c7ffc + 286a0b8c6681a7b79a5f85426fb60819b7bd349e93c1e11943c41c8dfd57d00f; + 49c9bcaf5434b8205a348df290a96da27df084f34047efbf91cf544057a50bde + 54269f8b23c05b94c9c1cd74eae88d9227af12a6fb337548eb93c2a88a75a8a5 + 1e3dcab5d7361f00e70198b52e8fb6d1fba74abd79d99480527f6711b21b833b; + c7212b43ff0d10fc52d3a2a7e15bb9cf53856f7ce33b79df89e1d13f95a94511 + 5213ec20bd054491a5485f665075835da609ad4f3dc005744f061480862ee602 + e136339fa01067c6c0823a7ebe2f84d917600e244a9883191900bb95a3b59764; + a15d62c6308448f1f70602f43be1825696d228eb170bed55322a21d2b1afa915 + 58afe65bc8adbe22807026aa1b0d68af330bca268fdc054406aa41c820960e71 + 996c597e0a024c2ff0ffe289d8c7270e4ea1d38e294d7210eddeb2c0725a1305; + ea233ae033ca93034bc8fd8ea3c0625bbb15b00265264fdd61f62cc808b505fe + 09e12ef9b3426991f7119c2e8349278342ed0072317c4ec7ab24561b61396b39 + 1fd72191fe9f00b0cca6a635d14904bd0c1e0ccd1275e84ad0b7bde134866e2a; + fe1fd3da05524af4ccd538f1d2ede09582101fc2a47c2e63758a546bed6cc5ec + 523346406ef3df6ac9a455d0d2e7f732ad01c88245b07a888cc21c40fce11a76 + de6cac8fdc878d9c4eda6b78f8b629ecf776980e49fed56bab936fc03924f44e; +} + +sqr { + 4dc87ad85abe8938634492ce2ea10bf48dc4aa2fa37a6a0fde3e5a01e90dde2b + 346297c756d77b8ea312ddc455223f4577aea2a2e8d944cb24adb487de143b46; + a33aed91165d7fd13a543760d3915c1b5d411a21895cde81b2f3640b87550cc9 + f923e49da9d5f45e23a782f911fac1a6c2e67cbbc7b03ee4675a0c5fbec03641; + a1406e453db6c92ef808295d1919c5bb2491c0ce6d5326d276c57e57943e0acc + 18425bbd92149991c4545e65495af0e6372289e9a32649c127122243204bdd73; + 70a56b8484234ed002a68c8d535a4e527a61c360f7d4ee5b4970c8ace1e12de0 + bed1c435139f151d0650577e03c32428f63cdb7bc1920110b5fdbb27f08cb200; + 86c3236a6194f1585b4dbeb09d0247fb7d3ae08063ccf9fc4fde70fc23735710 + 9b6509acadd5d31c4c843e91cec78be00aec4fdf7d86132965bd6f62f078f711; + 0e5e83c51a8c091673425c38785e0aa1d74c7261a035a891bf734cca8c8a3db3 + 4efe5029e99438835a8a45d0180b8de357eac82ba9dd8e127baac0d758dfca1d; + 31c39a5afc33f8a8d3c570b9dc5389bea6fdaaeeb44cd7524f6c445b9583a34c + 359b7c8791cc87aa8b3c799c9d4d7e1f4062a4f64436d8c01071597e8c262127; + 5a06203c15d6c1db6ea0fcf468a1e3f09d68c4a0f64c3d7e30e70290545e708e + 07f8b435ea322bedf613458b1eb420f1d6bc96741d69695899bd40f7eca6f72a; + b1973433646474167b791b326827641260fb33d6a2c96d4e963a6b091f824505 + 52a62396dbbceead62912f8adfed2707fc25edb51aac19b2a274c246f1209f26; + 96b157d71e107168fac206dd63007f10f4ae1225b0761c470ee18ecdb4fff20b + 84ba2a84b156437cdcc39ecc8b1ad01bba934d8dfb7d504b9f2442db0ad01167; + 6167d77675e19785960c28df551b458e41dfd0dfe5e56e785b7367997d5eb0a6 + 79a558436af4310d3eaee93f9a7cab72bd817fdf59d267e2c530ba8d667a005c; + fb4b2ab08db932ce10dd62639aabbb266f9ce775fb5a5e0875864f18a3ee506c + a6013661f2f0ab875739836039bff6bff4c993aefdda43dccecf4c404b9a612a; + ea26bcf0bfc3e15d4df7ddbcac8b5b9e9a7215d74196d77c9e0524e861c6aeb4 + 113e06861e3d338618695f4713464cea536c8019da77ff0dd9946f3b860df062; + b67b1ffffec757ef147e86cf60898c5dac515f080f43f4e120017caab4eee5ee + 7d9748a1ff6ca9cdebeb3dd6668522767afb9ebf1c680fe7fb3f5848f7222036; + c122016f17f12e53f9ac031719fe428c826217b2b757ae44f2f33d1b8b1bb553 + dc068d032b97c98eea6dd6c8bc42e9a0a6250c81a36a8f64615489e191287c72; + 1266ff0654c08e6183cd6c3cc760cc5c88d1f49b918f2f8f587b3edd4adea513 + f07f17fadb2e2dbd282d296fe5a15cc5c2baf6b45f5fd5dc4692f252039bc862; + 4c4dd9c6e46fb2c096388583f2a0c969245452b23f16ea43e15064c45357ccef + 5127fd4415a84115997b8b3ebe88a590efde9f8b53a43e86b1a4d60782c91025; + 095d149c2616212e5f770f133e719d285a8bc888414271296684e732feb5473f + d8038bee551cac14d1c610bbc0960d943b01526072c589bc913ab51fd47bc848; + ad55f9f916d62bb13ceb1461b99e82d76b860c084ca408dd4a18f154b907754a + e5c86ee6350c9ab736229f3f53e2ab852641bb8a62fdda609444df377e9ca04a; + 2e4976f45acc671e38add93aa6a89669cdc9874975879b3134dd37c633c9eb14 + c6513f1a11346d33b05ad65eca13f11125cf70bdcdd6e666a8a7b9c5b6098429; +} + +inv { + 344ff33a2bc801f861bd13575797da35352c8164808d27f51ea7ca46931b2b1f + 07f678dfeef071b038d692bf825cb4e75bc1934bd3a6e3dfea99e78bdaeecf3a; + 670097e2fbbecf6516e949d98f6d722d3f732c33bfb4e62512f9d28d94546b64 + b2f5cbcd7d8c996ed33a813b1aa39c36f5826b436867bd89b4ce43a156c1e96b; + 234f3bd8760cdc8ae7fac9457da418cd082e42906ffe59a13e52185a7845e61c + 7488ff9231e4819684148e16ebc25440732ffe7d211eb3d3263db1af9f1d0513; + 22a39a35b1cabc6dfecd1b898593070f0c81bbd04b4b20f43a9c7c14b5dfe967 + 892024fdbc8937b9b3aef794f57b8ed8517510b6de3f6088987b1310c1dc7932; + d850c118c69198268826d05650958c6bed55b181090cf1c96493ecc4635a913f + 7aeae1b0dc5211f374872e78da8f225325f3f33f9c42bfdeff25d425dff7710f; + 33ccdc0e39f7cd0a61782c87ac089eb7bbcabce4dd222f05088d3d6332efe959 + 668f87e23239e965c80740d4daa9c913609388f291c11f950684f45e689a061d; + 7c06f93b2f66001c74073b8ba7536049738e89b0a0d859a8a73710de3fadcf79 + 52df0aa24bb4a7474e6ed9852ba542c1f4c5219372a7b5d60988890e2955491e; + 9a0ce9ff4026bd58abcf8a128dbc55a19c446005834a143634452821d5f9f0ee + 79eaa93b6287705c0ba564c31bab1854c7f4354845be065c4bf94d38a7a8b421; + a13a127490bd3443fdba99baf61d113f3aadd6966c9732cd16348f9fbed32e0e + 67290b25e6bd7328dd804309a6203c02a75657f593c7e66e48c9192524e5b87b; + 35f7ac3eaa36c4ed338388fe7eda61cf237eff6ddc0847cb0c95255bbaef8ddc + 4950601d2cfe75eb9c4524e2f40f404fe82d0014c916dc43d9416cd2e274992d; + 8a058593035986da3a19e63dacaf40fb3db4cf382205fc4e4802e0ff0c61c8d6 + 6a23b743dc2fa4a4132780cf6edfed4c219760f3b31343b08710ea18d2c0fc43; + a10d03234e700ec92270b6105fdb787d8edd3f960fdd78ea2e169f49dc708d45 + 101dad1720604d0b093831baaa97fd0701ce773b9c1be1db0403533ad4376c0a; + 52d3b9e030a5ed574f7f9489c47e73f78ef88a39b570045f9ef73c9a3e65f425 + 83da3e088dbc9a3cbe06496c5a52edc8d7abcbb8a0f839361762137b119d961a; + 9928b5f14ec80a036479d7b08c128e601db880daebef36adde6fbd1f50f82629 + 8dcaa9a2439b77cbe44928ac2877e09755c505f58df4d0ffb9b23fc256b12779; + 5a956523956f8ace6166ea6cb3c7421f8d1142d26084a81d615f00d578e03448 + 658a3ee24ec5b5f1d2d12614e535247c9f2f259131b00f654e09b4517c63ea09; + 157492cd38a07802351ca0d1a5615d5c5f254692acce930fb6be351e9a52f088 + 32fec1eab85ddf8d27a437876e7f5a07f012b02e7971869811c3c31230204164; + 2d07c109c332847567621ef4338e149b99c5bc3bf6e50943695da7fed5d8141e + a81472d28a8235bff43ad5cc012d52e95ece3dbc5ba2962dc54cf9e7f81bee04; + 7e172f6c03415e83d4671e0191151d49e5257651ba73848254e68fad0a095c50 + 1a8efdc57ae9d39ffa64afbc5d2dc1f25a29d365ac3dee0d8529977c14aa4415; + 1934546e37b0bac80319fa7f2f86c18e83de8da830b9d0b5019c26768974eda4 + ffe55479771f20518d8fc6491055368b5dcce9b4e193d19897f1e116b9cdf946; + 90324ccfe0e34e9b08db1862be3bc17798a8581ee097be76545e3bf9d586a2e6 + 0f5ccf4920d617c9b4f7fe6e2c9e95e89b5516555d5e3f5f9900485451d7425f; +} + +sub-mulc-add-sub-mul { + ## A very basic smoke-test. + 03 01 5 01 07 02 37; + + ## Stress the bounds. + 000000feffff030000f0ffffdfffff7f000000ffffff03000008ffff1f000040 00 + 0xfffff + 000000feffff030000f0ffffdfffff7f000000ffffff03000008ffff1f000040 + 000000feffff030000f0ffffdfffff7f000000ffffff03000008ffff1f000040 00 + 00000cf2ffdf4400a09a000079fdffea040098edffa73b024072ffbfe0feff65; + + ## Random tests. + 6f0dc2547ad77f7b3da324876e63c724c1f69cb2bc70ff960f1594f90202e17f + 23a2d45c26880166516e6831ab5070dd2247992418ad7cc41155a9362744f763 + 4783 + e145b5f91b683918d1efe6193de27909bf9051b0157acfa458bc3f63e5272de5 + e694fc475a442bd72cbbf11206474e3169fa74ad453aec3034c1dcc1f1dfbce0 + f9b4f1ac9369241a302f874305cd2d2fba362b5735d36d3e072c7198d7dd4005 + 08518c9f52236272202869bf7cb549c667f6a27d1c4f18ab58605d63821a6c61; + f1b3d663fb9458ba8e4521e8e2c976fafc8936433b62e316c6dcff794b1ac8d0 + 2185dd6948ef2339c9a43518dd02eac6518852d1c35858dbdf5c1acf1d995884 + -128133 + d7a628aa59d409b755a1ff0d76f37ae8e4a73d47375940c61351074f3c42ad2a + 4a353e9e7e9a0700ab2c043cc05f9795cc267e001b64799ca1849f024a42dfd7 + dbef6513fba5ca19c545b162718eeba960b0930b8ce99c3edac5facf5fddf4f2 + 9dbfe69897a273171ab4a697f88a0d4b420a3ada7a7f9e67d0ba871e2097026d; + e8624dbe04a1b19102019d7fd9e06fe502764afb92330d03693279e872e84f81 + 8da8b04024ff708a8b5fa31a6a654fc31d3748f74a897898aa42c2bb65494586 + -118148 + 83f63df59ce4b3161427b932ecbd7796bd162b553b70bbfef8a990ba6f20c23e + 9d896d69ec478456fb6cad67774c6e3e0e313aadc0564d759380981fb12dc3b8 + 0e8799aa8de1b1ac4aeae4d1e8cd5936c66b8a0a60c8d17ec23bc9c9f94421ea + 2221a87687c127b352e077bf7e975c83c667c8ce603e9d9fa5e02b195c7ad873; + 2527a8ab98898038e9e132a1bc8349a33e2c2f6ec8b6c15b06410177297d30fb + 71cdbe058a0357295dd5134fb95a630055f39d9bd8fa0c759b1b70fbef6147d0 + -14780 + 0fe111e19ddf5d11d40b584a31888aea4a388540b86d14b962a1a7135114d5b3 + a6386eeff2764007a570ba7e82a214667c2949a793fa0a52b7665c09285c6661 + a7557aaa000c8840e08163b7eec08bf24ba275a35d85bfd7813aaafcc22715c8 + 867491333fc048b229ddb4740b1fe6dca1feebf1d40660814feeb1fce97fca12; + fe07ffcd6b14f75efd74ff983e8e5c20f379b0c78907f0e453753b457fb53159 + 617f2f7c1153eaea3ccae7dfb07fbe804b2cb5c8205d657c3106b5a60100cafc + -128312 + a0af624a3f54278e2f7cc77f5d0f10750179bb18e07da059feb0bc98fb1d688f + ce432859e6d7e691ba9b44d467e8a5f1a97979facd2d3c7deec46e121193133b + d60050fad0f5c4b2e042141e4f7a418d3e8a8a308acc96651728f363ee793628 + b9e4e8873e006d41104b1cd6a32282fa2b94ade000dd57dee6dd506cceb34c34; + b46e84be433f9e4f43ffa2f1f5200013953d820d292e3af978c54401a252111c + cb7cec0c48efa4e0e54f05ac18574d5f67bf1a1193d6d58f42183abc9f5bdd49 + -3288 + 18603c0562c4d893a699fa6258e2de57e1cafd3e8ad756806d09db0c1ac9f86c + 0854368c6b459ba6fcaa95db92605c9d791bf4dcc57bfda2ef318b48f7521d73 + 0c83d0b2d331ccc62bbcd2d2771b18bc0b909a3bbebcdcb99020952672f050d4 + d54073b190601f1cc36c25109a5b44f04b6c2240c268bc71cb983d3dfb25cd64; + 1f490365c3b82416fd3866a42a03f3bd643844e0d0cbebb66af54aa35de3efa4 + dc805758c7e714a2acbc1e467520ecc3b9d23ce5133fe79e377e6d4520cd5196 + -129141 + 1a53a3728ef3587085d0b530be8ae95e367927f9ac3b7126445d32b46f3abfad + 0af6ce5107ca06fdfe732d11161b197830a963ec667065c5293bbdf843e23672 + 5bbfb655a175dcf7a4f4fc886a03e3d99cddbd34524c3cd2fb13b793bb2132be + 7e5576020e7bac35c7f5b3dc6b60b2c1cca55f6c98141682efb4737963d0307e; + 09e03f5d67a92cf076740e4261833f244b60da6d02cb42278e2ed9183b6eba3a + 5a8b665434e91158559cfecf1d43d562c307518a3ce889165b406b1aae38a674 + -13299 + 16e0d405bc40f2075f5dd901769def19f15759faef9a5fb6356b94c5c48be3e0 + 101da726c0dfd102c8eb2178640054248225a6fcfac1bab28a09b59cf38703d6 + f268d102402d31a33e52a941a06c64a30a721ed55184edfcfcd85db6fe3db676 + 24fc9112ce6d3d70e8c245f9e2c883290384462065a5a1ab5821b9cf9d61953b; + f86ba345a367c676d99450d4893ef4033201c8ddbb631c4658e4ed24ece98749 + a3c51a17107327742f44cc2e265482522cf39b6bcbff919bcf99fc033d76348f + -88039 + 5dcbe0be5f9f80b46d2e651504b10a3e481257a4bfacf40d8ccacea9d73dc2ea + 2ac88350fbffcbc743d36fe45c80a70a8cbe88be36759943c8e124163e432dec + f54685e704df0e499e495b3b77e64aa0ec0c9fc0bb18ed456a412ac6079fe9d2 + 3e7c00ea9ee49825f5fbb9b704dc720b062d8e08ee8dd410f3d34d255602b56b; + 1cd810875f48ef4c71ecb0997acfa31249ebcc9280e7c2f741038a0529f653d3 + ba2e564cae8cdc2bea40d64b5598493ac2b11902d91b5aec7a0e4a8649aa7ab2 + -77488 + ca4929fb0fed99fd9dca7a21275dbd97185900d4eb362eb09bcb5ce7098aaaa2 + 89b395b440fab27d3b32e0740caac3ce38affbd80b2f8f7b6709294d6770873e + d2cd5ac1ea7b0e1f45bd3cc64b651518dcdf9f25078de7cc4970bf54976a0c8d + 06ea5e15853bb29db5cfccd62db2e64067c90a97eab43039cacf4ee639ffe979; + 9960c3131e513fba3241b1f9bdf44fdbe06c43d9cda9e05fdb451eda3e8df394 + 4c5f72392dddd1b8e84db5bfc2db04985fe3d1226a470c851a489f98d954e6ec + -33449 + 04eacb794c87c19087b793a0e2d187cc7e0b07d265fd6cc46946780d2ab83810 + 414db559bd22f12717b431310c1063534de5e627e983f4967b2fb38f2c5f87b0 + 599bee129778596ba77c27a896020fa55123b8635ff3d9efd56fe406de878146 + 26c7d591e0028dd611d6e84969da716eb740207e49f5ea02f17a3224612def12; + 2cfc4c9cb8ca259bd4ffcbbd498510f8a65cca584303a28ed475584bd1eed111 + 3223cab1344be86c0af9995a769276cddbf70eaf55caee24dd4ec0d06c2f4309 + -53728 + 3a48f8bd8a229fa9f5b9d7c4135ef6c3a251b0ce451f482471f017e963839405 + bec0e71a657e85e024d91c41b3a2f75f478528d8ae1afb70bdc5da255d20de06 + 53dc48614937bd5ea1bce51b166bc9947e248a98cd0fac55f632cded216e9cfd + 2999f607e716f0c769e9fdeff636568674f3601166d6b7f566a369b5d4e2797a; + c37946f39bde051e4a6f91a18fda21ec2bfcc375b703bbb81fdcc7ff15cfedf4 + f37e2f1bb6b71763263e742cfd1ce8a9e004a5e06faf0afc52cf49ea89ad9244 + -69223 + 0592777a8c77c1713a72b2ce18fa318fa93e9c6e6b870e54723a43408b502b74 + cf575224b26f636b4b4f28c90995b87d900f6ca5ec9de6c24e65b747d303d221 + 8ab7d691dd76d7cd63d3456c677987a9e53904a4895f053e3de08786d564c837 + 1d667ed7c9b93f2972de23f4d633bcb5f2e56e446cba103c3cf2cf0931687b13; + 2afb15ae6438b8eef9d07f2491e43e32a2b9287381c7bf036dce2a1b0a3542e5 + 18c9bb7c91d5e5ed96e6effc6609fe40e72c5634cac2cf6fb55b083c3339f40a + -33255 + c0c870e0ee2116aa1182b2b4ee4bfde446f721fa6f1bba3a2ed4fb703855c673 + 91fd148ae6d53a419983999020d845503cc923e012e6869201615fbb8e460360 + 5d02294b16e04ca4a7fa4fa5dc9f12c7156e2a68c02495ea171a8cf13d35af99 + f1f90c61c85c137e2e9d5ca7a6b388afd954f6c93548e5fc3b82e02fe875923e; + a7a7193315006fd85c97226f44708d61ecddad2e184d69b0db62f9fb2bc678bb + 87a64ef3c67f464441927c08fd1b57d17bd644a422b96f151aa2d7920989e4ac + 67151 + 6b31a270025004a878f214f100a764165033e170ea5e4c52832673cf88f9f0b3 + 0c0af67c56e0b2336ba06f447ea4f9cbf052147e35f16216b306ef1211ce4142 + 54fec5270095dfc21f44805ba17565386bc50d45bd1f797ebfa04da5040ad08c + 92fdbcabd91ab43baa827594672dede4f0c71d8302b4cc8cd79270aa0c96d915; + 5d72825031c526d9e1efe34fcbbbd81b91460112f9f950e78e12c777a2994a59 + 4334ab5246765e241435c6faac0a00add3a5adcdafe017121d03b9561d0d3516 + -66199 + ef3211bb4e1bac43430e48594577d1d271f76b46ea7e44542c23524b3d6236bc + 1bfb4fb8d8429f0dd7a6cbc23dc2357ce9c1c72be3acf0a49abc331ac751e68a + 7cf2ea00463d76152e230cb5465ab55d30f09afbd007d2b55393cfe8e348cf9a + c0a66f24e5800bec2b9ff65207afd4b29326525b910381d233bd91d51530d046; + 01a66900161a4ad98f028a8da8f10142b4855712a7eb8ffbd699b13b7cc1a471 + c6f8a0b42c3f2a1109201c7a97e3307f6cb56a1a27f13e4567456d9eac2a2f4e + -71473 + e2eef0e19682ebb44cfdf84c71cb556311432fd948951a661d4240e16bc68d1c + f69b44aaff31db04fa4e8f14f51051d353e50110b1562b27d4493ac5e44a5b07 + 29039f65da48a24be024cf8f1af95ea9b09a08a54eec810e97018add81eec1df + a827fcc698403d01f8c2597a582ac8ddd612e5c2a520b0a2428e8dd35c58002d; + bc607923665e36ea0f8206b358fb3a327ef9a70ddcadb2edc49660891e932726 + 33e905cd1b262eba320c61428704c1d7c354361043ea7233c69f1f691f7c7175 + -101999 + e1b804b366c155d0f1a09fdbede434c669bf7f188daa96b8542dcbc9bbe66b6b + 413d410937fbf0575bf34c0b70ee6cd5737e17b6058e774f889104b2f6cacc7b + 0784aeb89be679694477cf399b605473c2c48887b687a980d7d333632e060e6c + 16edc745e1596492fadffd303c31ce5827c3971c591c68f9e2b4f62c967aaa6f; + 904e1b99ee3f860409a2db3aa4e539b2b50cebc57359f23fca21b5bb1306dc50 + 133d76a1f768d77a6801d73a570b35f6dbbd98aabfa4258445684143cd30a39f + -50858 + 547270fdf41f753a4189d2a9630ab17b4c6034829770a40fd8346e75e92610cf + 076c4e34375857b9e6289287e03c744fd60d76c377a772f6edef3ad30e4afafb + f2e8a830fe203901f5f34ae6e606aa8f6cffb504b166cba4cdd84d011b783b90 + 659233daf8ebd47c77b9d46ac65237e1bfdeebbcc552c7db5cf8dd2d4bd1aa2d; + a6f4bce9adbe5a68585c205d20b05c6662748f5f1dd2458defbc126792837a63 + 8f1110929626587241fd403c8332b2ebbc398d4dd7042b993fc24a18b53e3847 + 81985 + 012aaf1e00a521e8e0639ed116e15476f0f8f1c6dff47d860ed60e0bb7c156f7 + 13aaf4f9f20d64110beef21bde5c7262c4bff80a85a08db1e098d45ec3d69457 + eb3230fdacce04518359bfa7cff0c9d6d015966545e40e1c1da18b760031b2aa + 607f7b77fc8d3ca15cdf6e2abba71853fa74739f83589a4458854ae2db0a4727; +} diff --git a/utils/curve25519.sage b/utils/curve25519.sage new file mode 100644 index 00000000..32ddb17c --- /dev/null +++ b/utils/curve25519.sage @@ -0,0 +1,43 @@ +#! /usr/local/bin/sage +### -*- mode: python; coding: utf-8 -*- + +###-------------------------------------------------------------------------- +### Define the field. + +p = 2^255 - 19; k = GF(p) + +###-------------------------------------------------------------------------- +### Arithmetic implementation. + +def sqrn(x, n): + for i in xrange(n): x = x*x + return x + +def inv(x): + t2 = sqrn(x, 1) # 1 | 2 + u = sqrn(t2, 2) # 3 | 8 + t = u*x # 4 | 9 + t11 = t*t2 # 5 | 11 + u = sqrn(t11, 1) # 6 | 22 + t = u*t # 7 | 2^5 - 1 = 31 + u = sqrn(t, 5) # 12 | 2^10 - 2^5 + t2p10m1 = u*t # 13 | 2^10 - 1 + u = sqrn(t2p10m1, 10) # 23 | 2^20 - 2^10 + t = u*t2p10m1 # 24 | 2^20 - 1 + u = sqrn(t, 20) # 44 | 2^40 - 2^20 + t = u*t # 45 | 2^40 - 1 + u = sqrn(t, 10) # 55 | 2^50 - 2^10 + t2p50m1 = u*t2p10m1 # 56 | 2^50 - 1 + u = sqrn(t2p50m1, 50) # 106 | 2^100 - 2^50 + t = u*t2p50m1 # 107 | 2^100 - 1 + u = sqrn(t, 100) # 207 | 2^200 - 2^100 + t = u*t # 208 | 2^200 - 1 + u = sqrn(t, 50) # 258 | 2^250 - 2^50 + t = u*t2p50m1 # 259 | 2^250 - 1 + u = sqrn(t, 5) # 264 | 2^255 - 2^5 + t = u*t11 # 265 | 2^255 - 21 + return t + +assert inv(k(9))*9 == 1 + +###----- That's all, folks -------------------------------------------------- diff --git a/utils/qfarith-test b/utils/qfarith-test new file mode 100755 index 00000000..3904d42c --- /dev/null +++ b/utils/qfarith-test @@ -0,0 +1,109 @@ +#! /usr/bin/python + +from sys import argv +import catacomb as C + +TESTS = {} +NTEST = 20 + +def test(arg): + def reg(fn, name): + TESTS[name] = fn + return fn + if isinstance(arg, str): return lambda fn: reg(fn, arg) + else: return reg(arg, arg.__name__.replace('_', '-')) + +FIELDS = {} + +class FieldElt (object): + def __init__(me, k, n): + me.k = k + me.v = (C.MP(n)%k.p).storel(k.len) + def __str__(me): return hex(me.v) + @property + def n(me): return C.MP.loadl(me.v) + def __pos__(me): return FieldElt(me.k, me.n) + def __neg__(me): return FieldElt(me.k, -me.n) + def __nonzero__(me): return me.n != 0 + def __add__(me, you): return FieldElt(me.k, me.n + me.k(you).n) + def __radd__(me, you): return FieldElt(me.k, me.k(you).n + me.n) + def __sub__(me, you): return FieldElt(me.k, me.n - me.k(you).n) + def __rsub__(me, you): return FieldElt(me.k, me.k(you).n - me.n) + def __mul__(me, you): return FieldElt(me.k, me.n*me.k(you).n) + def __rmul__(me, you): return FieldElt(me.k, me.k(you).n*me.n) + def __div__(me, you): return FieldElt(me.k, me.n*me.k(you).inv().n) + def __rdiv__(me, you): return FieldElt(me.k, me.k(you).n*me.inv().n) + def inv(me): return FieldElt(me.k, me.k.p.modinv(me.n)) + def sqrt(me): return FieldElt(me.k, me.k.p.modsqrt(me.n)) + @classmethod + def rand(cls, k): + me = cls(k, 0) + me.v = C.rand.block(k. len) + return me + +class Field (object): + def __init__(me, p, len = None): + me.p = C.MP(p) + me.len = len is None and me.p.noctets or len + @classmethod + def register(cls, name, *args, **kw): + FIELDS[name] = cls(*args, **kw) + def rand(me): return FieldElt.rand(me) + def __call__(me, n): + if isinstance(n, FieldElt): + assert n.k is me + return n + else: + return FieldElt(me, n) + +Field.register('f25519', C.MP(0).setbit(255) - 19) + +def binop(k, op): + x = k.rand(); y = k.rand() + print ' %s\n %s\n %s;' % (x, y, op(x, y)) + +def unop(k, op): + x = k.rand() + print ' %s\n %s;' % (x, op(x)) + +@test +def add(k): binop(k, lambda x, y: x + y) + +@test +def sub(k): binop(k, lambda x, y: x - y) + +@test +def mul(k): binop(k, lambda x, y: x*y) + +@test +def sqr(k): unop(k, lambda x: x*x) + +@test +def inv(k): unop(k, lambda x: x and x.inv() or k(0)) + +@test +def mulconst(k): + x = k.rand() + a = C.rand.range(1 << 20) - (1 << 19) + print ' %s %d\n %s;' % (x, a, a*x) + +def mask(): return C.rand.range(2)*0xffffffff + +@test +def condswap(k): + x = k.rand(); y = k.rand(); m = mask() + xx, yy = m and (+y, +x) or (+x, +y) + print ' %s\n %s\n 0x%08x\n %s\n %s;' % (x, y, m, xx, yy) + +@test +def sub_mulc_add_sub_mul(k): + u = k.rand(); v = k.rand(); w = k.rand(); x = k.rand(); y = k.rand(); + a = C.rand.range(1 << 20) - (1 << 19) + print ' %s\n %s %d\n %s\n %s\n %s\n %s;' % \ + (u, v, a, w, x, y, ((u - v)*a + w)*(x - y)) + +k = FIELDS[argv[1]] +for t in argv[2:]: + print '%s {' % t + for i in xrange(NTEST): TESTS[t](k) + print '}' -- 2.11.0