贴点ECC各种变形的型式：

Genus-1 curves over large-characteristic fields  一般域的：
Doubling-oriented Doche–Icart–Kohel curves: y^2=x^3+a*x^2+16*a*x
Tripling-oriented Doche–Icart–Kohel curves: y^2=x^3+3*a*(x+1)^2
Edwards curves: x^2+y^2=c^2*(1+d*x^2*y^2)
Hessian curves: x^3+y^3+1=3*d*x*y
Jacobi intersections: s^2+c^2=1, a*s^2+d^2=1
Jacobi quartics: y^2=x^4+2*a*x^2+1
Montgomery curves: b*y^2=x^3+a*x^2+x
Short Weierstrass curves: y^2=x^3+a*x+b
Twisted Edwards curves: a*x^2+y^2=1+d*x^2*y^2
Twisted Hessian curves: a*x^3+y^3+1=d*x*y==================================
1

NIST推荐的最常见的短式魏尔斯特拉斯，图就最常见了，为了运算快，可混和用各种雅格比座标，仿射座标：
Short Weierstrass curves: y^2=x^3+a*x+b
Jacobian coordinates with a4=0: x=X/Z^2, y=Y/Z^3, a=0
Jacobian coordinates with a4=-3: x=X/Z^2, y=Y/Z^3, a=-3
Jacobian coordinates: x=X/Z^2, y=Y/Z^3
Modified Jacobian coordinates: x=X/Z^2, y=Y/Z^3, T=a*Z^4
Projective coordinates with a4=-1: x=X/Z, y=Y/Z, a=-1
Projective coordinates with a4=-3: x=X/Z, y=Y/Z, a=-3
Projective coordinates: x=X/Z, y=Y/Z
W12 coordinates with a6=0: x=X/Z, y=Y/Z^2, b=0
XYZZ coordinates with a4=-3: x=X/ZZ, y=Y/ZZZ, ZZ^3=ZZZ^2, a=-3
XYZZ coordinates: x=X/ZZ, y=Y/ZZZ, ZZ^3=ZZZ^2
XZ coordinates: x=X/Z

y^2=x^3+a*x+b

Affine addition formulas: (x1,y1)+(x2,y2)=(x3,y3) where
  x3 = (y2-y1)^2/(x2-x1)^2-x1-x2
  y3 = (2*x1+x2)*(y2-y1)/(x2-x1)-(y2-y1)^3/(x2-x1)^3-y1

Affine doubling formulas: 2(x1,y1)=(x3,y3) where
  x3 = (3*x1^2+a)^2/(2*y1)^2-x1-x1
  y3 = (2*x1+x1)*(3*x1^2+a)/(2*y1)-(3*x1^2+a)^3/(2*y1)^3-y1

Affine negation formulas: -(x1,y1)=(x1,-y1).
The neutral element of the curve is the unique point at infinity, namely (0:1:0) in projective coordinates.

Jacobian coordinates with a4=0 for short Weierstrass curves
name Jacobian coordinates with a4=0
assume a = 0
variable X
variable Y
variable Z
satisfying x = X/Z^2
satisfying y = Y/Z^3

用最小2乘法时最佳组合：

Smallest multiplication counts assuming I=100M, S=1M, *param=0M, add=0M, *const=0M:
16M for addition: 11M+5S. 12M+4S. 12M+4S. 12M+4S.
11M for addition with Z2=1: 7M+4S. 8M+3S. 8M+3S. 8M+3S.
7M for addition with Z1=Z2: 5M+2S.
6M for addition with Z1=1 and Z2=1: 4M+2S.
14M for readdition: 10M+4S after 11M+5S. 11M+3S after 12M+4S. 11M+3S after 12M+4S. 11M+3S after 12M+4S.
11M for readdition with Z2=1: 7M+4S after 7M+4S. 8M+3S after 8M+3S. 8M+3S after 8M+3S. 8M+3S after 8M+3S.
7M for readdition with Z1=Z2: 5M+2S after 5M+2S.
6M for readdition with Z1=1 and Z2=1: 4M+2S after 4M+2S.
7M for doubling: 2M+5S.
6M for doubling with Z1=1: 1M+5S.
15M for tripling: 5M+10S. 8M+7S.
104M for scaling: 1I+3M+1S.
Smallest multiplication counts assuming I=100M, S=0.8M, *param=0M, add=0M, *const=0M:
15M for addition: 11M+5S.
10.2M for addition with Z2=1: 7M+4S.
6.6M for addition with Z1=Z2: 5M+2S.
5.6M for addition with Z1=1 and Z2=1: 4M+2S.
13.2M for readdition: 10M+4S after 11M+5S.
10.2M for readdition with Z2=1: 7M+4S after 7M+4S.
6.6M for readdition with Z1=Z2: 5M+2S after 5M+2S.
5.6M for readdition with Z1=1 and Z2=1: 4M+2S after 4M+2S.
6M for doubling: 2M+5S.
5M for doubling with Z1=1: 1M+5S. 13M for tripling: 5M+10S.
103.8M for scaling: 1I+3M+1S.
Smallest multiplication counts assuming I=100M, S=0.67M, *param=0M, add=0M, *const=0M:
14.35M for addition: 11M+5S.
9.68M for addition with Z2=1: 7M+4S.
6.34M for addition with Z1=Z2: 5M+2S.
5.34M for addition with Z1=1 and Z2=1: 4M+2S.
12.68M for readdition: 10M+4S after 11M+5S.
9.68M for readdition with Z2=1: 7M+4S after 7M+4S.
6.34M for readdition with Z1=Z2: 5M+2S after 5M+2S.
5.34M for readdition with Z1=1 and Z2=1: 4M+2S after 4M+2S.
5.35M for doubling: 2M+5S.
4.35M for doubling with Z1=1: 1M+5S.
11.7M for tripling: 5M+10S.
103.67M for scaling: 1I+3M+1S.

