From a0dcc763f8f4325e4ecf8be573f8fd5ad7076df6 Mon Sep 17 00:00:00 2001 From: MITSUNARI Shigeo Date: Sun, 12 Apr 2020 18:18:31 +0900 Subject: [PATCH] refactor mapto --- include/mcl/ec.hpp | 224 +++++++++++++++++++++---------------- include/mcl/mapto_wb19.hpp | 108 +++++++++--------- 2 files changed, 185 insertions(+), 147 deletions(-) diff --git a/include/mcl/ec.hpp b/include/mcl/ec.hpp index 0b9376b..20cc7ab 100644 --- a/include/mcl/ec.hpp +++ b/include/mcl/ec.hpp @@ -51,6 +51,11 @@ enum ModeCoeffA { namespace local { +/* + elliptic class E must have + member variables of type Fp x, y, z + static member a_, b_, specialA_ +*/ // x is negative <=> x < half(:=(p+1)/2) <=> a = 1 template bool get_a_flag(const F& x) @@ -68,61 +73,69 @@ bool get_a_flag(const mcl::Fp2T& x) } // mcl::ec::local -template -void normalizeJacobi(F& x, F& y, F& z) +template +void normalizeJacobi(E& P) { - if (z.isZero()) return; - F::inv(z, z); + typedef typename E::Fp F; + if (P.z.isZero()) return; + F::inv(P.z, P.z); F rz2; - F::sqr(rz2, z); - x *= rz2; - y *= rz2; - y *= z; - z = 1; + F::sqr(rz2, P.z); + P.x *= rz2; + P.y *= rz2; + P.y *= P.z; + P.z = 1; } // (x/z^2, y/z^3) -template -bool isEqualJacobi(const F& x1, const F& y1, const F& z1, const F& x2, const F& y2, const F& z2) +template +bool isEqualJacobi(const E& P1, const E& P2) { + typedef typename E::Fp F; F s1, s2, t1, t2; - F::sqr(s1, z1); - F::sqr(s2, z2); - F::mul(t1, x1, s2); - F::mul(t2, x2, s1); + F::sqr(s1, P1.z); + F::sqr(s2, P2.z); + F::mul(t1, P1.x, s2); + F::mul(t2, P2.x, s1); if (t1 != t2) return false; - F::mul(t1, y1, s2); - F::mul(t2, y2, s1); - t1 *= z2; - t2 *= z1; + F::mul(t1, P1.y, s2); + F::mul(t2, P2.y, s1); + t1 *= P2.z; + t2 *= P1.z; return t1 == t2; } // Y^2 == X(X^2 + aZ^4) + bZ^6 -template -bool isValidJacobi(const F& x, const F& y, const F& z, const F& a, const F& b) +template +bool isValidJacobi(const E& P) { + typedef typename E::Fp F; F y2, x2, z2, z4, t; - F::sqr(x2, x); - F::sqr(y2, y); - F::sqr(z2, z); + F::sqr(x2, P.x); + F::sqr(y2, P.y); + F::sqr(z2, P.z); F::sqr(z4, z2); - F::mul(t, z4, a); + F::mul(t, z4, E::a_); t += x2; - t *= x; + t *= P.x; z4 *= z2; - z4 *= b; + z4 *= E::b_; t += z4; return y2 == t; } /* + M > S + A a = 0 2M + 5S + 14A a = -3 2M + 7S + 15A generic 3M + 7S + 15A + M == S + a = 0 3M + 4S + 13A + a = -3 3M + 6S + 14A + generic 4M + 6S + 14A */ template -void dblJacobi(E& R, const E& P, int specialA, const typename E::Fp& a) +void dblJacobi(E& R, const E& P) { typedef typename E::Fp F; if (P.isZero()) { @@ -133,13 +146,19 @@ void dblJacobi(E& R, const E& P, int specialA, const typename E::Fp& a) F x2, y2, xy, t; F::sqr(x2, P.x); F::sqr(y2, P.y); - F::add(xy, P.x, y2); - F::sqr(y2, y2); - F::sqr(xy, xy); - xy -= x2; - xy -= y2; - xy += xy; - switch (specialA) { + if (sizeof(F) <= 32) { // M == S + F::mul(xy, P.x, y2); + xy += xy; + F::sqr(y2, y2); + } else { // M > S + A + F::add(xy, P.x, y2); + F::sqr(y2, y2); + F::sqr(xy, xy); + xy -= x2; + xy -= y2; + } + xy += xy; // 4xy^2 + switch (E::specialA_) { case Zero: F::add(t, x2, x2); x2 += t; @@ -158,11 +177,11 @@ void dblJacobi(E& R, const E& P, int specialA, const typename E::Fp& a) case GenericA: default: if (isPzOne) { - t = a; + t = E::a_; } else { F::sqr(t, P.z); F::sqr(t, t); - t *= a; + t *= E::a_; } t += x2; x2 += x2; @@ -187,9 +206,10 @@ void dblJacobi(E& R, const E& P, int specialA, const typename E::Fp& a) } // 7M + 4S + 7A -template -void addJacobi(E& R, const E& P, const E& Q, int specialA, const F& a) +template +void addJacobi(E& R, const E& P, const E& Q) { + typedef typename E::Fp F; if (P.isZero()) { R = Q; return; } if (Q.isZero()) { R = P; return; } bool isPzOne = P.z.isOne(); @@ -230,7 +250,7 @@ void addJacobi(E& R, const E& P, const E& Q, int specialA, const F& a) r -= S1; if (H.isZero()) { if (r.isZero()) { - ec::dblJacobi(R, P, specialA, a); + ec::dblJacobi(R, P); } else { R.clear(); } @@ -263,43 +283,46 @@ void addJacobi(E& R, const E& P, const E& Q, int specialA, const F& a) F::sub(R.y, U1, H3); } -template -void normalizeProj(F& x, F& y, F& z) +template +void normalizeProj(E& P) { - if (z.isZero()) return; - F::inv(z, z); - x *= z; - y *= z; - z = 1; + typedef typename E::Fp F; + if (P.z.isZero()) return; + F::inv(P.z, P.z); + P.x *= P.z; + P.y *= P.z; + P.z = 1; } // (Y^2 - bZ^2)Z = X(X^2 + aZ^2) -template -bool isValidProj(const F& x, const F& y, const F& z, const F& a, const F& b) +template +bool isValidProj(const E& P) { + typedef typename E::Fp F; F y2, x2, z2, t; - F::sqr(x2, x); - F::sqr(y2, y); - F::sqr(z2, z); - F::mul(t, a, z2); + F::sqr(x2, P.x); + F::sqr(y2, P.y); + F::sqr(z2, P.z); + F::mul(t, E::a_, z2); t += x2; - t *= x; - z2 *= b; + t *= P.x; + z2 *= E::b_; y2 -= z2; - y2 *= z; + y2 *= P.z; return y2 == t; } // (x/z, y/z) -template -bool isEqualProj(const F& x1, const F& y1, const F& z1, const F& x2, const F& y2, const F& z2) +template +bool isEqualProj(const E& P1, const E& P2) { + typedef typename E::Fp F; F t1, t2; - F::mul(t1, x1, z2); - F::mul(t2, x2, z1); + F::mul(t1, P1.x, P2.z); + F::mul(t2, P2.x, P1.z); if (t1 != t2) return false; - F::mul(t1, y1, z2); - F::mul(t2, y2, z1); + F::mul(t1, P1.y, P2.z); + F::mul(t2, P2.y, P1.z); return t1 == t2; } @@ -309,16 +332,17 @@ bool isEqualProj(const F& x1, const F& y1, const F& z1, const F& x2, const F& y2 sqr| 4| 5| 5 add| 11|12|12 */ -template -void dblProj(E& R, const E& P, int specialA, const F& a) +template +void dblProj(E& R, const E& P) { + typedef typename E::Fp F; if (P.isZero()) { R.clear(); return; } const bool isPzOne = P.z.isOne(); F w, t, h; - switch (specialA) { + switch (E::specialA_) { case Zero: F::sqr(w, P.x); F::add(t, w, w); @@ -338,10 +362,10 @@ void dblProj(E& R, const E& P, int specialA, const F& a) case GenericA: default: if (isPzOne) { - w = a; + w = E::a_; } else { F::sqr(w, P.z); - w *= a; + w *= E::a_; } F::sqr(t, P.x); w += t; @@ -379,9 +403,10 @@ void dblProj(E& R, const E& P, int specialA, const F& a) sqr| 2 add| 7 */ -template -void addProj(E& R, const E& P, const E& Q, int specialA, const F& a) +template +void addProj(E& R, const E& P, const E& Q) { + typedef typename E::Fp F; if (P.isZero()) { R = Q; return; } if (Q.isZero()) { R = P; return; } bool isPzOne = P.z.isOne(); @@ -404,7 +429,7 @@ void addProj(E& R, const E& P, const E& Q, int specialA, const F& a) v -= r; if (v.isZero()) { if (A == PyQz) { - dblProj(R, P, specialA, a); + dblProj(R, P); } else { R.clear(); } @@ -442,22 +467,25 @@ void addProj(E& R, const E& P, const E& Q, int specialA, const F& a) } // y^2 == (x^2 + a)x + b -template -bool isValidAffine(const F& x, const F& y, const F& a, const F& b) +template +bool isValidAffine(const E& P) { + typedef typename E::Fp F; + assert(!P.z.isZero()); F y2, t; - F::sqr(y2, y); - F::sqr(t, x); - t += a; - t *= x; - t += b; + F::sqr(y2, P.y); + F::sqr(t, P.x); + t += E::a_; + t *= P.x; + t += E::b_; return y2 == t; } // y^2 = x^3 + ax + b -template -static inline void dblAffine(E& R, const E& P, const F& a) +template +static inline void dblAffine(E& R, const E& P) { + typedef typename E::Fp F; if (P.isZero()) { R.clear(); return; @@ -470,7 +498,7 @@ static inline void dblAffine(E& R, const E& P, const F& a) F::sqr(t, P.x); F::add(s, t, t); t += s; - t += a; + t += E::a_; F::add(s, P.y, P.y); t /= s; F::sqr(s, t); @@ -484,16 +512,17 @@ static inline void dblAffine(E& R, const E& P, const F& a) R.z = 1; } -template -void addAffine(E& R, const E& P, const E& Q, const F& a) +template +void addAffine(E& R, const E& P, const E& Q) { + typedef typename E::Fp F; if (P.isZero()) { R = Q; return; } if (Q.isZero()) { R = P; return; } F t; F::sub(t, Q.x, P.x); if (t.isZero()) { if (P.y == Q.y) { - dblAffine(R, P, a); + dblAffine(R, P); } else { R.clear(); } @@ -552,17 +581,17 @@ public: private: bool isValidAffine() const { - return ec::isValidAffine(x, y, a_, b_); + return ec::isValidAffine(*this); } public: void normalize() { switch (mode_) { case ec::Jacobi: - ec::normalizeJacobi(x, y, z); + ec::normalizeJacobi(*this); break; case ec::Proj: - ec::normalizeProj(x, y, z); + ec::normalizeProj(*this); break; } } @@ -629,10 +658,10 @@ public: { switch (mode_) { case ec::Jacobi: - if (!ec::isValidJacobi(x, y, z, a_, b_)) return false; + if (!ec::isValidJacobi(*this)) return false; break; case ec::Proj: - if (!ec::isValidProj(x, y, z, a_, b_)) return false; + if (!ec::isValidProj(*this)) return false; break; case ec::Affine: if (z.isZero()) return true; @@ -664,26 +693,27 @@ public: { switch (mode_) { case ec::Jacobi: - ec::dblJacobi(R, P, specialA_, a_); + ec::dblJacobi(R, P); break; case ec::Proj: - ec::dblProj(R, P, specialA_, a_); + ec::dblProj(R, P); break; case ec::Affine: - ec::dblAffine(R, P, a_); + ec::dblAffine(R, P); break; } } - static inline void add(EcT& R, const EcT& P, const EcT& Q) { + static inline void add(EcT& R, const EcT& P, const EcT& Q) + { switch (mode_) { case ec::Jacobi: - ec::addJacobi(R, P, Q, specialA_, a_); + ec::addJacobi(R, P, Q); break; case ec::Proj: - ec::addProj(R, P, Q, specialA_, a_); + ec::addProj(R, P, Q); break; case ec::Affine: - ec::addAffine(R, P, Q, a_); + ec::addAffine(R, P, Q); break; } } @@ -1075,9 +1105,9 @@ public: { switch (mode_) { case ec::Jacobi: - return ec::isEqualJacobi(x, y, z, rhs.x, rhs.y, rhs.z); + return ec::isEqualJacobi(*this, rhs); case ec::Proj: - return ec::isEqualProj(x, y, z, rhs.x, rhs.y, rhs.z); + return ec::isEqualProj(*this, rhs); case ec::Affine: default: return x == rhs.x && y == rhs.y && z == rhs.z; diff --git a/include/mcl/mapto_wb19.hpp b/include/mcl/mapto_wb19.hpp index 3c96d41..92ab315 100644 --- a/include/mcl/mapto_wb19.hpp +++ b/include/mcl/mapto_wb19.hpp @@ -43,6 +43,9 @@ template struct PointT { typedef F Fp; F x, y, z; + static F a_; + static F b_; + static int specialA_; bool isZero() const { return z.isZero(); @@ -53,15 +56,24 @@ struct PointT { y.clear(); z.clear(); } +#if 0 + bool isEqual(const PointT& rhs) const + { + return ec::isEqualJacobi(*this, rhs); + } +#endif }; +template F PointT::a_; +template F PointT::b_; +template int PointT::specialA_; + } // mcl::local template struct MapToG2_WB19 { + typedef local::PointT Point; Fp2 xi; - Fp2 Ell2p_a; - Fp2 Ell2p_b; Fp half; mpz_class sqrtConst; // (p^2 - 9) / 16 Fp2 root4[4]; @@ -71,7 +83,6 @@ struct MapToG2_WB19 { Fp2 ynum[4]; Fp2 yden[4]; int draftVersion_; - typedef local::PointT Point; void setDraftVersion(int version) { draftVersion_ = version; @@ -174,16 +185,12 @@ struct MapToG2_WB19 { void dbl(Point& Q, const Point& P) const { dblT(Q, P); - } - void dbl(G2& Q, const G2& P) const - { - dblT(Q, P); -// G2::dbl(Q, P); +// ec::dblJacobi(Q, P); } // P is on y^2 = x^3 + Ell2p_a x + Ell2p_b bool isValidPoint(const Point& P) const { - return ec::isValidJacobi(P.x, P.y, P.z, Ell2p_a, Ell2p_b); + return ec::isValidJacobi(P); } bool isValidPoint(const G2& P) const { @@ -194,10 +201,11 @@ struct MapToG2_WB19 { bool b; xi.a = -2; xi.b = -1; - Ell2p_a.a = 0; - Ell2p_a.b = 240; - Ell2p_b.a = 1012; - Ell2p_b.b = 1012; + Point::a_.a = 0; + Point::a_.b = 240; + Point::b_.a = 1012; + Point::b_.b = 1012; + Point::specialA_ = ec::GenericA; half = -1; half /= 2; sqrtConst = Fp::getOp().mp; @@ -386,11 +394,11 @@ struct MapToG2_WB19 { // (t^2 * xi)^2 + (t^2 * xi) Fp2::add(deno, t2xi2, t2xi); Fp2::add(nume, deno, 1); - nume *= Ell2p_b; + nume *= Point::b_; if (deno.isZero()) { - Fp2::mul(deno, Ell2p_a, xi); + Fp2::mul(deno, Point::a_, xi); } else { - deno *= -Ell2p_a; + deno *= -Point::a_; } Fp2 u, v; { @@ -398,8 +406,8 @@ struct MapToG2_WB19 { Fp2::sqr(deno2, deno); Fp2::mul(v, deno2, deno); - Fp2::mul(u, Ell2p_b, v); - Fp2::mul(tmp, Ell2p_a, nume); + Fp2::mul(u, Point::b_, v); + Fp2::mul(tmp, Point::a_, nume); tmp *= deno2; u += tmp; Fp2::sqr(tmp, nume); @@ -506,20 +514,20 @@ struct MapToG2_WB19 { den += den2; Fp2 x0_num, x0_den; Fp2::add(x0_num, den, 1); - x0_num *= Ell2p_b; + x0_num *= Point::b_; if (den.isZero()) { - Fp2::mul(x0_den, Ell2p_a, xi); + Fp2::mul(x0_den, Point::a_, xi); } else { - Fp2::mul(x0_den, -Ell2p_a, den); + Fp2::mul(x0_den, -Point::a_, den); } Fp2 x0_den2, x0_den3, gx0_den, gx0_num; Fp2::sqr(x0_den2, x0_den); Fp2::mul(x0_den3, x0_den2, x0_den); gx0_den = x0_den3; - Fp2::mul(gx0_num, Ell2p_b, gx0_den); + Fp2::mul(gx0_num, Point::b_, gx0_den); Fp2 tmp, tmp1, tmp2; - Fp2::mul(tmp, Ell2p_a, x0_num); + Fp2::mul(tmp, Point::a_, x0_num); tmp *= x0_den2; gx0_num += tmp; Fp2::sqr(tmp, x0_num); @@ -581,22 +589,22 @@ struct MapToG2_WB19 { { G2 t[16]; t[0] = P; - dbl(t[1], t[0]); - add(t[4], t[1], t[0]); - add(t[2], t[4], t[1]); - add(t[3], t[2], t[1]); - add(t[11], t[3], t[1]); - add(t[9], t[11], t[1]); - add(t[10], t[9], t[1]); - add(t[5], t[10], t[1]); - add(t[7], t[5], t[1]); - add(t[15], t[7], t[1]); - add(t[13], t[15], t[1]); - add(t[6], t[13], t[1]); - add(t[14], t[6], t[1]); - add(t[12], t[14], t[1]); - add(t[8], t[12], t[1]); - dbl(t[1], t[6]); + G2::dbl(t[1], t[0]); + G2::add(t[4], t[1], t[0]); + G2::add(t[2], t[4], t[1]); + G2::add(t[3], t[2], t[1]); + G2::add(t[11], t[3], t[1]); + G2::add(t[9], t[11], t[1]); + G2::add(t[10], t[9], t[1]); + G2::add(t[5], t[10], t[1]); + G2::add(t[7], t[5], t[1]); + G2::add(t[15], t[7], t[1]); + G2::add(t[13], t[15], t[1]); + G2::add(t[6], t[13], t[1]); + G2::add(t[14], t[6], t[1]); + G2::add(t[12], t[14], t[1]); + G2::add(t[8], t[12], t[1]); + G2::dbl(t[1], t[6]); const struct { uint32_t n; @@ -618,21 +626,21 @@ struct MapToG2_WB19 { }; for (size_t j = 0; j < CYBOZU_NUM_OF_ARRAY(tbl); j++) { const uint32_t n = tbl[j].n; - for (size_t i = 0; i < n; i++) dbl(t[1], t[1]); - add(t[1], t[1], t[tbl[j].idx]); + for (size_t i = 0; i < n; i++) G2::dbl(t[1], t[1]); + G2::add(t[1], t[1], t[tbl[j].idx]); } - for (size_t i = 0; i < 5; i++) dbl(t[1], t[1]); - add(out, t[1], t[2]); + for (size_t i = 0; i < 5; i++) G2::dbl(t[1], t[1]); + G2::add(out, t[1], t[2]); } void mx_chain(G2& Q, const G2& P) const { G2 T; - dbl(T, P); + G2::dbl(T, P); const size_t tbl[] = { 2, 3, 9, 32, 16 }; for (size_t i = 0; i < CYBOZU_NUM_OF_ARRAY(tbl); i++) { - add(T, T, P); + G2::add(T, T, P); for (size_t j = 0; j < tbl[i]; j++) { - dbl(T, T); + G2::dbl(T, T); } } Q = T; @@ -644,12 +652,12 @@ struct MapToG2_WB19 { #else G2 work, work2; h2_chain(work, P); - dbl(work2, work); - add(work2, work, work2); + G2::dbl(work2, work); + G2::add(work2, work, work2); mx_chain(work, work2); mx_chain(work, work); - neg(work2, work2); - add(Q, work, work2); + G2::neg(work2, work2); + G2::add(Q, work, work2); #endif } template