add fp_mont for gmp

5 files changed, 46 insertions, 15 deletions

diff --git a/include/mcl/fp.hpp b/include/mcl/fp.hpp
index e052584..8cfd6cf 100644
--- a/include/mcl/fp.hpp
+++ b/include/mcl/fp.hpp
@@ -119,10 +119,10 @@ public:
 #else
 		if (mode == fp::FP_LLVM || mode == fp::FP_LLVM_MONT) mode = fp::FP_AUTO;
 #endif
-		if (mode == fp::FP_AUTO) mode = fp::FP_GMP;
+		if (mode == fp::FP_AUTO) mode = fp::FP_GMP_MONT;
 
-		op_.isMont = mode == fp::FP_LLVM_MONT || mode == fp::FP_XBYAK;
-		if (mode == fp::FP_LLVM_MONT) {
+		op_.isMont = mode == fp::FP_GMP_MONT || mode == fp::FP_LLVM_MONT || mode == fp::FP_XBYAK;
+		if (mode == fp::FP_GMP_MONT || mode == fp::FP_LLVM_MONT) {
 			op_.fp_mul = fp_montW;
 			op_.fp_sqr = fp_montSqrW;
 			op_.fpDbl_mod = fp_montRedW;
diff --git a/sample/bench.cpp b/sample/bench.cpp
index ee21fe9..4bc2a9b 100644
--- a/sample/bench.cpp
+++ b/sample/bench.cpp
@@ -62,13 +62,14 @@ void benchFp(size_t bitSize, int mode)
 	};
 	for (size_t i = 0; i < CYBOZU_NUM_OF_ARRAY(tbl); i++) {
 		if (bitSize != 0 && tbl[i].bitSize != bitSize) continue;
-		if (mode & (1 << 0)) benchFpSub(tbl[i].p, tbl[i].x, tbl[i].y, mcl::fp::FP_GMP);
+		if (mode & 1) benchFpSub(tbl[i].p, tbl[i].x, tbl[i].y, mcl::fp::FP_GMP);
+		if (mode & 2) benchFpSub(tbl[i].p, tbl[i].x, tbl[i].y, mcl::fp::FP_GMP_MONT);
 #ifdef MCL_USE_LLVM
-		if (mode & (1 << 1)) benchFpSub(tbl[i].p, tbl[i].x, tbl[i].y, mcl::fp::FP_LLVM);
-		if (mode & (1 << 2)) benchFpSub(tbl[i].p, tbl[i].x, tbl[i].y, mcl::fp::FP_LLVM_MONT);
+		if (mode & 4) benchFpSub(tbl[i].p, tbl[i].x, tbl[i].y, mcl::fp::FP_LLVM);
+		if (mode & 8) benchFpSub(tbl[i].p, tbl[i].x, tbl[i].y, mcl::fp::FP_LLVM_MONT);
 #endif
 #ifdef MCL_USE_XBYAK
-		if (mode & (1 << 3)) benchFpSub(tbl[i].p, tbl[i].x, tbl[i].y, mcl::fp::FP_XBYAK);
+		if (mode & 16) benchFpSub(tbl[i].p, tbl[i].x, tbl[i].y, mcl::fp::FP_XBYAK);
 #endif
 	}
 }
@@ -113,13 +114,14 @@ void benchEc(size_t bitSize, int mode, mcl::ec::Mode ecMode)
 	for (size_t i = 0; i < CYBOZU_NUM_OF_ARRAY(tbl); i++) {
 		if (bitSize != 0 && tbl[i].bitSize != bitSize) continue;
 		benchEcSub(tbl[i], mcl::fp::FP_AUTO, ecMode);
-		if (mode & (1 << 0)) benchEcSub(tbl[i], mcl::fp::FP_GMP, ecMode);
+		if (mode & 1) benchEcSub(tbl[i], mcl::fp::FP_GMP, ecMode);
+		if (mode & 2) benchEcSub(tbl[i], mcl::fp::FP_GMP_MONT, ecMode);
 #ifdef MCL_USE_LLVM
-		if (mode & (1 << 1)) benchEcSub(tbl[i], mcl::fp::FP_LLVM, ecMode);
-		if (mode & (1 << 2)) benchEcSub(tbl[i], mcl::fp::FP_LLVM_MONT, ecMode);
+		if (mode & 4) benchEcSub(tbl[i], mcl::fp::FP_LLVM, ecMode);
+		if (mode & 8) benchEcSub(tbl[i], mcl::fp::FP_LLVM_MONT, ecMode);
 #endif
 #ifdef MCL_USE_XBYAK
-		if (mode & (1 << 3)) benchEcSub(tbl[i], mcl::fp::FP_XBYAK, ecMode);
+		if (mode & 16) benchEcSub(tbl[i], mcl::fp::FP_XBYAK, ecMode);
 #endif
 	}
 }