挑一个用SAGE运行如下：
Assumptions: Z1=1.
Cost: 1M + 5S + 7add + 1*8 + 3*2 + 1*3.
Source: 2007 Bernstein–Lange.
Explicit formulas:
      XX = X1^2
      YY = Y1^2
      YYYY = YY^2
      S = 2*((X1+YY)^2-XX-YYYY)
      M = 3*XX+a
      T = M^2-2*S
      X3 = T
      Y3 = M*(S-T)-8*YYYY
      Z3 = 2*Y1

def mynumerator(x):
  if parent(x) == R:
    return x
  return numerator(x)

class fastfrac:
  def __init__(self,top,bot=1):
    if parent(top) == ZZ or parent(top) == R:
      self.top = R(top)
      self.bot = R(bot)
    elif top.__class__ == fastfrac:
      self.top = top.top
      self.bot = top.bot * bot
    else:
      self.top = R(numerator(top))
      self.bot = R(denominator(top)) * bot
  def reduce(self):
    return fastfrac(self.top / self.bot)
  def sreduce(self):
    return fastfrac(I.reduce(self.top),I.reduce(self.bot))
  def iszero(self):
    return self.top in I and not (self.bot in I)
  def isdoublingzero(self):
    return self.top in J and not (self.bot in J)
  def __add__(self,other):
    if parent(other) == ZZ:
      return fastfrac(self.top + self.bot * other,self.bot)
    if other.__class__ == fastfrac:
      return fastfrac(self.top * other.bot + self.bot * other.top,self.bot * other.bot)
    return NotImplemented
  def __sub__(self,other):
    if parent(other) == ZZ:
      return fastfrac(self.top - self.bot * other,self.bot)
    if other.__class__ == fastfrac:
      return fastfrac(self.top * other.bot - self.bot * other.top,self.bot * other.bot)
    return NotImplemented
  def __neg__(self):
    return fastfrac(-self.top,self.bot)
  def __mul__(self,other):
    if parent(other) == ZZ:
      return fastfrac(self.top * other,self.bot)
    if other.__class__ == fastfrac:
      return fastfrac(self.top * other.top,self.bot * other.bot)
    return NotImplemented
  def __rmul__(self,other):
    return self.__mul__(other)
  def __div__(self,other):
    if parent(other) == ZZ:
      return fastfrac(self.top,self.bot * other)
    if other.__class__ == fastfrac:
      return fastfrac(self.top * other.bot,self.bot * other.top)
    return NotImplemented
  def __pow__(self,other):
    if parent(other) == ZZ:
      return fastfrac(self.top ^ other,self.bot ^ other)
    return NotImplemented

def isidentity(x):
  return x.iszero()

def isdoublingidentity(x):
  return x.isdoublingzero()

R.<ua,ub,ux1,uy1,uX1,uY1,uZ1> = PolynomialRing(QQ,7,order='invlex')
I = R.ideal([
  mynumerator((uy1^2)-(ux1^3+ua*ux1+ub))
, mynumerator((ux1)-(uX1/uZ1^2))
, mynumerator((uy1)-(uY1/uZ1^3))
, mynumerator((ua)-(0))
, mynumerator((uZ1)-(1))
])

ua = fastfrac(ua)
ub = fastfrac(ub)
ux1 = fastfrac(ux1)
uy1 = fastfrac(uy1)
uX1 = fastfrac(uX1)
uY1 = fastfrac(uY1)
uZ1 = fastfrac(uZ1)

