add fast mod for SECP256k1

include/mcl/op.hpp

@@ -146,7 +146,7 @@ enum Mode {
 enum PrimeMode {
 	PM_GENERIC = 0,
 	PM_NIST_P192,
-	PM_NIST_P256,
+	PM_SECP256K1,
 	PM_NIST_P521
 };
 

include/mcl/vint.hpp

@@ -874,6 +874,55 @@ public:
 	T& operator[](size_t n) { verify(n); return v_[n]; }
 };
 
+#if MCL_SIZEOF_UNIT == 8
+/*
+	M = 1 << 256
+	a = M mod p = (1 << 32) + 0x3d1
+	[H:L] mod p = H * a + L
+
+	if H = L = M - 1, t = H * a + L = aM + (M - a - 1)
+	H' = a, L' = M - a - 1
+	t' = H' * a + L' = M + (a^2 - a - 1)
+	H'' = 1, L'' = a^2 - a - 1
+	t'' = H'' * a + L'' = a^2 - 1
+*/
+inline void mcl_fpDbl_mod_SECP256K1(Unit *z, const Unit *x, const Unit *p)
+{
+	const Unit a = (uint64_t(1) << 32) + 0x3d1;
+	Unit buf[5];
+	buf[4] = mulu1(buf, x + 4, 4, a); // H * a
+	buf[4] += addN(buf, buf, x, 4); // t = H * a + L
+	Unit x2[2];
+	x2[0] = mulUnit(&x2[1], buf[4], a);
+	Unit x3 = addN(buf, buf, x2, 2);
+	if (x3) {
+		x3 = addu1(buf + 2, buf + 2, 2, Unit(1)); // t' = H' * a + L'
+		if (x3) {
+			x3 = addu1(buf, buf, 4, a);
+			assert(x3 == 0);
+		}
+	}
+	if (fp::isGreaterOrEqualArray(buf, p, 4)) {
+		subN(z, buf, p, 4);
+	} else {
+		fp::copyArray(z, buf, 4);
+	}
+}
+
+inline void mcl_fp_mul_SECP256K1(Unit *z, const Unit *x, const Unit *y, const Unit *p)
+{
+	Unit xy[8];
+	mulNM(xy, x, 4, y, 4);
+	mcl_fpDbl_mod_SECP256K1(z, xy, p);
+}
+inline void mcl_fp_sqr_SECP256K1(Unit *y, const Unit *x, const Unit *p)
+{
+	Unit xx[8];
+	sqrN(xx, x, 4);
+	mcl_fpDbl_mod_SECP256K1(y, xx, p);
+}
+#endif
+
 } // vint
 
 /**

src/fp.cpp

@@ -500,6 +500,13 @@ void Op::init(const std::string& mstr, size_t maxBitSize, Mode mode, size_t mclM
 		isMont = false;
 		isFastMod = true;
 	}
+#endif
+#if defined(MCL_USE_VINT) && MCL_SIZEOF_UNIT == 8
+	if (mp == mpz_class("0xfffffffffffffffffffffffffffffffffffffffffffffffffffffffefffffc2f")) {
+		primeMode = PM_SECP256K1;
+		isMont = false;
+		isFastMod = true;
+	}
 #endif
 	N = (bitSize + UnitBitSize - 1) / UnitBitSize;
 	switch (N) {
@@ -548,6 +555,13 @@ void Op::init(const std::string& mstr, size_t maxBitSize, Mode mode, size_t mclM
 	if (primeMode == PM_NIST_P521) {
 		fpDbl_mod = &mcl_fpDbl_mod_NIST_P521L;
 	}
+#endif
+#if defined(MCL_USE_VINT) && MCL_SIZEOF_UNIT == 8
+	if (primeMode == PM_SECP256K1) {
+		fp_mul = &mcl::vint::mcl_fp_mul_SECP256K1;
+		fp_sqr = &mcl::vint::mcl_fp_sqr_SECP256K1;
+		fpDbl_mod = &mcl::vint::mcl_fpDbl_mod_SECP256K1;
+	}
 #endif
 	fp::initForMont(*this, p, mode);
 	sq.set(mp);

test/vint_test.cpp

@@ -5,6 +5,7 @@
 #include <set>
 #include <cybozu/benchmark.hpp>
 #include <cybozu/test.hpp>
+#include <cybozu/xorshift.hpp>
 
 #define PUT(x) std::cout << #x "=" << x << std::endl;
 
@@ -1220,4 +1221,35 @@ CYBOZU_TEST_AUTO(divUnit)
 		}
 	}
 }
+
+void compareMod(const uint64_t *x, const uint64_t *p)
+{
+	uint64_t y1[4] = {};
+	uint64_t y2[4] = {};
+	mcl::vint::divNM((uint64_t*)0, 0, y1, x, 8, p, 4);
+	mcl::vint::mcl_fpDbl_mod_SECP256K1(y2, x, p);
+	CYBOZU_TEST_EQUAL_ARRAY(y1, y2, 4);
+}
+CYBOZU_TEST_AUTO(SECP256k1)
+{
+	const uint64_t F = uint64_t(-1);
+	const uint64_t p[4] = { uint64_t(0xfffffffefffffc2full), F, F, F };
+	const uint64_t tbl[][8] = {
+		{ 0, 0, 0, 0, 0, 0, 0, 0 },
+		{ F, F, F, F, F, F, F, F },
+		{ F, F, F, F, 1, 0, 0, 0 },
+	};
+	for (size_t i = 0; i < CYBOZU_NUM_OF_ARRAY(tbl); i++) {
+		const uint64_t *x = tbl[i];
+		compareMod(x, p);
+	}
+	cybozu::XorShift rg;
+	for (size_t i = 0; i < 100; i++) {
+		uint64_t x[8];
+		for (int j = 0; j < 8; j++) {
+			x[j] = rg();
+		}
+		compareMod(x, p);
+	}
+}
 #endif