@@ -183,7 +185,7 @@ int main(int argc, char *argv[])
 	std::string ecModeStr;
 	cybozu::Option opt;
 	opt.appendOpt(&bitSize, 0, "b", ": bitSize");
-	opt.appendOpt(&mode, 0, "m", ": mode(0:all, sum of 1:gmp, 2:llvm, 8:llvm+mont, 8:xbyak");
+	opt.appendOpt(&mode, 0, "m", ": mode(0:all, sum of 1:gmp, 2:gmp+mont, 4:llvm, 8:llvm+mont, 16:xbyak");
 	opt.appendBoolOpt(&ecOnly, "ec", ": ec only");
 	opt.appendBoolOpt(&fpOnly, "fp", ": fp only");
 	opt.appendBoolOpt(&misc, "misc", ": other benchmark");
@@ -202,12 +204,12 @@ int main(int argc, char *argv[])
 		opt.usage();
 		return 1;
 	}
-	if (mode < 0 || mode > 15) {
+	if (mode < 0 || mode > 31) {
 		printf("bad mode %d\n", mode);
 		opt.usage();
 		return 1;
 	}
-	if (mode == 0) mode = 15;
+	if (mode == 0) mode = 31;
 	if (misc) {
 		benchToStr16();
 		benchFromStr16();
diff --git a/sample/rawbench.cpp b/sample/rawbench.cpp
index 721a4ca..10e9e4a 100644
--- a/sample/rawbench.cpp
+++ b/sample/rawbench.cpp
@@ -116,6 +116,7 @@ int main(int argc, char *argv[])
 		}
 		printf("prime=%s\n", p);
 		benchRaw(tbl[i], mcl::fp::FP_GMP);
+		benchRaw(tbl[i], mcl::fp::FP_GMP_MONT);
 #ifdef MCL_USE_LLVM
 		benchRaw(tbl[i], mcl::fp::FP_LLVM);
 		benchRaw(tbl[i], mcl::fp::FP_LLVM_MONT);
diff --git a/src/fp.cpp b/src/fp.cpp
index ed6b8ab..df9aab2 100644
--- a/src/fp.cpp
+++ b/src/fp.cpp
@@ -277,6 +277,32 @@ struct OpeFunc {
 			}
 		}
 	}
+	// z[N] <- montRed(xy[N * 2])
+	static inline void fp_montRedPUC(Unit *z, const Unit *xy, const Unit *p, Unit rp)
+	{
+		Unit t[N * 2];
+		Unit buf[N * 2 + 1];
+		clearArray(t, N + 1, N * 2);
+		Unit *c = buf;
+		Unit q = xy[0] * rp;
+		low_mul_Unit<N>(t, p, q);
+		buf[N * 2] = low_add<N * 2>(buf, xy, t);
+		c++;
+		for (size_t i = 1; i < N; i++) {
+			q = c[0] * rp;
+			low_mul_Unit<N>(t, p, q);
+			// QQQ
+			mpn_add_n((mp_limb_t*)c, (const mp_limb_t*)c, (const mp_limb_t*)t, N * 2 + 1 - i);
+			c++;
+		}
+		if (c[N]) {
+			low_sub<N>(z, c, p);
+		} else {
+			if (low_sub<N>(z, c, p)) {
+				memcpy(z, c, N * sizeof(Unit));
+			}
+		}
+	}
 	static inline void fp_invOpC(Unit *y, const Unit *x, const Op& op)
 	{
 		mpz_class my;
@@ -376,6 +402,8 @@ struct OpeFunc {
 		fpDbl_mulPre = OpeFunc<n>::fpDbl_mulPreC; \
 		fpDbl_sqrPre = OpeFunc<n>::fpDbl_sqrPreC; \
 		fpDbl_modP = OpeFunc<n>::fpDbl_modPC; \
+		montPU = OpeFunc<n>::fp_montPUC; \
+		montRedPU = OpeFunc<n>::fp_montRedPUC; \
 		SET_OP_LLVM(n)
 
 #ifdef MCL_USE_XBYAK
diff --git a/test/bn_test.cpp b/test/bn_test.cpp
index 3430bd7..a48f9a0 100644
--- a/test/bn_test.cpp
+++ b/test/bn_test.cpp
@@ -100,7 +100,7 @@ int main(int argc, char *argv[])
 {
 	cybozu::Option opt;
 	std::string mode;
-	opt.appendOpt(&mode, "auto", "m", ": mode(gmp/llvm/llvm_mont/xbyak)");
+	opt.appendOpt(&mode, "auto", "m", ": mode(gmp/gmp_mont/llvm/llvm_mont/xbyak)");
 	if (!opt.parse(argc, argv)) {
 		opt.usage();
 		return 1;