uXX = ((uX1^2))
uYY = ((uY1^2))
uYYYY = ((uYY^2))
uS = ((fastfrac(2)*((uX1+uYY)^2-uXX-uYYYY)))
uM = ((fastfrac(3)*uXX+ua))
uT = ((uM^2-fastfrac(2)*uS))
uX3 = ((uT))
uY3 = ((uM*(uS-uT)-fastfrac(8)*uYYYY))
uZ3 = ((fastfrac(2)*uY1))

ux3 = (((fastfrac(3)*ux1^2+ua)^2/(fastfrac(2)*uy1)^2-ux1-ux1)).reduce()
uy3 = (((fastfrac(2)*ux1+ux1)*(fastfrac(3)*ux1^2+ua)/(fastfrac(2)*uy1)-(fastfrac(3)*ux1^2+ua)^3/(fastfrac(2)*uy1)^3-uy1)).reduce()

print isidentity((uy3^2)-(ux3^3+ua*ux3+ub))
print isidentity((ux3)-(uX3/uZ3^2))
print isidentity((uy3)-(uY3/uZ3^3))
True
True
True

上传的附件：

• half-weierstrass.png （1.48kb，28次下载）
• half-weierstrass_turtle.png （57.25kb，27次下载）
• 1.JPG （20.71kb，16次下载）

2

爱德沃德式的图很4脚海星：

elliptic curve in Edwards form 爱德沃德式：
x^2+y^2=c^2*(1+d*x^2*y^2)

今年国产英文版论文集就说了爱德沃德式

求点公式：

Affine addition formulas: (x1,y1)+(x2,y2)=(x3,y3) where
  x3 = (x1*y2+y1*x2)/(c*(1+d*x1*x2*y1*y2))
  y3 = (y1*y2-x1*x2)/(c*(1-d*x1*x2*y1*y2))

Affine doubling formulas: 2(x1,y1)=(x3,y3) where
  x3 = (x1*y1+y1*x1)/(c*(1+d*x1*x1*y1*y1))
  y3 = (y1*y1-x1*x1)/(c*(1-d*x1*x1*y1*y1))

Affine negation formulas: -(x1,y1)=(-x1,y1).

x^2+y^2=c^2*(1+d*x^2*y^2)是超奇异的，不算EC，但特点是 point (0,c). The point (0,-c) has order 2. The points (c,0) and (-c,0) have order 4.         

如用仿射座标：
x=Z/X
  y=Z/Y

Smallest multiplication counts assuming I=100M, S=1M, *param=0M, add=0M, *const=0M:
10M for addition: 9M+1S.
9M for addition with Z2=1: 8M+1S. 9M.
9M for addition with X2=1: 8M+1S.
7M for addition with Z1=1 and Z2=1: 7M.
10M for readdition: 9M+1S after 9M+1S.
9M for readdition with Z2=1: 8M+1S after 8M+1S. 9M after 9M.
9M for readdition with X2=1: 8M+1S after 8M+1S.
7M for readdition with Z1=1 and Z2=1: 7M after 7M.
7M for doubling: 3M+4S.
6M for doubling with Z1=1: 3M+3S.
13M for tripling: 9M+4S.
102M for scaling: 1I+2M.
Smallest multiplication counts assuming I=100M, S=0.8M, *param=0M, add=0M, *const=0M:
9.8M for addition: 9M+1S.
8.8M for addition with Z2=1: 8M+1S.
8.8M for addition with X2=1: 8M+1S.
7M for addition with Z1=1 and Z2=1: 7M.
9.8M for readdition: 9M+1S after 9M+1S.
8.8M for readdition with Z2=1: 8M+1S after 8M+1S.
8.8M for readdition with X2=1: 8M+1S after 8M+1S.
7M for readdition with Z1=1 and Z2=1: 7M after 7M.
6.2M for doubling: 3M+4S.
5.4M for doubling with Z1=1: 3M+3S.
12.2M for tripling: 9M+4S.
102M for scaling: 1I+2M.
Smallest multiplication counts assuming I=100M, S=0.67M, *param=0M, add=0M, *const=0M:
9.67M for addition: 9M+1S.
8.67M for addition with Z2=1: 8M+1S.
8.67M for addition with X2=1: 8M+1S.
7M for addition with Z1=1 and Z2=1: 7M.
9.67M for readdition: 9M+1S after 9M+1S.
8.67M for readdition with Z2=1: 8M+1S after 8M+1S.
8.67M for readdition with X2=1: 8M+1S after 8M+1S.
7M for readdition with Z1=1 and Z2=1: 7M after 7M.
5.68M for doubling: 3M+4S.
5.01M for doubling with Z1=1: 3M+3S. 11.68M for tripling: 9M+4S.
102M for scaling: 1I+2M.

