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 <limits.h>],
+[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])
$(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 --------------------------------------------------
--- /dev/null
+/* -*-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 <stdio.h>
+
+#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<i} x_j 2^{o_i}, so
+ * x = X_10. We see, inductively, that |X_i| < 2^{31+o_{i-1}}: X_0 = 0;
+ * |x_i| <= 2^30; and |X_{i+1}| = |X_i + x_i 2^{o_i}| <= |X_i| + 2^{30+o_i}
+ * < 2^{31+o_i}. Then x = X_9 + 2^230 x_9, and we have better bounds for
+ * x_9, so
+ *
+ * -2^235 < x + 19 c_9 < 2^255 + 2^235
+ *
+ * Here, the x_i are signed, so we must be cautious about bithacking them.
+ */
+ c = ASR(piece, x9, 25); x9 = (upiece)x9&M25;
+ x0 += 19*c; c = ASR(piece, x0, 26); x0 = (upiece)x0&M26;
+ x1 += c; c = ASR(piece, x1, 25); x1 = (upiece)x1&M25;
+ x2 += c; c = ASR(piece, x2, 26); x2 = (upiece)x2&M26;
+ x3 += c; c = ASR(piece, x3, 25); x3 = (upiece)x3&M25;
+ x4 += c; c = ASR(piece, x4, 26); x4 = (upiece)x4&M26;
+ x5 += c; c = ASR(piece, x5, 25); x5 = (upiece)x5&M25;
+ x6 += c; c = ASR(piece, x6, 26); x6 = (upiece)x6&M26;
+ x7 += c; c = ASR(piece, x7, 25); x7 = (upiece)x7&M25;
+ x8 += c; c = ASR(piece, x8, 26); x8 = (upiece)x8&M26;
+ x9 += c; c = ASR(piece, x9, 25); x9 = (upiece)x9&M25;
+
+ /* Now we have a slightly fiddly job to do. If c = +1, or if c = 0 and
+ * x >= 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<i} y_j 2^{o_i}, so
+ * y = Y_26. We see, inductively, that |Y_i| < 2^{31+o_{i-1}}: Y_0 = 0;
+ * |y_i| <= 2^14; and |Y_{i+1}| = |Y_i + y_i 2^{o_i}| <= |Y_i| + 2^{14+o_i}
+ * < 2^{15+o_i}. Then x = Y_25 + 2^246 y_25, and we have better bounds for
+ * y_25, so
+ *
+ * -2^251 < y + 19 c_25 < 2^255 + 2^251
+ *
+ * Here, the y_i are signed, so we must be cautious about bithacking them.
+ *
+ * (Rather closer than the 10-piece case above, but still doable in one
+ * pass.)
+ */
+ c = 19*ASR(piece, y[NPIECE - 1], 9);
+ y[NPIECE - 1] = (upiece)y[NPIECE - 1]&M9;
+ for (i = 0; i < NPIECE; i++) {
+ wd = PIECEWD(i);
+ y[i] += c;
+ c = ASR(piece, y[i], wd);
+ y[i] = (upiece)y[i]&MASK(wd);
+ }
+
+ /* Now the addition or subtraction. */
+ m = SIGN(c);
+ d = m&1;
+
+ d += y[0] + (19 ^ (M10&m));
+ yy[0] = d&M10;
+ d >>= 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 <mLib/report.h>
+#include <mLib/testrig.h>
+
+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 -------------------------------------------------*/
--- /dev/null
+/* -*-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 <mLib/bits.h>
+
+#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
--- /dev/null
+/* -*-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 <limits.h>
+
+#include <mLib/bits.h>
+
+/*----- 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 <stdio.h>, <catacomb/mp.h>, and
+ * <catacomb/mptext.h> 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
--- /dev/null
+### 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;
+}
--- /dev/null
+#! /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 --------------------------------------------------
--- /dev/null
+#! /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 '}'
~mdw / catacomb / commitdiff