add mulT

misc/low_test.cpp

@@ -0,0 +1,42 @@
+#include "../src/low_funct.hpp"
+
+#define MCL_VINT_FIXED_BUFFER
+#define MCL_SIZEOF_UNIT 4
+#define MCL_MAX_BIT_SIZE 384
+#include <mcl/vint.hpp>
+#include <cybozu/test.hpp>
+#include <cybozu/xorshift.hpp>
+
+void mul3(uint32_t z[6], const uint32_t x[3], uint32_t y[3])
+{
+	return mcl::mulT<3>(z, x, y);
+}
+
+template<class RG>
+void setRand(uint32_t *x, size_t n, RG& rg)
+{
+	for (size_t i = 0; i < n; i++) {
+		x[i] = rg.get32();
+	}
+}
+
+CYBOZU_TEST_AUTO(mul3)
+{
+	cybozu::XorShift rg;
+	uint32_t x[3];
+	uint32_t y[3];
+	uint32_t z[6];
+	for (size_t i = 0; i < 1000; i++) {
+		setRand(x, 3, rg);
+		setRand(y, 3, rg);
+		mcl::Vint vx, vy;
+		vx.setArray(x, 3);
+		vy.setArray(y, 3);
+		printf("vx=%s\n", vx.getStr(16).c_str());
+		printf("vy=%s\n", vy.getStr(16).c_str());
+		vx *= vy;
+		printf("xy=%s\n", vx.getStr(16).c_str());
+		mul3(z, x, y);
+		CYBOZU_TEST_EQUAL_ARRAY(z, vx.getUnit(), 6);
+	}
+}

src/low_funct.hpp

@@ -7,13 +7,14 @@
 */
 #include <stdint.h>
 #include <stdlib.h>
+#include <assert.h>
 
-// for 32bit not full version
+// only for 32bit not full bit prime version
 
 namespace mcl {
 
 template<size_t N>
-void copyT(uint32_t *y, const uint32_t *x)
+void copyT(uint32_t y[N], const uint32_t x[N])
 {
 	for (size_t i = 0; i < N; i++) {
 		y[i] = x[i];
@@ -21,7 +22,7 @@ void copyT(uint32_t *y, const uint32_t *x)
 }
 
 template<size_t N>
-void addT(uint32_t *z, const uint32_t *x, const uint32_t *y)
+void addT(uint32_t z[N], const uint32_t x[N], const uint32_t y[N])
 {
 	bool c = false;
 	for (size_t i = 0; i < N; i++) {
@@ -29,10 +30,11 @@ void addT(uint32_t *z, const uint32_t *x, const uint32_t *y)
 		z[i] = uint32_t(v);
 		c = (v >> 32) != 0;
 	}
+	assert(!c);
 }
 
 template<size_t N>
-bool subT(uint32_t *z, const uint32_t *x, const uint32_t *y)
+bool subT(uint32_t z[N], const uint32_t x[N], const uint32_t y[N])
 {
 	bool c = false;
 	for (size_t i = 0; i < N; i++) {
@@ -43,8 +45,47 @@ bool subT(uint32_t *z, const uint32_t *x, const uint32_t *y)
 	return c;
 }
 
+// [return:z[N]] = x[N] * y
 template<size_t N>
-void addModT(uint32_t *z, const uint32_t *x, const uint32_t *y, const uint32_t *p)
+uint32_t mulUnitT(uint32_t z[N], const uint32_t x[N], uint32_t y)
+{
+	uint32_t H = 0;
+	for (size_t i = 0; i < N; i++) {
+		uint64_t v = uint64_t(x[i]) * y;
+		v += H;
+		z[i] = uint32_t(v);
+		H = uint32_t(v >> 32);
+	}
+	return H;
+}
+
+// [return:z[N]] = z[N] + x[N] * z
+template<size_t N>
+uint32_t addMulUnitT(uint32_t z[N], const uint32_t x[N], uint32_t y)
+{
+	uint32_t H = 0;
+	for (size_t i = 0; i < N; i++) {
+		uint64_t v = uint64_t(x[i]) * y;
+		v += H;
+		v += z[i];
+		z[i] = uint32_t(v);
+		H = uint32_t(v >> 32);
+	}
+	return H;
+}
+
+// z[N * 2] = x[N] * y[N]
+template<size_t N>
+void mulT(uint32_t z[N * 2], const uint32_t x[N], const uint32_t y[N])
+{
+	z[N] = mulUnitT<N>(z, x, y[0]);
+	for (size_t i = 1; i < N; i++) {
+		z[N + i] = addMulUnitT<N>(&z[i], x, y[i]);
+	}
+}
+
+template<size_t N>
+void addModT(uint32_t z[N], const uint32_t x[N], const uint32_t y[N], const uint32_t p[N])
 {
 	uint32_t t[N];
 	addT<N>(z, x, y);