挑一个5.01M for doubling with Z1=1: 3M+3SCost: 9M + 4S + 1*c + 1*d + 7add + 1*4 + 1*2.
Source: 2007 Bernstein–Lange.
Explicit formulas:
      XX = X1^2
      YY = Y1^2
      ZZ = (c*Z1)^2
      D = XX+YY
      DD = D^2
      E = 4*(D-d*ZZ)
      H = 2*D*(YY-XX)
      P = DD-XX*E
      Q = DD-YY*E
      X3 = (H+Q)*Q*X1
      Y3 = (H-P)*P*Y1
      Z3 = P*Q*Z1

SAGE运行：

def mynumerator(x):
  if parent(x) == R:
    return x
  return numerator(x)

class fastfrac:
  def __init__(self,top,bot=1):
    if parent(top) == ZZ or parent(top) == R:
      self.top = R(top)
      self.bot = R(bot)
    elif top.__class__ == fastfrac:
      self.top = top.top
      self.bot = top.bot * bot
    else:
      self.top = R(numerator(top))
      self.bot = R(denominator(top)) * bot
  def reduce(self):
    return fastfrac(self.top / self.bot)
  def sreduce(self):
    return fastfrac(I.reduce(self.top),I.reduce(self.bot))
  def iszero(self):
    return self.top in I and not (self.bot in I)
  def isdoublingzero(self):
    return self.top in J and not (self.bot in J)
  def __add__(self,other):
    if parent(other) == ZZ:
      return fastfrac(self.top + self.bot * other,self.bot)
    if other.__class__ == fastfrac:
      return fastfrac(self.top * other.bot + self.bot * other.top,self.bot * other.bot)
    return NotImplemented
  def __sub__(self,other):
    if parent(other) == ZZ:
      return fastfrac(self.top - self.bot * other,self.bot)
    if other.__class__ == fastfrac:
      return fastfrac(self.top * other.bot - self.bot * other.top,self.bot * other.bot)
    return NotImplemented
  def __neg__(self):
    return fastfrac(-self.top,self.bot)
  def __mul__(self,other):
    if parent(other) == ZZ:
      return fastfrac(self.top * other,self.bot)
    if other.__class__ == fastfrac:
      return fastfrac(self.top * other.top,self.bot * other.bot)
    return NotImplemented
  def __rmul__(self,other):
    return self.__mul__(other)
  def __div__(self,other):
    if parent(other) == ZZ:
      return fastfrac(self.top,self.bot * other)
    if other.__class__ == fastfrac:
      return fastfrac(self.top * other.bot,self.bot * other.top)
    return NotImplemented
  def __pow__(self,other):
    if parent(other) == ZZ:
      return fastfrac(self.top ^ other,self.bot ^ other)
    return NotImplemented

def isidentity(x):
  return x.iszero()

def isdoublingidentity(x):
  return x.isdoublingzero()

R.<uc,ud,ux1,uy1,uX1,uY1,uZ1> = PolynomialRing(QQ,7,order='invlex')
I = R.ideal([
  mynumerator((ux1^2+uy1^2)-(uc^2*(1+ud*ux1^2*uy1^2)))
, mynumerator((ux1)-(uZ1/uX1))
, mynumerator((uy1)-(uZ1/uY1))
])

uc = fastfrac(uc)
ud = fastfrac(ud)
ux1 = fastfrac(ux1)
uy1 = fastfrac(uy1)
uX1 = fastfrac(uX1)
uY1 = fastfrac(uY1)
uZ1 = fastfrac(uZ1)

uXX = ((uX1^2))
uYY = ((uY1^2))
uZZ = (((uc*uZ1)^2))
uD = ((uXX+uYY))
uDD = ((uD^2))
uE = ((fastfrac(4)*(uD-ud*uZZ)))
uH = ((fastfrac(2)*uD*(uYY-uXX)))
uP = ((uDD-uXX*uE))
uQ = ((uDD-uYY*uE))
uX3 = (((uH+uQ)*uQ*uX1))
uY3 = (((uH-uP)*uP*uY1))
uZ3 = ((uP*uQ*uZ1))

ux2 = (((ux1*uy1+uy1*ux1)/(uc*(fastfrac(1)+ud*ux1*ux1*uy1*uy1)))).reduce()
uy2 = (((uy1*uy1-ux1*ux1)/(uc*(fastfrac(1)-ud*ux1*ux1*uy1*uy1)))).reduce()
ux3 = (((ux1*uy2+uy1*ux2)/(uc*(fastfrac(1)+ud*ux1*ux2*uy1*uy2)))).reduce()
uy3 = (((uy1*uy2-ux1*ux2)/(uc*(fastfrac(1)-ud*ux1*ux2*uy1*uy2)))).reduce()

print isidentity((ux3^2+uy3^2)-(uc^2*(fastfrac(1)+ud*ux3^2*uy3^2)))
print isidentity((ux3)-(uZ3/uX3))
print isidentity((uy3)-(uZ3/uY3))

True
True
True

上传的附件：

