From 56811bbb125cd1f4f7351241c4615ec6042c90ef Mon Sep 17 00:00:00 2001
From: MITSUNARI Shigeo <herumi@nifty.com>
Date: Fri, 17 Mar 2017 15:17:57 +0900
Subject: [PATCH] add constant time scalar mul for EC

---
 include/mcl/ec.hpp   | 15 ++++++++++---
 include/mcl/util.hpp | 50 +++++++++++++++++++++++++++++++++++++-------
 test/ec_test.cpp     | 38 +++++++--------------------------
 3 files changed, 63 insertions(+), 40 deletions(-)

diff --git a/include/mcl/ec.hpp b/include/mcl/ec.hpp
index f4268dc..16fb43e 100644
--- a/include/mcl/ec.hpp
+++ b/include/mcl/ec.hpp
@@ -589,13 +589,18 @@ public:
 		mulArray(z, x, gmp::getUnit(y), abs(y.get_mpz_t()->_mp_size), y < 0);
 	}
 	template<class tag, size_t maxBitSize, template<class _tag, size_t _maxBitSize>class FpT>
-	static inline void mul_s(EcT& z, const EcT& x, const FpT<tag, maxBitSize>& y)
+	static inline void mulCT(EcT& z, const EcT& x, const FpT<tag, maxBitSize>& y)
 	{
 		fp::Block b;
 		y.getBlock(b);
 		mulArray(z, x, b.p, b.n, false, true);
 	}
-	static inline void mul_s(EcT& z, const EcT& x, const mpz_class& y)
+	static inline void mulCT(EcT& z, const EcT& x, int y)
+	{
+		const fp::Unit u = abs(y);
+		mulArray(z, x, &u, 1, y < 0, true);
+	}
+	static inline void mulCT(EcT& z, const EcT& x, const mpz_class& y)
 	{
 		mulArray(z, x, gmp::getUnit(y), abs(y.get_mpz_t()->_mp_size), y < 0, true);
 	}
@@ -776,7 +781,11 @@ private:
 			px = &tmp;
 		}
 		z.clear();
-		fp::powGeneric(z, *px, y, yn, EcT::add, EcT::dbl, constTime);
+		if (constTime) {
+			fp::powGenericCT(z, *px, y, yn, EcT::add, EcT::dbl);
+		} else {
+			fp::powGeneric(z, *px, y, yn, EcT::add, EcT::dbl);
+		}
 		if (isNegative) {
 			neg(z, z);
 		}
diff --git a/include/mcl/util.hpp b/include/mcl/util.hpp
index 4cdf229..565fea1 100644
--- a/include/mcl/util.hpp
+++ b/include/mcl/util.hpp
@@ -195,7 +195,7 @@ void getRandVal(T *out, RG& rg, const T *in, size_t bitSize)
 	@note &out != x and out = the unit element of G
 */
 template<class G, class T>
