start to impl naive optimal ate

Makefile

@@ -3,7 +3,7 @@ OBJ_DIR=obj
 EXE_DIR=bin
 
 SRC_SRC=fp.cpp
-TEST_SRC=fp_test.cpp ec_test.cpp fp_util_test.cpp window_method_test.cpp elgamal_test.cpp fp_tower_test.cpp gmp_test.cpp
+TEST_SRC=fp_test.cpp ec_test.cpp fp_util_test.cpp window_method_test.cpp elgamal_test.cpp fp_tower_test.cpp gmp_test.cpp bn_test.cpp
 ifeq ($(CPU),x86-64)
   TEST_SRC+=fp_generator_test.cpp mont_fp_test.cpp
 endif

include/mcl/bn.hpp

@@ -135,13 +135,14 @@ struct ParamT {
 	void init(const CurveParam& cp = CurveFp254BNb, fp::Mode mode = fp::FP_AUTO)
 	{
 		z = cp.z;
-		Fp::setModulo(z.get_str(), 10, mode);
 		const int pCoff[] = { 1, 6, 24, 36, 36 };
 		const int rCoff[] = { 1, 6, 18, 36, 36 };
 		const int tCoff[] = { 1, 0,  6,  0,  0 };
 		eval(p, z, pCoff);
 		eval(r, z, rCoff);
 		eval(t, z, tCoff);
+		Fp::setModulo(p.get_str(), 10, mode);
+		Fp2::init(cp.xi_a);
 		b = cp.b; // set b before calling Fp::setModulo
 		half = Fp(1) / Fp(2);
 		Fp2 xi(cp.xi_a, 1);
@@ -163,6 +164,7 @@ struct ParamT {
 		Fp2::power(tmp, xi, (p * p - 1) / 6);
 		assert(tmp.b.isZero());
 		Fp::sqr(Z, tmp.a);
+
 		const mpz_class largest_c = abs(6 * z + 2);
 		useNAF = getGoodRepl(siTbl, largest_c);
 		getGoodRepl(zReplTbl, abs(z)); // QQQ : snark
@@ -173,5 +175,78 @@ struct ParamT {
 	}
 };
 
+namespace naive {
+/*
+	v is the line arising in the addition of Q1 and Q2 in G2 evaluated at point P in G1
+*/
+template<class Fp12, class G2, class G1>
+void evalLine(Fp12& v, const G2& Q1, const G2& Q2, const G1& P)
+{
+	typedef typename G2::Fp Fp2;
+	Q1.normalize();
+	Q2.normalize();
+	P.normalize();
+	Fp2 lambda;
+	if (Q1.x == Q2.x) {
+		lambda = Q1.x * Q1.x / (Q1.y + Q1.y);
+	} else {
+		lambda = (Q1.y - Q2.y) / (Q1.x - Q1.x);
+	}
+	v = lambda * (P.x - Q1.x) + Q1.y;
+}
+
+template<class G>
+void Frobenius(G& y, const G& x)
+{
+	typedef typename G::Fp Fp;
+	mpz_class& mp = Fp::getOp().mp;
+	Fp::power(y.x, x.x, mp);
+	Fp::power(y.y, x.y, mp);
+	Fp::power(y.z, x.z, mp);
+}
+
+template<class Fp12, class G2, class G1>
+void optimalAtePairing(Fp12& f, const G2& Q, const G1& P, const mpz_class& u)
+{
+	mpz_class r = abs(6 * u + 2);
+	mpz_class p = (((36 * u + 36) * u + 24) * u + 6) * u + 1;
+	mpz_class n = (((36 * u + 36) * u + 18) * u + 6) * u + 1;
+	G2 T = Q;
+	Fp12 t;
+	f = 1;
+	const int c = (int)mcl::gmp::getBitSize(r);
+	for (int i = c - 1; i >= 0; i--) {
+		Fp12::sqr(f, f);
+		evalLine(t, T, T, P);
+		f *= t;
+		G2::dbl(T, T);
+		if (mcl::gmp::testBit(r, i)) {
+			evalLine(t, T, Q, P);
+			f *= t;
+			T += Q;
+		}
+	}
+	G2 Q1, Q2;
+	Frobenius(Q1, Q);
+	Frobenius(Q2, Q1);
+	if (u < 0) {
+		G2::neg(T, T);
+		Fp12::inv(f, f);
+	}
+	evalLine(t, T, Q1, P);
+	f *= t;
+	T += Q1;
+	evalLine(t, T, -Q2, P);
+	f *= t;
+	mpz_class a = p * p * p;
+	a *= a;
+	a *= a;
+	a = (a - 1) / n;
+	Fp12::power(f, f, a);
+}
+
+} // mcl::bn::naive
+
+
 } } // mcl::bn
 

include/mcl/ec.hpp

@@ -469,21 +469,6 @@ public:
 		R.z = P.z;
 #endif
 	}
-	static inline void mulArray(EcT& z, const EcT& x, const fp::Unit *y, size_t yn, bool isNegative)
-	{
-		x.normalize();
-		EcT tmp;
-		const EcT *px = &x;
-		if (&z == &x) {
-			tmp = x;
-			px = &tmp;
-		}
-		z.clear();
-		fp::powerGeneric(z, *px, y, yn, EcT::add, EcT::dbl);
-		if (isNegative) {
-			neg(z, z);
-		}
-	}
 	template<class tag, size_t maxBitSize, template<class _tag, size_t _maxBitSize>class FpT>
 	static inline void mul(EcT& z, const EcT& x, const FpT<tag, maxBitSize>& y)
 	{
@@ -642,6 +627,22 @@ public:
 		return R.isZero();
 	}
 	bool operator!=(const EcT& rhs) const { return !operator==(rhs); }
+private:
+	static inline void mulArray(EcT& z, const EcT& x, const fp::Unit *y, size_t yn, bool isNegative)
+	{
+		x.normalize();
+		EcT tmp;
+		const EcT *px = &x;
+		if (&z == &x) {
+			tmp = x;
+			px = &tmp;
+		}
+		z.clear();
+		fp::powerGeneric(z, *px, y, yn, EcT::add, EcT::dbl);
+		if (isNegative) {
+			neg(z, z);
+		}
+	}
 };
 
 template<class Fp> Fp EcT<Fp>::a_;

test/bn_test.cpp

@@ -0,0 +1,19 @@
+#include <cybozu/test.hpp>
+#include <mcl/gmp_util.hpp>
+#include <mcl/bn.hpp>
+#include <mcl/ec.hpp>
+
+typedef mcl::FpT<mcl::FpTag, 256> Fp;
+typedef mcl::Fp2T<Fp> Fp2;
+typedef mcl::Fp6T<Fp> Fp6;
+typedef mcl::Fp12T<Fp> Fp12;
+typedef mcl::EcT<Fp> G1;
+typedef mcl::EcT<Fp2> G2;
+
+typedef mcl::bn::ParamT<Fp> Param;
+
+CYBOZU_TEST_AUTO(naive)
+{
+	Param param;
+	param.init();
+}