• half-edwards.png （1.48kb，26次下载）
• half-edwards_star.png （13.65kb，26次下载）

3
海神式：Hessian curves: x^3+y^3+1=3*d*x*y

上传的附件：

• half-hessian.png （1.30kb，26次下载）
• half-hessian-ray.png （8.21kb，27次下载）

4

雅可比4次式
Jacobi quartics: y^2=x^4+2*a*x^2+1

上传的附件：

• half-quartic.png （1.86kb，27次下载）
• half-Jac-squid.png （11.88kb，26次下载）

对了解ECC各种变形的速度比拟：

The Edwards–Weierstrass race

上传的附件：

• half-race_3_subtitles.png （78.62kb，24次下载）
• half-race_5.png （63.01kb，25次下载）

Edwards最好

genus-1 curves over binary fields 迹1   2进制域

Binary Edwards curves: d1*(x+y)+d2*(x^2+y^2)=(x+x^2)*(y+y^2)
Hessian curves: x^3+y^3+1=3*d*x*y
Short Weierstrass curves: y^2+x*y=x^3+a2*x^2+a6

1
NIST 推荐就一般书上的了，Lopez-Dahab 两人用了一种新座标，   包扩K氏曲线
Short Weierstrass curves: y^2+x*y=x^3+a2*x^2+a6     

y^2+x*y=x^3+a2*x^2+a6

Affine addition formulas: (x1,y1)+(x2,y2)=(x3,y3) where
  x3 = ((y1+y2)/(x1+x2))^2+((y1+y2)/(x1+x2))+x1+x2+a2
  y3 = ((y1+y2)/(x1+x2))^3+(x2+a2+1)*((y1+y2)/(x1+x2))+x1+x2+a2+y1

Affine doubling formulas: 2(x1,y1)=(x3,y3) where
  x3 = (x1+y1/x1)^2+(x1+y1/x1)+a2
  y3 = (x1+y1/x1)^3+(x1+a2+1)*(x1+y1/x1)+a2+y1

Affine negation formulas: -(x1,y1)=(x1,x1+y1).

试a2=1时

SAGE：

def mynumerator(x):
  if parent(x) == R:
    return x
  return numerator(x)

class fastfrac:
  def __init__(self,top,bot=1):
    if parent(top) == ZZ or parent(top) == R:
      self.top = R(top)
      self.bot = R(bot)
    elif top.__class__ == fastfrac:
      self.top = top.top
      self.bot = top.bot * bot
    else:
      self.top = R(numerator(top))
      self.bot = R(denominator(top)) * bot
  def reduce(self):
    return fastfrac(self.top / self.bot)
  def sreduce(self):
    return fastfrac(I.reduce(self.top),I.reduce(self.bot))
  def iszero(self):
    return self.top in I and not (self.bot in I)
  def isdoublingzero(self):
    return self.top in J and not (self.bot in J)
  def __add__(self,other):
    if parent(other) == ZZ:
      return fastfrac(self.top + self.bot * other,self.bot)
    if other.__class__ == fastfrac:
      return fastfrac(self.top * other.bot + self.bot * other.top,self.bot * other.bot)
    return NotImplemented
  def __sub__(self,other):
    if parent(other) == ZZ:
      return fastfrac(self.top - self.bot * other,self.bot)
    if other.__class__ == fastfrac:
      return fastfrac(self.top * other.bot - self.bot * other.top,self.bot * other.bot)
    return NotImplemented
  def __neg__(self):
    return fastfrac(-self.top,self.bot)
  def __mul__(self,other):
    if parent(other) == ZZ:
      return fastfrac(self.top * other,self.bot)
    if other.__class__ == fastfrac:
      return fastfrac(self.top * other.top,self.bot * other.bot)
    return NotImplemented
  def __rmul__(self,other):
    return self.__mul__(other)
  def __div__(self,other):
    if parent(other) == ZZ:
      return fastfrac(self.top,self.bot * other)
    if other.__class__ == fastfrac:
      return fastfrac(self.top * other.bot,self.bot * other.top)
    return NotImplemented
  def __pow__(self,other):
    if parent(other) == ZZ:
      return fastfrac(self.top ^ other,self.bot ^ other)
    return NotImplemented

def isidentity(x):
  return x.iszero()

def isdoublingidentity(x):
  return x.isdoublingzero()

R.<ua2,ua6,ux1,uy1,ux2,uy2> = PolynomialRing(GF(2),6,order='invlex')
I = R.ideal([
  mynumerator((uy1^2+ux1*uy1)-(ux1^3+ua2*ux1^2+ua6))
, mynumerator((uy2^2+ux2*uy2)-(ux2^3+ua2*ux2^2+ua6))
])

J = I + R.ideal([0
, ux1-ux2
, uy1-uy2
])