-void powGeneric(G& out, const G& x, const T *y, size_t n, void mul(G&, const G&, const G&) , void sqr(G&, const G&), bool constTime = false)
+void powGeneric(G& out, const G& x, const T *y, size_t n, void mul(G&, const G&, const G&) , void sqr(G&, const G&))
 {
 #if 0
 	assert(&out != &x);
@@ -262,10 +262,6 @@ void powGeneric(G& out, const G& x, const T *y, size_t n, void mul(G&, const G&,
 	} else {
 		out = x;
 	}
-	G dummy;
-	if (constTime) {
-		dummy = x;
-	}
 	for (int i = (int)n - 1; i >= 0; i--) {
 		T v = y[i];
 		for (int j = m - 2; j >= 0; j -= 2) {
@@ -274,8 +270,6 @@ void powGeneric(G& out, const G& x, const T *y, size_t n, void mul(G&, const G&,
 			T idx = (v >> j) & 3;
 			if (idx > 0) {
 				mul(out, out, tbl[idx - 1]);
-			} else if (constTime) {
-				mul(dummy, dummy, tbl[0]);
 			}
 		}
 		m = (int)sizeof(T) * 8;
@@ -283,6 +277,48 @@ void powGeneric(G& out, const G& x, const T *y, size_t n, void mul(G&, const G&,
 #endif
 }
 
+/*
+	constant time pow
+	@note depends on bit length of y[n]
+*/
+template<class G, class T>
+void powGenericCT(G& out, const G& x, const T *y, size_t n, void mul(G&, const G&, const G&) , void sqr(G&, const G&))
+{
+	assert(&out != &x);
+	while (n > 0) {
+		if (y[n - 1]) break;
+		n--;
+	}
+	if (n == 0) return;
+	G tbl[4]; // tbl = { x, x^2, x^3 }
+	tbl[0] = x; // discard
+	tbl[1] = x;
+	mul(tbl[2], x, x); tbl[2].normalize();
+	mul(tbl[3], tbl[2], x); tbl[3].normalize();
+	T v = y[n - 1];
+	int m = cybozu::bsr<T>(v);
+	if (m & 1) {
+		m--;
+		T idx = (v >> m) & 3;
+		assert(idx > 0);
+		out = tbl[idx];
+	} else {
+		out = x;
+	}
+	G *pTbl[] = { &tbl[0], &out, &out, &out };
+
+	for (int i = (int)n - 1; i >= 0; i--) {
+		T v = y[i];
+		for (int j = m - 2; j >= 0; j -= 2) {
+			sqr(out, out);
+			sqr(out, out);
+			T idx = (v >> j) & 3;
+			mul(*pTbl[idx], *pTbl[idx], tbl[idx]);
+		}
+		m = (int)sizeof(T) * 8;
+	}
+}
+
 /*
 	shortcut of multiplication by Unit
 */
diff --git a/test/ec_test.cpp b/test/ec_test.cpp
index 1cfa377..210503b 100644
--- a/test/ec_test.cpp
+++ b/test/ec_test.cpp
@@ -2,7 +2,6 @@
 #define CYBOZU_TEST_DISABLE_AUTO_RUN
 #include <cybozu/test.hpp>
 #include <cybozu/benchmark.hpp>
-#include <cybozu/xorshift.hpp>
 #include <mcl/gmp_util.hpp>
 
 #include <mcl/fp.hpp>
@@ -308,38 +307,17 @@ struct Test {
 		}
 		Fp::setIoMode(mcl::IoAuto);
 	}
-	void bench_mul_s(void f(Ec&, const Ec&, const Zn&)) const
+	void mulCT() const
 	{
 		Fp x(para.gx);
 		Fp y(para.gy);
-		Ec P(x, y);
-		cybozu::XorShift rg;
-		std::vector<double> tv(100);
-		for (size_t i = 0; i < tv.size(); i++) {
-			Zn r;
-			r.setRand(rg);
-			CYBOZU_BENCH_C("", 30, f, P, P, r);
-			tv[i] = cybozu::bench::g_clk.getClock();
-printf("%.1f\n", tv[i]);
-		}
-		double ave = 0;
-		for (size_t i = 0; i < tv.size(); i++) {
-			ave += tv[i];
-		}
-		ave /= tv.size();
-		double v = 0;
-		for (size_t i = 0; i < tv.size(); i++) {
-			double t = tv[i] - ave;
-			v += t * t;
+		Ec P(x, y), Q1, Q2;
+		for (int i = 0; i < 100; i++) {
+			Zn r = i;
+			Ec::mul(Q1, P, r);
+			Ec::mulCT(Q2, P, r);
+			CYBOZU_TEST_EQUAL(Q1, Q2);
 		}
-		v /= tv.size();
-		v = sqrt(v);
-		printf("ave %.2f v %.2f\n", ave, v);
-	}
-	void mul_s() const
-	{
-		bench_mul_s(Ec::mul);
-		bench_mul_s(Ec::mul_s);
 	}
 
 	template<class F>
@@ -379,7 +357,7 @@ mul 499.00usec
 		squareRoot();
 		str();
 		ioMode();
-//		mul_s();
+		mulCT();
 	}
 private:
 	Test(const Test&);