github.com/platonnetwork/platon-go@v0.7.6/cases/tool/win/bls_win/include/mcl/ec.hpp (about)

     1  #pragma once
     2  /**
     3  	@file
     4  	@brief elliptic curve
     5  	@author MITSUNARI Shigeo(@herumi)
     6  	@license modified new BSD license
     7  	http://opensource.org/licenses/BSD-3-Clause
     8  */
     9  #include <stdlib.h>
    10  #include <cybozu/exception.hpp>
    11  #include <mcl/op.hpp>
    12  #include <mcl/util.hpp>
    13  
    14  //#define MCL_EC_USE_AFFINE
    15  
    16  #ifdef _MSC_VER
    17  	#pragma warning(push)
    18  	#pragma warning(disable : 4458)
    19  #endif
    20  
    21  namespace mcl {
    22  
    23  template<class _Fp> class Fp2T;
    24  
    25  namespace ec {
    26  
    27  enum Mode {
    28  	Jacobi = 0,
    29  	Proj = 1
    30  };
    31  
    32  namespace local {
    33  
    34  // x is negative <=> x < half(:=(p+1)/2) <=> a = 1
    35  template<class Fp>
    36  bool get_a_flag(const Fp& x)
    37  {
    38  	return x.isNegative();
    39  }
    40  
    41  // Im(x) is negative <=> Im(x)  < half(:=(p+1)/2) <=> a = 1
    42  
    43  template<class Fp>
    44  bool get_a_flag(const mcl::Fp2T<Fp>& x)
    45  {
    46  	return get_a_flag(x.b); // x = a + bi
    47  }
    48  
    49  } // mcl::ec::local
    50  
    51  } // mcl::ec
    52  
    53  /*
    54  	elliptic curve
    55  	y^2 = x^3 + ax + b (affine)
    56  	y^2 = x^3 + az^4 + bz^6 (Jacobi) x = X/Z^2, y = Y/Z^3
    57  */
    58  template<class _Fp>
    59  class EcT : public fp::Serializable<EcT<_Fp> > {
    60  	enum {
    61  		zero,
    62  		minus3,
    63  		generic
    64  	};
    65  public:
    66  	typedef _Fp Fp;
    67  	typedef _Fp BaseFp;
    68  #ifdef MCL_EC_USE_AFFINE
    69  	Fp x, y;
    70  	bool inf_;
    71  #else
    72  	Fp x, y, z;
    73  	static int mode_;
    74  #endif
    75  	static Fp a_;
    76  	static Fp b_;
    77  	static int specialA_;
    78  	static int ioMode_;
    79  	/*
    80  		order_ is the order of G2 which is the subgroup of EcT<Fp2>.
    81  		check the order of the elements if verifyOrder_ is true
    82  	*/
    83  	static bool verifyOrder_;
    84  	static mpz_class order_;
    85  	static void (*mulArrayGLV)(EcT& z, const EcT& x, const fp::Unit *y, size_t yn, bool isNegative, bool constTime);
    86  	/* default constructor is undefined value */
    87  	EcT() {}
    88  	EcT(const Fp& _x, const Fp& _y)
    89  	{
    90  		set(_x, _y);
    91  	}
    92  	bool isNormalized() const
    93  	{
    94  #ifdef MCL_EC_USE_AFFINE
    95  		return true;
    96  #else
    97  		return isZero() || z.isOne();
    98  #endif
    99  	}
   100  #ifndef MCL_EC_USE_AFFINE
   101  private:
   102  	void normalizeJacobi()
   103  	{
   104  		assert(!z.isZero());
   105  		Fp rz2;
   106  		Fp::inv(z, z);
   107  		Fp::sqr(rz2, z);
   108  		x *= rz2;
   109  		y *= rz2;
   110  		y *= z;
   111  		z = 1;
   112  	}
   113  	void normalizeProj()
   114  	{
   115  		assert(!z.isZero());
   116  		Fp::inv(z, z);
   117  		x *= z;
   118  		y *= z;
   119  		z = 1;
   120  	}
   121  	// Y^2 == X(X^2 + aZ^4) + bZ^6
   122  	bool isValidJacobi() const
   123  	{
   124  		Fp y2, x2, z2, z4, t;
   125  		Fp::sqr(x2, x);
   126  		Fp::sqr(y2, y);
   127  		Fp::sqr(z2, z);
   128  		Fp::sqr(z4, z2);
   129  		Fp::mul(t, z4, a_);
   130  		t += x2;
   131  		t *= x;
   132  		z4 *= z2;
   133  		z4 *= b_;
   134  		t += z4;
   135  		return y2 == t;
   136  	}
   137  	// (Y^2 - bZ^2)Z = X(X^2 + aZ^2)
   138  	bool isValidProj() const
   139  	{
   140  		Fp y2, x2, z2, t;
   141  		Fp::sqr(x2, x);
   142  		Fp::sqr(y2, y);
   143  		Fp::sqr(z2, z);
   144  		Fp::mul(t, a_, z2);
   145  		t += x2;
   146  		t *= x;
   147  		z2 *= b_;
   148  		y2 -= z2;
   149  		y2 *= z;
   150  		return y2 == t;
   151  	}
   152  #endif
   153  	// y^2 == (x^2 + a)x + b
   154  	static inline bool isValid(const Fp& _x, const Fp& _y)
   155  	{
   156  		Fp y2, t;
   157  		Fp::sqr(y2, _y);
   158  		Fp::sqr(t, _x);
   159  		t += a_;
   160  		t *= _x;
   161  		t += b_;
   162  		return y2 == t;
   163  	}
   164  public:
   165  	void normalize()
   166  	{
   167  #ifndef MCL_EC_USE_AFFINE
   168  		if (isNormalized()) return;
   169  		switch (mode_) {
   170  		case ec::Jacobi:
   171  			normalizeJacobi();
   172  			break;
   173  		case ec::Proj:
   174  			normalizeProj();
   175  			break;
   176  		}
   177  #endif
   178  	}
   179  	static void normalize(EcT& y, const EcT& x)
   180  	{
   181  		y = x;
   182  		y.normalize();
   183  	}
   184  	static inline void init(const Fp& a, const Fp& b, int mode = ec::Jacobi)
   185  	{
   186  		a_ = a;
   187  		b_ = b;
   188  		if (a_.isZero()) {
   189  			specialA_ = zero;
   190  		} else if (a_ == -3) {
   191  			specialA_ = minus3;
   192  		} else {
   193  			specialA_ = generic;
   194  		}
   195  		ioMode_ = 0;
   196  		verifyOrder_ = false;
   197  		order_ = 0;
   198  		mulArrayGLV = 0;
   199  #ifdef MCL_EC_USE_AFFINE
   200  		cybozu::disable_warning_unused_variable(mode);
   201  #else
   202  		assert(mode == ec::Jacobi || mode == ec::Proj);
   203  		mode_ = mode;
   204  #endif
   205  	}
   206  	/*
   207  		verify the order of *this is equal to order if order != 0
   208  		in constructor, set, setStr, operator<<().
   209  	*/
   210  	static void setOrder(const mpz_class& order)
   211  	{
   212  		if (order != 0) {
   213  			verifyOrder_ = true;
   214  			order_ = order;
   215  		} else {
   216  			verifyOrder_ = false;
   217  			// don't clear order_ because it is used for isValidOrder()
   218  		}
   219  	}
   220  	static void setMulArrayGLV(void f(EcT& z, const EcT& x, const fp::Unit *y, size_t yn, bool isNegative, bool constTime))
   221  	{
   222  		mulArrayGLV = f;
   223  	}
   224  	static inline void init(bool *pb, const char *astr, const char *bstr, int mode = ec::Jacobi)
   225  	{
   226  		Fp a, b;
   227  		a.setStr(pb, astr);
   228  		if (!*pb) return;
   229  		b.setStr(pb, bstr);
   230  		if (!*pb) return;
   231  		init(a, b, mode);
   232  	}
   233  	// verify the order
   234  	bool isValidOrder() const
   235  	{
   236  		EcT Q;
   237  		EcT::mulGeneric(Q, *this, order_);
   238  		return Q.isZero();
   239  	}
   240  	bool isValid() const
   241  	{
   242  		if (isZero()) return true;
   243  		bool isOK = false;
   244  #ifndef MCL_EC_USE_AFFINE
   245  		if (!z.isOne()) {
   246  			switch (mode_) {
   247  			case ec::Jacobi:
   248  				isOK = isValidJacobi();
   249  				break;
   250  			case ec::Proj:
   251  				isOK = isValidProj();
   252  				break;
   253  			}
   254  		} else
   255  #endif
   256  		{
   257  			isOK = isValid(x, y);
   258  		}
   259  		if (!isOK) return false;
   260  		if (verifyOrder_) return isValidOrder();
   261  		return true;
   262  	}
   263  	void set(bool *pb, const Fp& _x, const Fp& _y, bool verify = true)
   264  	{
   265  		if (verify && !isValid(_x, _y)) {
   266  			*pb = false;
   267  			return;
   268  		}
   269  		x = _x; y = _y;
   270  #ifdef MCL_EC_USE_AFFINE
   271  		inf_ = false;
   272  #else
   273  		z = 1;
   274  #endif
   275  		if (verify && verifyOrder_ && !isValidOrder()) {
   276  			*pb = false;
   277  		} else {
   278  			*pb = true;
   279  		}
   280  	}
   281  	void clear()
   282  	{
   283  #ifdef MCL_EC_USE_AFFINE
   284  		inf_ = true;
   285  #else
   286  		z.clear();
   287  #endif
   288  		x.clear();
   289  		y.clear();
   290  	}
   291  #ifndef MCL_EC_USE_AFFINE
   292  	static inline void dblNoVerifyInfJacobi(EcT& R, const EcT& P)
   293  	{
   294  		Fp S, M, t, y2;
   295  		Fp::sqr(y2, P.y);
   296  		Fp::mul(S, P.x, y2);
   297  		const bool isPzOne = P.z.isOne();
   298  		S += S;
   299  		S += S;
   300  		Fp::sqr(M, P.x);
   301  		switch (specialA_) {
   302  		case zero:
   303  			Fp::add(t, M, M);
   304  			M += t;
   305  			break;
   306  		case minus3:
   307  			if (isPzOne) {
   308  				M -= P.z;
   309  			} else {
   310  				Fp::sqr(t, P.z);
   311  				Fp::sqr(t, t);
   312  				M -= t;
   313  			}
   314  			Fp::add(t, M, M);
   315  			M += t;
   316  			break;
   317  		case generic:
   318  		default:
   319  			if (isPzOne) {
   320  				t = a_;
   321  			} else {
   322  				Fp::sqr(t, P.z);
   323  				Fp::sqr(t, t);
   324  				t *= a_;
   325  			}
   326  			t += M;
   327  			M += M;
   328  			M += t;
   329  			break;
   330  		}
   331  		Fp::sqr(R.x, M);
   332  		R.x -= S;
   333  		R.x -= S;
   334  		if (isPzOne) {
   335  			R.z = P.y;
   336  		} else {
   337  			Fp::mul(R.z, P.y, P.z);
   338  		}
   339  		R.z += R.z;
   340  		Fp::sqr(y2, y2);
   341  		y2 += y2;
   342  		y2 += y2;
   343  		y2 += y2;
   344  		Fp::sub(R.y, S, R.x);
   345  		R.y *= M;
   346  		R.y -= y2;
   347  	}
   348  	static inline void dblNoVerifyInfProj(EcT& R, const EcT& P)
   349  	{
   350  		const bool isPzOne = P.z.isOne();
   351  		Fp w, t, h;
   352  		switch (specialA_) {
   353  		case zero:
   354  			Fp::sqr(w, P.x);
   355  			Fp::add(t, w, w);
   356  			w += t;
   357  			break;
   358  		case minus3:
   359  			Fp::sqr(w, P.x);
   360  			if (isPzOne) {
   361  				w -= P.z;
   362  			} else {
   363  				Fp::sqr(t, P.z);
   364  				w -= t;
   365  			}
   366  			Fp::add(t, w, w);
   367  			w += t;
   368  			break;
   369  		case generic:
   370  		default:
   371  			if (isPzOne) {
   372  				w = a_;
   373  			} else {
   374  				Fp::sqr(w, P.z);
   375  				w *= a_;
   376  			}
   377  			Fp::sqr(t, P.x);
   378  			w += t;
   379  			w += t;
   380  			w += t; // w = a z^2 + 3x^2
   381  			break;
   382  		}
   383  		if (isPzOne) {
   384  			R.z = P.y;
   385  		} else {
   386  			Fp::mul(R.z, P.y, P.z); // s = yz
   387  		}
   388  		Fp::mul(t, R.z, P.x);
   389  		t *= P.y; // xys
   390  		t += t;
   391  		t += t; // 4(xys) ; 4B
   392  		Fp::sqr(h, w);
   393  		h -= t;
   394  		h -= t; // w^2 - 8B
   395  		Fp::mul(R.x, h, R.z);
   396  		t -= h; // h is free
   397  		t *= w;
   398  		Fp::sqr(w, P.y);
   399  		R.x += R.x;
   400  		R.z += R.z;
   401  		Fp::sqr(h, R.z);
   402  		w *= h;
   403  		R.z *= h;
   404  		Fp::sub(R.y, t, w);
   405  		R.y -= w;
   406  	}
   407  #endif
   408  	static inline void dblNoVerifyInf(EcT& R, const EcT& P)
   409  	{
   410  #ifdef MCL_EC_USE_AFFINE
   411  		Fp t, s;
   412  		Fp::sqr(t, P.x);
   413  		Fp::add(s, t, t);
   414  		t += s;
   415  		t += a_;
   416  		Fp::add(s, P.y, P.y);
   417  		t /= s;
   418  		Fp::sqr(s, t);
   419  		s -= P.x;
   420  		Fp x3;
   421  		Fp::sub(x3, s, P.x);
   422  		Fp::sub(s, P.x, x3);
   423  		s *= t;
   424  		Fp::sub(R.y, s, P.y);
   425  		R.x = x3;
   426  		R.inf_ = false;
   427  #else
   428  		switch (mode_) {
   429  		case ec::Jacobi:
   430  			dblNoVerifyInfJacobi(R, P);
   431  			break;
   432  		case ec::Proj:
   433  			dblNoVerifyInfProj(R, P);
   434  			break;
   435  		}
   436  #endif
   437  	}
   438  	static inline void dbl(EcT& R, const EcT& P)
   439  	{
   440  		if (P.isZero()) {
   441  			R.clear();
   442  			return;
   443  		}
   444  		dblNoVerifyInf(R, P);
   445  	}
   446  #ifndef MCL_EC_USE_AFFINE
   447  	static inline void addJacobi(EcT& R, const EcT& P, const EcT& Q, bool isPzOne, bool isQzOne)
   448  	{
   449  		Fp r, U1, S1, H, H3;
   450  		if (isPzOne) {
   451  			// r = 1;
   452  		} else {
   453  			Fp::sqr(r, P.z);
   454  		}
   455  		if (isQzOne) {
   456  			U1 = P.x;
   457  			if (isPzOne) {
   458  				H = Q.x;
   459  			} else {
   460  				Fp::mul(H, Q.x, r);
   461  			}
   462  			H -= U1;
   463  			S1 = P.y;
   464  		} else {
   465  			Fp::sqr(S1, Q.z);
   466  			Fp::mul(U1, P.x, S1);
   467  			if (isPzOne) {
   468  				H = Q.x;
   469  			} else {
   470  				Fp::mul(H, Q.x, r);
   471  			}
   472  			H -= U1;
   473  			S1 *= Q.z;
   474  			S1 *= P.y;
   475  		}
   476  		if (isPzOne) {
   477  			r = Q.y;
   478  		} else {
   479  			r *= P.z;
   480  			r *= Q.y;
   481  		}
   482  		r -= S1;
   483  		if (H.isZero()) {
   484  			if (r.isZero()) {
   485  				dblNoVerifyInf(R, P);
   486  			} else {
   487  				R.clear();
   488  			}
   489  			return;
   490  		}
   491  		if (isPzOne) {
   492  			R.z = H;
   493  		} else {
   494  			Fp::mul(R.z, P.z, H);
   495  		}
   496  		if (!isQzOne) {
   497  			R.z *= Q.z;
   498  		}
   499  		Fp::sqr(H3, H); // H^2
   500  		Fp::sqr(R.y, r); // r^2
   501  		U1 *= H3; // U1 H^2
   502  		H3 *= H; // H^3
   503  		R.y -= U1;
   504  		R.y -= U1;
   505  		Fp::sub(R.x, R.y, H3);
   506  		U1 -= R.x;
   507  		U1 *= r;
   508  		H3 *= S1;
   509  		Fp::sub(R.y, U1, H3);
   510  	}
   511  	static inline void addProj(EcT& R, const EcT& P, const EcT& Q, bool isPzOne, bool isQzOne)
   512  	{
   513  		Fp r, PyQz, v, A, vv;
   514  		if (isQzOne) {
   515  			r = P.x;
   516  			PyQz = P.y;
   517  		} else {
   518  			Fp::mul(r, P.x, Q.z);
   519  			Fp::mul(PyQz, P.y, Q.z);
   520  		}
   521  		if (isPzOne) {
   522  			A = Q.y;
   523  			v = Q.x;
   524  		} else {
   525  			Fp::mul(A, Q.y, P.z);
   526  			Fp::mul(v, Q.x, P.z);
   527  		}
   528  		v -= r;
   529  		if (v.isZero()) {
   530  			if (A == PyQz) {
   531  				dblNoVerifyInf(R, P);
   532  			} else {
   533  				R.clear();
   534  			}
   535  			return;
   536  		}
   537  		Fp::sub(R.y, A, PyQz);
   538  		Fp::sqr(A, R.y);
   539  		Fp::sqr(vv, v);
   540  		r *= vv;
   541  		vv *= v;
   542  		if (isQzOne) {
   543  			R.z = P.z;
   544  		} else {
   545  			if (isPzOne) {
   546  				R.z = Q.z;
   547  			} else {
   548  				Fp::mul(R.z, P.z, Q.z);
   549  			}
   550  		}
   551  		// R.z = 1 if isPzOne && isQzOne
   552  		if (isPzOne && isQzOne) {
   553  			R.z = vv;
   554  		} else {
   555  			A *= R.z;
   556  			R.z *= vv;
   557  		}
   558  		A -= vv;
   559  		vv *= PyQz;
   560  		A -= r;
   561  		A -= r;
   562  		Fp::mul(R.x, v, A);
   563  		r -= A;
   564  		R.y *= r;
   565  		R.y -= vv;
   566  	}
   567  #endif
   568  	static inline void add(EcT& R, const EcT& P, const EcT& Q) {
   569  		if (P.isZero()) { R = Q; return; }
   570  		if (Q.isZero()) { R = P; return; }
   571  		if (&P == &Q) {
   572  			dblNoVerifyInf(R, P);
   573  			return;
   574  		}
   575  #ifdef MCL_EC_USE_AFFINE
   576  		Fp t;
   577  		Fp::neg(t, Q.y);
   578  		if (P.y == t) { R.clear(); return; }
   579  		Fp::sub(t, Q.x, P.x);
   580  		if (t.isZero()) {
   581  			dblNoVerifyInf(R, P);
   582  			return;
   583  		}
   584  		Fp s;
   585  		Fp::sub(s, Q.y, P.y);
   586  		Fp::div(t, s, t);
   587  		R.inf_ = false;
   588  		Fp x3;
   589  		Fp::sqr(x3, t);
   590  		x3 -= P.x;
   591  		x3 -= Q.x;
   592  		Fp::sub(s, P.x, x3);
   593  		s *= t;
   594  		Fp::sub(R.y, s, P.y);
   595  		R.x = x3;
   596  #else
   597  		bool isPzOne = P.z.isOne();
   598  		bool isQzOne = Q.z.isOne();
   599  		switch (mode_) {
   600  		case ec::Jacobi:
   601  			addJacobi(R, P, Q, isPzOne, isQzOne);
   602  			break;
   603  		case ec::Proj:
   604  			addProj(R, P, Q, isPzOne, isQzOne);
   605  			break;
   606  		}
   607  #endif
   608  	}
   609  	static inline void sub(EcT& R, const EcT& P, const EcT& Q)
   610  	{
   611  		EcT nQ;
   612  		neg(nQ, Q);
   613  		add(R, P, nQ);
   614  	}
   615  	static inline void neg(EcT& R, const EcT& P)
   616  	{
   617  		if (P.isZero()) {
   618  			R.clear();
   619  			return;
   620  		}
   621  		R.x = P.x;
   622  		Fp::neg(R.y, P.y);
   623  #ifdef MCL_EC_USE_AFFINE
   624  		R.inf_ = false;
   625  #else
   626  		R.z = P.z;
   627  #endif
   628  	}
   629  	template<class tag, size_t maxBitSize, template<class _tag, size_t _maxBitSize>class FpT>
   630  	static inline void mul(EcT& z, const EcT& x, const FpT<tag, maxBitSize>& y)
   631  	{
   632  		fp::Block b;
   633  		y.getBlock(b);
   634  		mulArray(z, x, b.p, b.n, false);
   635  	}
   636  	static inline void mul(EcT& z, const EcT& x, int64_t y)
   637  	{
   638  		const uint64_t u = fp::abs_(y);
   639  #if MCL_SIZEOF_UNIT == 8
   640  		mulArray(z, x, &u, 1, y < 0);
   641  #else
   642  		uint32_t ua[2] = { uint32_t(u), uint32_t(u >> 32) };
   643  		size_t un = ua[1] ? 2 : 1;
   644  		mulArray(z, x, ua, un, y < 0);
   645  #endif
   646  	}
   647  	static inline void mul(EcT& z, const EcT& x, const mpz_class& y)
   648  	{
   649  		mulArray(z, x, gmp::getUnit(y), gmp::getUnitSize(y), y < 0);
   650  	}
   651  	template<class tag, size_t maxBitSize, template<class _tag, size_t _maxBitSize>class FpT>
   652  	static inline void mulCT(EcT& z, const EcT& x, const FpT<tag, maxBitSize>& y)
   653  	{
   654  		fp::Block b;
   655  		y.getBlock(b);
   656  		mulArray(z, x, b.p, b.n, false, true);
   657  	}
   658  	static inline void mulCT(EcT& z, const EcT& x, const mpz_class& y)
   659  	{
   660  		mulArray(z, x, gmp::getUnit(y), gmp::getUnitSize(y), y < 0, true);
   661  	}
   662  	/*
   663  		0 <= P for any P
   664  		(Px, Py) <= (P'x, P'y) iff Px < P'x or Px == P'x and Py <= P'y
   665  		@note compare function calls normalize()
   666  	*/
   667  	template<class F>
   668  	static inline int compareFunc(const EcT& P_, const EcT& Q_, F comp)
   669  	{
   670  		const bool QisZero = Q_.isZero();
   671  		if (P_.isZero()) {
   672  			if (QisZero) return 0;
   673  			return -1;
   674  		}
   675  		if (QisZero) return 1;
   676  		EcT P(P_), Q(Q_);
   677  		P.normalize();
   678  		Q.normalize();
   679  		int c = comp(P.x, Q.x);
   680  		if (c > 0) return 1;
   681  		if (c < 0) return -1;
   682  		return comp(P.y, Q.y);
   683  	}
   684  	static inline int compare(const EcT& P, const EcT& Q)
   685  	{
   686  		return compareFunc(P, Q, Fp::compare);
   687  	}
   688  	static inline int compareRaw(const EcT& P, const EcT& Q)
   689  	{
   690  		return compareFunc(P, Q, Fp::compareRaw);
   691  	}
   692  	bool isZero() const
   693  	{
   694  #ifdef MCL_EC_USE_AFFINE
   695  		return inf_;
   696  #else
   697  		return z.isZero();
   698  #endif
   699  	}
   700  	static inline bool isMSBserialize()
   701  	{
   702  		return !b_.isZero() && (Fp::BaseFp::getBitSize() & 7) != 0;
   703  	}
   704  	template<class OutputStream>
   705  	void save(bool *pb, OutputStream& os, int ioMode) const
   706  	{
   707  		const char sep = *fp::getIoSeparator(ioMode);
   708  		if (ioMode & IoEcProj) {
   709  			cybozu::writeChar(pb, os, '4'); if (!*pb) return;
   710  			if (sep) {
   711  				cybozu::writeChar(pb, os, sep);
   712  				if (!*pb) return;
   713  			}
   714  			x.save(pb, os, ioMode); if (!*pb) return;
   715  			if (sep) {
   716  				cybozu::writeChar(pb, os, sep);
   717  				if (!*pb) return;
   718  			}
   719  			y.save(pb, os, ioMode); if (!*pb) return;
   720  			if (sep) {
   721  				cybozu::writeChar(pb, os, sep);
   722  				if (!*pb) return;
   723  			}
   724  #ifndef MCL_EC_USE_AFFINE
   725  			z.save(pb, os, ioMode);
   726  #endif
   727  			return;
   728  		}
   729  		EcT P(*this);
   730  		P.normalize();
   731  		if (ioMode & (IoSerialize | IoSerializeHexStr)) {
   732  			const size_t n = Fp::getByteSize();
   733  			const size_t adj = isMSBserialize() ? 0 : 1;
   734  			uint8_t buf[sizeof(Fp) + 1];
   735  			if (Fp::BaseFp::isETHserialization()) {
   736  				const uint8_t c_flag = 0x80;
   737  				const uint8_t b_flag = 0x40;
   738  				const uint8_t a_flag = 0x20;
   739  				if (P.isZero()) {
   740  					buf[0] = c_flag | b_flag;
   741  					memset(buf + 1, 0, n - 1);
   742  				} else {
   743  					cybozu::MemoryOutputStream mos(buf, n);
   744  					P.x.save(pb, mos, IoSerialize); if (!*pb) return;
   745  					uint8_t cba = c_flag;
   746  					if (ec::local::get_a_flag(P.y)) cba |= a_flag;
   747  					buf[0] |= cba;
   748  				}
   749  			} else {
   750  				/*
   751  					if (isMSBserialize()) {
   752  					  // n bytes
   753  					  x | (y.isOdd ? 0x80 : 0)
   754  					} else {
   755  					  // n + 1 bytes
   756  					  (y.isOdd ? 3 : 2), x
   757  					}
   758  				*/
   759  				if (isZero()) {
   760  					memset(buf, 0, n + adj);
   761  				} else {
   762  					cybozu::MemoryOutputStream mos(buf + adj, n);
   763  					P.x.save(pb, mos, IoSerialize); if (!*pb) return;
   764  					if (adj) {
   765  						buf[0] = P.y.isOdd() ? 3 : 2;
   766  					} else {
   767  						if (P.y.isOdd()) {
   768  							buf[n - 1] |= 0x80;
   769  						}
   770  					}
   771  				}
   772  			}
   773  			if (ioMode & IoSerializeHexStr) {
   774  				mcl::fp::writeHexStr(pb, os, buf, n + adj);
   775  			} else {
   776  				cybozu::write(pb, os, buf, n + adj);
   777  			}
   778  			return;
   779  		}
   780  		if (isZero()) {
   781  			cybozu::writeChar(pb, os, '0');
   782  			return;
   783  		}
   784  		if (ioMode & IoEcCompY) {
   785  			cybozu::writeChar(pb, os, P.y.isOdd() ? '3' : '2');
   786  			if (!*pb) return;
   787  			if (sep) {
   788  				cybozu::writeChar(pb, os, sep);
   789  				if (!*pb) return;
   790  			}
   791  			P.x.save(pb, os, ioMode);
   792  		} else {
   793  			cybozu::writeChar(pb, os, '1'); if (!*pb) return;
   794  			if (sep) {
   795  				cybozu::writeChar(pb, os, sep);
   796  				if (!*pb) return;
   797  			}
   798  			P.x.save(pb, os, ioMode); if (!*pb) return;
   799  			if (sep) {
   800  				cybozu::writeChar(pb, os, sep);
   801  				if (!*pb) return;
   802  			}
   803  			P.y.save(pb, os, ioMode);
   804  		}
   805  	}
   806  	template<class InputStream>
   807  	void load(bool *pb, InputStream& is, int ioMode)
   808  	{
   809  #ifdef MCL_EC_USE_AFFINE
   810  		inf_ = false;
   811  #else
   812  		z = 1;
   813  #endif
   814  		if (ioMode & (IoSerialize | IoSerializeHexStr)) {
   815  			const size_t n = Fp::getByteSize();
   816  			const size_t adj = isMSBserialize() ? 0 : 1;
   817  			const size_t n1 = n + adj;
   818  			uint8_t buf[sizeof(Fp) + 1];
   819  			size_t readSize;
   820  			if (ioMode & IoSerializeHexStr) {
   821  				readSize = mcl::fp::readHexStr(buf, n1, is);
   822  			} else {
   823  				readSize = cybozu::readSome(buf, n1, is);
   824  			}
   825  			if (readSize != n1) {
   826  				*pb = false;
   827  				return;
   828  			}
   829  			if (Fp::BaseFp::isETHserialization()) {
   830  				const uint8_t c_flag = 0x80;
   831  				const uint8_t b_flag = 0x40;
   832  				const uint8_t a_flag = 0x20;
   833  				*pb = false;
   834  				if ((buf[0] & c_flag) == 0) { // assume compressed
   835  					return;
   836  				}
   837  				if (buf[0] & b_flag) { // infinity
   838  					if (buf[0] != (c_flag | b_flag)) return;
   839  					for (size_t i = 1; i < n - 1; i++) {
   840  						if (buf[i]) return;
   841  					}
   842  					clear();
   843  					*pb = true;
   844  					return;
   845  				}
   846  				bool a = (buf[0] & a_flag) != 0;
   847  				buf[0] &= ~(c_flag | b_flag | a_flag);
   848  				mcl::fp::local::byteSwap(buf, n);
   849  				x.setArray(pb, buf, n);
   850  				if (!*pb) return;
   851  				getWeierstrass(y, x);
   852  				if (!Fp::squareRoot(y, y)) {
   853  					*pb = false;
   854  					return;
   855  				}
   856  				if (ec::local::get_a_flag(y) ^ a) {
   857  					Fp::neg(y, y);
   858  				}
   859  				return;
   860  			}
   861  			if (fp::isZeroArray(buf, n1)) {
   862  				clear();
   863  				*pb = true;
   864  				return;
   865  			}
   866  			bool isYodd;
   867  			if (adj) {
   868  				char c = buf[0];
   869  				if (c != 2 && c != 3) {
   870  					*pb = false;
   871  					return;
   872  				}
   873  				isYodd = c == 3;
   874  			} else {
   875  				isYodd = (buf[n - 1] >> 7) != 0;
   876  				buf[n - 1] &= 0x7f;
   877  			}
   878  			x.setArray(pb, buf + adj, n);
   879  			if (!*pb) return;
   880  			*pb = getYfromX(y, x, isYodd);
   881  			if (!*pb) return;
   882  		} else {
   883  			char c = 0;
   884  			if (!fp::local::skipSpace(&c, is)) {
   885  				*pb = false;
   886  				return;
   887  			}
   888  			if (c == '0') {
   889  				clear();
   890  				*pb = true;
   891  				return;
   892  			}
   893  			x.load(pb, is, ioMode); if (!*pb) return;
   894  			if (c == '1') {
   895  				y.load(pb, is, ioMode); if (!*pb) return;
   896  				if (!isValid(x, y)) {
   897  					*pb = false;
   898  					return;
   899  				}
   900  			} else if (c == '2' || c == '3') {
   901  				bool isYodd = c == '3';
   902  				*pb = getYfromX(y, x, isYodd);
   903  				if (!*pb) return;
   904  			} else if (c == '4') {
   905  				y.load(pb, is, ioMode); if (!*pb) return;
   906  #ifndef MCL_EC_USE_AFFINE
   907  				z.load(pb, is, ioMode); if (!*pb) return;
   908  #endif
   909  			} else {
   910  				*pb = false;
   911  				return;
   912  			}
   913  		}
   914  		if (verifyOrder_ && !isValidOrder()) {
   915  			*pb = false;
   916  		} else {
   917  			*pb = true;
   918  		}
   919  	}
   920  	// deplicated
   921  	static void setCompressedExpression(bool compressedExpression = true)
   922  	{
   923  		if (compressedExpression) {
   924  			ioMode_ |= IoEcCompY;
   925  		} else {
   926  			ioMode_ &= ~IoEcCompY;
   927  		}
   928  	}
   929  	/*
   930  		set IoMode for operator<<(), or operator>>()
   931  	*/
   932  	static void setIoMode(int ioMode)
   933  	{
   934  		assert(!(ioMode & 0xff));
   935  		ioMode_ = ioMode;
   936  	}
   937  	static inline int getIoMode() { return Fp::BaseFp::getIoMode() | ioMode_; }
   938  	static inline void getWeierstrass(Fp& yy, const Fp& x)
   939  	{
   940  		Fp t;
   941  		Fp::sqr(t, x);
   942  		t += a_;
   943  		t *= x;
   944  		Fp::add(yy, t, b_);
   945  	}
   946  	static inline bool getYfromX(Fp& y, const Fp& x, bool isYodd)
   947  	{
   948  		getWeierstrass(y, x);
   949  		if (!Fp::squareRoot(y, y)) {
   950  			return false;
   951  		}
   952  		if (y.isOdd() ^ isYodd) {
   953  			Fp::neg(y, y);
   954  		}
   955  		return true;
   956  	}
   957  	inline friend EcT operator+(const EcT& x, const EcT& y) { EcT z; add(z, x, y); return z; }
   958  	inline friend EcT operator-(const EcT& x, const EcT& y) { EcT z; sub(z, x, y); return z; }
   959  	template<class INT>
   960  	inline friend EcT operator*(const EcT& x, const INT& y) { EcT z; mul(z, x, y); return z; }
   961  	EcT& operator+=(const EcT& x) { add(*this, *this, x); return *this; }
   962  	EcT& operator-=(const EcT& x) { sub(*this, *this, x); return *this; }
   963  	template<class INT>
   964  	EcT& operator*=(const INT& x) { mul(*this, *this, x); return *this; }
   965  	EcT operator-() const { EcT x; neg(x, *this); return x; }
   966  	bool operator==(const EcT& rhs) const
   967  	{
   968  		EcT R;
   969  		sub(R, *this, rhs); // QQQ : optimized later
   970  		return R.isZero();
   971  	}
   972  	bool operator!=(const EcT& rhs) const { return !operator==(rhs); }
   973  	bool operator<(const EcT& rhs) const
   974  	{
   975  		return compare(*this, rhs) < 0;
   976  	}
   977  	bool operator>=(const EcT& rhs) const { return !operator<(rhs); }
   978  	bool operator>(const EcT& rhs) const { return rhs < *this; }
   979  	bool operator<=(const EcT& rhs) const { return !operator>(rhs); }
   980  	static inline void mulArray(EcT& z, const EcT& x, const fp::Unit *y, size_t yn, bool isNegative, bool constTime = false)
   981  	{
   982  		if (!constTime && x.isZero()) {
   983  			z.clear();
   984  			return;
   985  		}
   986  		if (mulArrayGLV && (constTime || yn > 1)) {
   987  			mulArrayGLV(z, x, y, yn, isNegative, constTime);
   988  			return;
   989  		}
   990  		mulArrayBase(z, x, y, yn, isNegative, constTime);
   991  	}
   992  	static inline void mulArrayBase(EcT& z, const EcT& x, const fp::Unit *y, size_t yn, bool isNegative, bool constTime)
   993  	{
   994  		EcT tmp;
   995  		const EcT *px = &x;
   996  		if (&z == &x) {
   997  			tmp = x;
   998  			px = &tmp;
   999  		}
  1000  		z.clear();
  1001  		fp::powGeneric(z, *px, y, yn, EcT::add, EcT::dbl, EcT::normalize, constTime ? Fp::BaseFp::getBitSize() : 0);
  1002  		if (isNegative) {
  1003  			neg(z, z);
  1004  		}
  1005  	}
  1006  	/*
  1007  		generic mul
  1008  	*/
  1009  	static inline void mulGeneric(EcT& z, const EcT& x, const mpz_class& y, bool constTime = false)
  1010  	{
  1011  		mulArrayBase(z, x, gmp::getUnit(y), gmp::getUnitSize(y), y < 0, constTime);
  1012  	}
  1013  #ifndef CYBOZU_DONT_USE_EXCEPTION
  1014  	static inline void init(const std::string& astr, const std::string& bstr, int mode = ec::Jacobi)
  1015  	{
  1016  		bool b;
  1017  		init(&b, astr.c_str(), bstr.c_str(), mode);
  1018  		if (!b) throw cybozu::Exception("mcl:EcT:init");
  1019  	}
  1020  	void set(const Fp& _x, const Fp& _y, bool verify = true)
  1021  	{
  1022  		bool b;
  1023  		set(&b, _x, _y, verify);
  1024  		if (!b) throw cybozu::Exception("ec:EcT:set") << _x << _y;
  1025  	}
  1026  	template<class OutputStream>
  1027  	void save(OutputStream& os, int ioMode = IoSerialize) const
  1028  	{
  1029  		bool b;
  1030  		save(&b, os, ioMode);
  1031  		if (!b) throw cybozu::Exception("EcT:save");
  1032  	}
  1033  	template<class InputStream>
  1034  	void load(InputStream& is, int ioMode = IoSerialize)
  1035  	{
  1036  		bool b;
  1037  		load(&b, is, ioMode);
  1038  		if (!b) throw cybozu::Exception("EcT:load");
  1039  	}
  1040  #endif
  1041  #ifndef CYBOZU_DONT_USE_STRING
  1042  	// backward compatilibity
  1043  	static inline void setParam(const std::string& astr, const std::string& bstr, int mode = ec::Jacobi)
  1044  	{
  1045  		init(astr, bstr, mode);
  1046  	}
  1047  	friend inline std::istream& operator>>(std::istream& is, EcT& self)
  1048  	{
  1049  		self.load(is, fp::detectIoMode(getIoMode(), is));
  1050  		return is;
  1051  	}
  1052  	friend inline std::ostream& operator<<(std::ostream& os, const EcT& self)
  1053  	{
  1054  		self.save(os, fp::detectIoMode(getIoMode(), os));
  1055  		return os;
  1056  	}
  1057  #endif
  1058  };
  1059  
  1060  template<class Fp> Fp EcT<Fp>::a_;
  1061  template<class Fp> Fp EcT<Fp>::b_;
  1062  template<class Fp> int EcT<Fp>::specialA_;
  1063  template<class Fp> int EcT<Fp>::ioMode_;
  1064  template<class Fp> bool EcT<Fp>::verifyOrder_;
  1065  template<class Fp> mpz_class EcT<Fp>::order_;
  1066  template<class Fp> void (*EcT<Fp>::mulArrayGLV)(EcT& z, const EcT& x, const fp::Unit *y, size_t yn, bool isNegative, bool constTime);
  1067  #ifndef MCL_EC_USE_AFFINE
  1068  template<class Fp> int EcT<Fp>::mode_;
  1069  #endif
  1070  
  1071  struct EcParam {
  1072  	const char *name;
  1073  	const char *p;
  1074  	const char *a;
  1075  	const char *b;
  1076  	const char *gx;
  1077  	const char *gy;
  1078  	const char *n;
  1079  	size_t bitSize; // bit length of p
  1080  	int curveType;
  1081  };
  1082  
  1083  } // mcl
  1084  
  1085  #ifdef CYBOZU_USE_BOOST
  1086  namespace mcl {
  1087  template<class Fp>
  1088  size_t hash_value(const mcl::EcT<Fp>& P_)
  1089  {
  1090  	if (P_.isZero()) return 0;
  1091  	mcl::EcT<Fp> P(P_); P.normalize();
  1092  	return mcl::hash_value(P.y, mcl::hash_value(P.x));
  1093  }
  1094  
  1095  }
  1096  #else
  1097  namespace std { CYBOZU_NAMESPACE_TR1_BEGIN
  1098  
  1099  template<class Fp>
  1100  struct hash<mcl::EcT<Fp> > {
  1101  	size_t operator()(const mcl::EcT<Fp>& P_) const
  1102  	{
  1103  		if (P_.isZero()) return 0;
  1104  		mcl::EcT<Fp> P(P_); P.normalize();
  1105  		return hash<Fp>()(P.y, hash<Fp>()(P.x));
  1106  	}
  1107  };
  1108  
  1109  CYBOZU_NAMESPACE_TR1_END } // std
  1110  #endif
  1111  
  1112  #ifdef _MSC_VER
  1113  	#pragma warning(pop)
  1114  #endif