ua2 = fastfrac(ua2)
ua6 = fastfrac(ua6)
ux1 = fastfrac(ux1)
uy1 = fastfrac(uy1)
ux2 = fastfrac(ux2)
uy2 = fastfrac(uy2)

ux3 = ((((uy1+uy2)/(ux1+ux2))^2+((uy1+uy2)/(ux1+ux2))+ux1+ux2+ua2)).reduce()
uy3 = ((((uy1+uy2)/(ux1+ux2))^3+(ux2+ua2+fastfrac(1))*((uy1+uy2)/(ux1+ux2))+ux1+ux2+ua2+uy1)).reduce()
ux4 = (((ux1+uy1/ux1)^2+(ux1+uy1/ux1)+ua2)).reduce()
uy4 = (((ux1+uy1/ux1)^3+(ux1+ua2+fastfrac(1))*(ux1+uy1/ux1)+ua2+uy1)).reduce()
a0 = fastfrac((fastfrac(1)))
a1 = fastfrac((fastfrac(1)))
a2 = fastfrac((ua2))
a3 = fastfrac((fastfrac(0)))
a4 = fastfrac((fastfrac(0)))
a6 = fastfrac((ua6))
wweierx1 = ((ux1)).reduce().sreduce()
wweiery1 = ((uy1)).reduce().sreduce()
print isidentity(a0*(wweiery1^2)+a1*(wweierx1*wweiery1)+a3*wweiery1-(((wweierx1+a2)*wweierx1+a4)*wweierx1+a6))
print isidentity(ux1-(wweierx1))
print isidentity(uy1-(wweiery1))
wweierx2 = ((ux2)).reduce().sreduce()
wweiery2 = ((uy2)).reduce().sreduce()
wweierx3 = ((ux3)).reduce().sreduce()
wweiery3 = ((uy3)).reduce().sreduce()
wweierx4 = ((ux4)).reduce().sreduce()
wweiery4 = ((uy4)).reduce().sreduce()
slope = ((wweiery2-wweiery1)/(wweierx2-wweierx1)).reduce().sreduce()
print isidentity(a0*slope^2+a1*slope-a2-wweierx1-wweierx2-wweierx3)
print isidentity(slope*(wweierx1-wweierx3)-wweiery1-a1*wweierx3-a3-wweiery3)
slope = ((fastfrac(3)*wweierx1^2+fastfrac(2)*a2*wweierx1+a4-a1*wweiery1)/(fastfrac(2)*a0*wweiery1+a1*wweierx1+a3)).reduce().sreduce()
print isidentity(a0*slope^2+a1*slope-a2-wweierx1-wweierx1-wweierx4)
print isidentity(slope*(wweierx1-wweierx4)-wweiery1-a1*wweierx4-a3-wweiery4)

True
True
True
True
True
True
True

2

Binary Edwards curves

d1*(x+y)+d2*(x^2+y^2)=(x+x^2)*(y+y^2)

Affine addition formulas: (x1,y1)+(x2,y2)=(x3,y3) where
  x3 = (d1*(x1+x2)+d2*(x1+y1)*(x2+y2)+(x1+x1^2)*(x2*(y1+y2+1)+y1*y2))/(d1+(x1+x1^2)*(x2+y2))
  y3 = (d1*(y1+y2)+d2*(x1+y1)*(x2+y2)+(y1+y1^2)*(y2*(x1+x2+1)+x1*x2))/(d1+(y1+y1^2)*(x2+y2))

Affine doubling formulas: 2(x1,y1)=(x3,y3) where
  x3 = (d1*(x1+x1)+d2*(x1+y1)*(x1+y1)+(x1+x1^2)*(x1*(y1+y1+1)+y1*y1))/(d1+(x1+x1^2)*(x1+y1))
  y3 = (d1*(y1+y1)+d2*(x1+y1)*(x1+y1)+(y1+y1^2)*(y1*(x1+x1+1)+x1*x1))/(d1+(y1+y1^2)*(x1+y1))

Affine negation formulas: -(x1,y1)=(y1,x1).

可以x+y=w/d1=d2/ x+y=W/Z等多种限制

Best operation counts
Smallest multiplication counts assuming I=10M, S=0M, *param=0M, add=0M, *const=0M:
11M for doubling: 1I+1M+2S.
13M for differential addition: 1I+3M+1S.
24M for differential addition and doubling: 2I+4M+3S.
0M for scaling: 0M.
Smallest multiplication counts assuming I=10M, S=0.2M, *param=0M, add=0M, *const=0M:
11.4M for doubling: 1I+1M+2S.
13.2M for differential addition: 1I+3M+1S.
24.6M for differential addition and doubling: 2I+4M+3S.
0M for scaling: 0M.

挑一种看：
24.6M for differential addition and doubling: 2I+4M+3S.
Assumptions: d2overd1plus1=d2/d1+1.
Cost: 2I + 4M + 3S + 2*d2overd1plus1 + 9add.
Source: 2008 Bernstein–Lange–Rezaeian Farashahi.
Explicit formulas:
      R = w2*w3
      S = R^2
      T = R*(1+w2+w3)+S
      w5 = T/(d1+T+d2overd1plus1*S)+w1
      A = w2^2
      J = A^2
      K = A+J
      w4 = K/(d1+K+d2overd1plus1*J)

SAGE：

def mynumerator(x):
  if parent(x) == R:
    return x
  return numerator(x)

class fastfrac:
  def __init__(self,top,bot=1):
    if parent(top) == ZZ or parent(top) == R:
      self.top = R(top)
      self.bot = R(bot)
    elif top.__class__ == fastfrac:
      self.top = top.top
      self.bot = top.bot * bot
    else:
      self.top = R(numerator(top))
      self.bot = R(denominator(top)) * bot
  def reduce(self):
    return fastfrac(self.top / self.bot)
  def sreduce(self):
    return fastfrac(I.reduce(self.top),I.reduce(self.bot))
  def iszero(self):
    return self.top in I and not (self.bot in I)
  def isdoublingzero(self):
    return self.top in J and not (self.bot in J)
  def __add__(self,other):
    if parent(other) == ZZ:
      return fastfrac(self.top + self.bot * other,self.bot)
    if other.__class__ == fastfrac:
      return fastfrac(self.top * other.bot + self.bot * other.top,self.bot * other.bot)
    return NotImplemented
  def __sub__(self,other):
    if parent(other) == ZZ:
      return fastfrac(self.top - self.bot * other,self.bot)
    if other.__class__ == fastfrac:
      return fastfrac(self.top * other.bot - self.bot * other.top,self.bot * other.bot)
    return NotImplemented
  def __neg__(self):
    return fastfrac(-self.top,self.bot)
  def __mul__(self,other):
    if parent(other) == ZZ:
      return fastfrac(self.top * other,self.bot)
    if other.__class__ == fastfrac:
      return fastfrac(self.top * other.top,self.bot * other.bot)
    return NotImplemented
  def __rmul__(self,other):
    return self.__mul__(other)
  def __div__(self,other):
    if parent(other) == ZZ:
      return fastfrac(self.top,self.bot * other)
    if other.__class__ == fastfrac:
      return fastfrac(self.top * other.bot,self.bot * other.top)
    return NotImplemented
  def __pow__(self,other):
    if parent(other) == ZZ:
      return fastfrac(self.top ^ other,self.bot ^ other)
    return NotImplemented

def isidentity(x):
  return x.iszero()

def isdoublingidentity(x):
  return x.isdoublingzero()

R.<ud2overd1plus1,ud1,ud2,ux3,uy3,ux2,uy2,uw1,uw2,uw3> = PolynomialRing(GF(2),10,order='invlex')
ux_2 = (uy2)
uy_2 = (ux2)
ux1 = ((ud1*(ux3+ux_2)+ud2*(ux3+uy3)*(ux_2+uy_2)+(ux3+ux3^2)*(ux_2*(uy3+uy_2+1)+uy3*uy_2))/(ud1+(ux3+ux3^2)*(ux_2+uy_2)))
uy1 = ((ud1*(uy3+uy_2)+ud2*(ux3+uy3)*(ux_2+uy_2)+(uy3+uy3^2)*(uy_2*(ux3+ux_2+1)+ux3*ux_2))/(ud1+(uy3+uy3^2)*(ux_2+uy_2)))
I = R.ideal([
  mynumerator((ud1*(ux1+uy1)+ud2*(ux1^2+uy1^2))-((ux1+ux1^2)*(uy1+uy1^2)))
, mynumerator((ux1+uy1)-(uw1))
, mynumerator((ud1*(ux2+uy2)+ud2*(ux2^2+uy2^2))-((ux2+ux2^2)*(uy2+uy2^2)))
, mynumerator((ux2+uy2)-(uw2))
, mynumerator((ud1*(ux3+uy3)+ud2*(ux3^2+uy3^2))-((ux3+ux3^2)*(uy3+uy3^2)))
, mynumerator((ux3+uy3)-(uw3))
, mynumerator((ud2overd1plus1)-(ud2/ud1+1))
])

ud2overd1plus1 = fastfrac(ud2overd1plus1)
ud1 = fastfrac(ud1)
ud2 = fastfrac(ud2)
ux3 = fastfrac(ux3)
uy3 = fastfrac(uy3)
ux2 = fastfrac(ux2)
uy2 = fastfrac(uy2)
uw1 = fastfrac(uw1)
uw2 = fastfrac(uw2)
uw3 = fastfrac(uw3)
ux_2 = fastfrac(ux_2)
uy_2 = fastfrac(uy_2)
ux1 = fastfrac(ux1)
uy1 = fastfrac(uy1)

uR = ((uw2*uw3))
uS = ((uR^2))
uT = ((uR*(fastfrac(1)+uw2+uw3)+uS))
uw5 = ((uT/(ud1+uT+ud2overd1plus1*uS)+uw1))
uA = ((uw2^2))
uJ = ((uA^2))
uK = ((uA+uJ))
uw4 = ((uK/(ud1+uK+ud2overd1plus1*uJ)))

ux4 = (((ud1*(ux2+ux2)+ud2*(ux2+uy2)*(ux2+uy2)+(ux2+ux2^2)*(ux2*(uy2+uy2+fastfrac(1))+uy2*uy2))/(ud1+(ux2+ux2^2)*(ux2+uy2)))).reduce()
ux5 = (((ud1*(ux3+ux2)+ud2*(ux3+uy3)*(ux2+uy2)+(ux3+ux3^2)*(ux2*(uy3+uy2+fastfrac(1))+uy3*uy2))/(ud1+(ux3+ux3^2)*(ux2+uy2)))).reduce()
uy4 = (((ud1*(uy2+uy2)+ud2*(ux2+uy2)*(ux2+uy2)+(uy2+uy2^2)*(uy2*(ux2+ux2+fastfrac(1))+ux2*ux2))/(ud1+(uy2+uy2^2)*(ux2+uy2)))).reduce()
uy5 = (((ud1*(uy3+uy2)+ud2*(ux3+uy3)*(ux2+uy2)+(uy3+uy3^2)*(uy2*(ux3+ux2+fastfrac(1))+ux3*ux2))/(ud1+(uy3+uy3^2)*(ux2+uy2)))).reduce()

print isidentity((ud1*(ux5+uy5)+ud2*(ux5^2+uy5^2))-((ux5+ux5^2)*(uy5+uy5^2))) or isidentity(uy1*uy1*uy2*uy3*ux1*ux1*ux2*ux3*((ud1*(ux5+uy5)+ud2*(ux5^2+uy5^2))-((ux5+ux5^2)*(uy5+uy5^2))))
print isidentity((ux5+uy5)-(uw5)) or isidentity(uy1*uy1*uy2*uy3*ux1*ux1*ux2*ux3*((ux5+uy5)-(uw5)))
print isidentity((ud1*(ux5+uy5)+ud2*(ux5^2+uy5^2))-((ux5+ux5^2)*(uy5+uy5^2))) or isidentity(uy1*uy1*uy2*uy3*ux1*ux1*ux2*ux3*((ud1*(ux5+uy5)+ud2*(ux5^2+uy5^2))-((ux5+ux5^2)*(uy5+uy5^2))))
print isidentity((ux5+uy5)-(uw5)) or isidentity(uy1*uy1*uy2*uy3*ux1*ux1*ux2*ux3*((ux5+uy5)-(uw5)))

True
True
True
True

3
Hessian curves   就一种仿射变换

x^3+y^3+1=3*d*x*y

Affine addition formulas: (x1,y1)+(x2,y2)=(x3,y3) where
  x3 = (y1^2*x2-y2^2*x1)/(x2*y2-x1*y1)
  y3 = (x1^2*y2-x2^2*y1)/(x2*y2-x1*y1)

Affine doubling formulas: 2(x1,y1)=(x3,y3) where
  x3 = y1*(1-x1^3)/(x1^3-y1^3)
  y3 = x1*(y1^3-1)/(x1^3-y1^3)

Affine negation formulas: -(x1,y1)=(y1,x1).

x=X/Z
  y=Y/Z

Best operation counts
Smallest multiplication counts assuming I=10M, S=0M, *param=0M, add=0M, *const=0M:
12M for addition: 12M. 12M. 12M+6S.
10M for addition with Z2=1: 10M.
8M for addition with Z1=1 and Z2=1: 8M.
9M for readdition: 9M+3S after 12M+6S.
9M for readdition with Z2=1: 9M+3S after 12M+6S.
7M for readdition with Z1=1 and Z2=1: 7M after 8M.
6M for doubling: 6M+3S.
11M for tripling: 11M+4S.
12M for scaling: 1I+2M.
Smallest multiplication counts assuming I=10M, S=0.2M, *param=0M, add=0M, *const=0M:
12M for addition: 12M. 12M.
10M for addition with Z2=1: 10M.
8M for addition with Z1=1 and Z2=1: 8M.
9.6M for readdition: 9M+3S after 12M+6S.
9.6M for readdition with Z2=1: 9M+3S after 12M+6S.
7M for readdition with Z1=1 and Z2=1: 7M after 8M.
6.6M for doubling: 6M+3S.
11.8M for tripling: 11M+4S.
12M for scaling: 1I+2M.