Implement binary EEA inversion for faster BNADD precompile

core/vm/contracts.go

@@ -423,7 +423,7 @@ func runBn256Add(input []byte) ([]byte, error) {
 	}
 	res := new(bn256.G1)
 	res.Add(x, y)
-	return res.Marshal(), nil
+	return res.MarshalVariableTime(), nil
 }
 
 // bn256Add implements a native elliptic curve point addition conforming to
@@ -461,7 +461,7 @@ func runBn256ScalarMul(input []byte) ([]byte, error) {
 	}
 	res := new(bn256.G1)
 	res.ScalarMult(p, new(big.Int).SetBytes(getData(input, 64, 32)))
-	return res.Marshal(), nil
+	return res.MarshalVariableTime(), nil
 }
 
 // bn256ScalarMulIstanbul implements a native elliptic curve scalar

crypto/bn256/cloudflare/bn256.go

@@ -124,6 +124,32 @@ func (e *G1) Marshal() []byte {
 	return ret
 }
 
+// Marshal converts e to a byte slice.
+// Uses variable time algorithms for inversion
+// for transformation to affine coordinates
+func (e *G1) MarshalVariableTime() []byte {
+	// Each value is a 256-bit number.
+	const numBytes = 256 / 8
+
+	if e.p == nil {
+		e.p = &curvePoint{}
+	}
+
+	e.p.MakeAffineVariableTime()
+	ret := make([]byte, numBytes*2)
+	if e.p.IsInfinity() {
+		return ret
+	}
+	temp := &gfP{}
+
+	montDecode(temp, &e.p.x)
+	temp.Marshal(ret)
+	montDecode(temp, &e.p.y)
+	temp.Marshal(ret[numBytes:])
+
+	return ret
+}
+
 // Unmarshal sets e to the result of converting the output of Marshal back into
 // a group element and then returns e.
 func (e *G1) Unmarshal(m []byte) ([]byte, error) {

crypto/bn256/cloudflare/bn256_test.go

@@ -92,6 +92,71 @@ func TestTripartiteDiffieHellman(t *testing.T) {
 	}
 }
 
+func TestBinaryEAA(t *testing.T) {
+	for i := 0; i < 10000; i++ {
+		_, Ga, err := RandomG1(rand.Reader)
+		if err != nil {
+			t.Fatal(err)
+		}
+		tmpLittleFermat := &gfP{}
+		tmpLittleFermat.Invert(&Ga.p.x)
+
+		tmpBinaryEAA := &gfP{}
+		tmpBinaryEAA.InvertVariableTime(&Ga.p.x)
+
+		tmpBinaryEAASelfSet := &gfP{}
+		tmpBinaryEAASelfSet.Set(&Ga.p.x)
+		tmpBinaryEAASelfSet.InvertVariableTime(tmpBinaryEAASelfSet)
+
+		if *tmpLittleFermat != *tmpBinaryEAA {
+			t.Fatalf("results of different inversion do not agree")
+		}
+
+		if *tmpLittleFermat != *tmpBinaryEAASelfSet {
+			t.Fatalf("self-assigned inversion is invalid")
+		}
+	}
+}
+
+func BenchmarkLittleFermatInversion(b *testing.B) {
+	el := gfP{0x0, 0x97816a916871ca8d, 0xb85045b68181585d, 0x30644e72e131a029}
+
+	b.ResetTimer()
+
+	tmp := &gfP{}
+	for i := 0; i < b.N; i++ {
+		tmp.Invert(&el)
+	}
+}
+
+func BenchmarkBinaryEEAInversion(b *testing.B) {
+	el := gfP{0x0, 0x97816a916871ca8d, 0xb85045b68181585d, 0x30644e72e131a029}
+
+	b.ResetTimer()
+
+	tmp := &gfP{}
+	for i := 0; i < b.N; i++ {
+		tmp.InvertVariableTime(&el)
+	}
+}
+
+func BenchmarkG1AddAndMakeAffine(b *testing.B) {
+	_, Ga, err := RandomG1(rand.Reader)
+	if err != nil {
+		b.Fatal(err)
+	}
+	_, Gb, err := RandomG1(rand.Reader)
+	if err != nil {
+		b.Fatal(err)
+	}
+	b.ResetTimer()
+
+	for i := 0; i < b.N; i++ {
+		e := new(G1).Add(Ga, Gb)
+		e.p.MakeAffine()
+	}
+}
+
 func TestG2SelfAddition(t *testing.T) {
 	s, _ := rand.Int(rand.Reader, Order)
 	p := new(G2).ScalarBaseMult(s)

crypto/bn256/cloudflare/curve.go

@@ -206,6 +206,30 @@ func (c *curvePoint) Mul(a *curvePoint, scalar *big.Int) {
 	c.Set(sum)
 }
 
+func (c *curvePoint) MakeAffineVariableTime() {
+	if c.z == *newGFp(1) {
+		return
+	} else if c.z == *newGFp(0) {
+		c.x = gfP{0}
+		c.y = *newGFp(1)
+		c.t = gfP{0}
+		return
+	}
+
+	zInv := &gfP{}
+	zInv.InvertVariableTime(&c.z)
+
+	t, zInv2 := &gfP{}, &gfP{}
+	gfpMul(t, &c.y, zInv)
+	gfpMul(zInv2, zInv, zInv)
+
+	gfpMul(&c.x, &c.x, zInv2)
+	gfpMul(&c.y, t, zInv2)
+
+	c.z = *newGFp(1)
+	c.t = *newGFp(1)
+}
+
 func (c *curvePoint) MakeAffine() {
 	if c.z == *newGFp(1) {
 		return

crypto/bn256/cloudflare/gfp.go

@@ -3,6 +3,7 @@ package bn256
 import (
 	"errors"
 	"fmt"
+	"math/bits"
 )
 
 type gfP [4]uint64
@@ -80,3 +81,110 @@ func (e *gfP) Unmarshal(in []byte) error {
 
 func montEncode(c, a *gfP) { gfpMul(c, a, r2) }
 func montDecode(c, a *gfP) { gfpMul(c, a, &gfP{1}) }
+
+func isZero(a *gfP) bool {
+	return (a[0] | a[1] | a[2] | a[3]) == 0
+}
+
+func isEven(a *gfP) bool {
+	return a[0]&1 == 0
+}
+
+func div2(a *gfP) {
+	a[0] = a[0]>>1 | a[1]<<63
+	a[1] = a[1]>>1 | a[2]<<63
+	a[2] = a[2]>>1 | a[3]<<63
+	a[3] = a[3] >> 1
+}
+
+func (e *gfP) addNocarry(f *gfP) {
+	carry := uint64(0)
+	e[0], carry = bits.Add64(e[0], f[0], carry)
+	e[1], carry = bits.Add64(e[1], f[1], carry)
+	e[2], carry = bits.Add64(e[2], f[2], carry)
+	e[3], _ = bits.Add64(e[3], f[3], carry)
+}
+
+func (e *gfP) subNoborrow(f *gfP) {
+	borrow := uint64(0)
+	e[0], borrow = bits.Sub64(e[0], f[0], borrow)
+	e[1], borrow = bits.Sub64(e[1], f[1], borrow)
+	e[2], borrow = bits.Sub64(e[2], f[2], borrow)
+	e[3], _ = bits.Sub64(e[3], f[3], borrow)
+}
+
+func gte(a, b *gfP) bool {
+	// subtract b from a. If no borrow occures then a >= b
+	borrow := uint64(0)
+	_, borrow = bits.Sub64(a[0], b[0], borrow)
+	_, borrow = bits.Sub64(a[1], b[1], borrow)
+	_, borrow = bits.Sub64(a[2], b[2], borrow)
+	_, borrow = bits.Sub64(a[3], b[3], borrow)
+
+	return borrow == 0
+}
+
+// Performs inversion of the field element using binary EEA.
+// If element is zero (no inverse exists) then set `e` to zero
+func (e *gfP) InvertVariableTime(f *gfP) {
+	if isZero(f) {
+		e.Set(&gfP{0, 0, 0, 0})
+		return
+	}
+
+	// Guajardo Kumar Paar Pelzl
+	// Efficient Software-Implementation of Finite Fields with Applications to Cryptography
+	// Algorithm 16 (BEA for Inversion in Fp)
+
+	one := gfP{1, 0, 0, 0}
+
+	u, b := gfP{}, gfP{}
+	u.Set(f)
+	b.Set(r2)
+
+	v := gfP(p2)
+	c := gfP{0, 0, 0, 0}
+	modulus := gfP(p2)
+
+	for u != one && v != one {
+		// while u is even
+		for isEven(&u) {
+			div2(&u)
+			if !isEven(&b) {
+				// we will not overflow a modulus here,
+				// so we can use specialized function
+				// do perform addition without reduction
+				b.addNocarry(&modulus)
+			}
+			div2(&b)
+		}
+
+		// while v is even
+		for isEven(&v) {
+			div2(&v)
+			if !isEven(&c) {
+				// we will not overflow a modulus here,
+				// so we can use specialized function
+				// do perform addition without reduction
+				c.addNocarry(&modulus)
+			} else {
+
+			}
+			div2(&c)
+		}
+
+		if gte(&u, &v) {
+			u.subNoborrow(&v)
+			gfpSub(&b, &b, &c)
+		} else {
+			v.subNoborrow(&u)
+			gfpSub(&c, &c, &b)
+		}
+	}
+
+	if u == one {
+		e.Set(&b)
+	} else {
+		e.Set(&c)
+	}
+}

crypto/bn256/cloudflare/gfp12.go

@@ -143,6 +143,22 @@ func (e *gfP12) Square(a *gfP12) *gfP12 {
 	return e
 }
 
+func (e *gfP12) InvertVariableTime(a *gfP12) *gfP12 {
+	// See "Implementing cryptographic pairings", M. Scott, section 3.2.
+	// ftp://136.206.11.249/pub/crypto/pairings.pdf
+	t1, t2 := &gfP6{}, &gfP6{}
+
+	t1.Square(&a.x)
+	t2.Square(&a.y)
+	t1.MulTau(t1).Sub(t2, t1)
+	t2.InvertVariableTime(t1)
+
+	e.x.Neg(&a.x)
+	e.y.Set(&a.y)
+	e.MulScalar(e, t2)
+	return e
+}
+
 func (e *gfP12) Invert(a *gfP12) *gfP12 {
 	// See "Implementing cryptographic pairings", M. Scott, section 3.2.
 	// ftp://136.206.11.249/pub/crypto/pairings.pdf

crypto/bn256/cloudflare/gfp2.go

@@ -137,6 +137,24 @@ func (e *gfP2) Square(a *gfP2) *gfP2 {
 	return e
 }
 
+func (e *gfP2) InvertVariableTime(a *gfP2) *gfP2 {
+	// See "Implementing cryptographic pairings", M. Scott, section 3.2.
+	// ftp://136.206.11.249/pub/crypto/pairings.pdf
+	t1, t2 := &gfP{}, &gfP{}
+	gfpMul(t1, &a.x, &a.x)
+	gfpMul(t2, &a.y, &a.y)
+	gfpAdd(t1, t1, t2)
+
+	inv := &gfP{}
+	inv.InvertVariableTime(t1)
+
+	gfpNeg(t1, &a.x)
+
+	gfpMul(&e.x, t1, inv)
+	gfpMul(&e.y, &a.y, inv)
+	return e
+}
+
 func (e *gfP2) Invert(a *gfP2) *gfP2 {
 	// See "Implementing cryptographic pairings", M. Scott, section 3.2.
 	// ftp://136.206.11.249/pub/crypto/pairings.pdf

crypto/bn256/cloudflare/gfp6.go

@@ -211,3 +211,47 @@ func (e *gfP6) Invert(a *gfP6) *gfP6 {
 	e.z.Mul(A, F)
 	return e
 }
+
+func (e *gfP6) InvertVariableTime(a *gfP6) *gfP6 {
+	// See "Implementing cryptographic pairings", M. Scott, section 3.2.
+	// ftp://136.206.11.249/pub/crypto/pairings.pdf
+
+	// Here we can give a short explanation of how it works: let j be a cubic root of
+	// unity in GF(p²) so that 1+j+j²=0.
+	// Then (xτ² + yτ + z)(xj²τ² + yjτ + z)(xjτ² + yj²τ + z)
+	// = (xτ² + yτ + z)(Cτ²+Bτ+A)
+	// = (x³ξ²+y³ξ+z³-3ξxyz) = F is an element of the base field (the norm).
+	//
+	// On the other hand (xj²τ² + yjτ + z)(xjτ² + yj²τ + z)
+	// = τ²(y²-ξxz) + τ(ξx²-yz) + (z²-ξxy)
+	//
+	// So that's why A = (z²-ξxy), B = (ξx²-yz), C = (y²-ξxz)
+	t1 := (&gfP2{}).Mul(&a.x, &a.y)
+	t1.MulXi(t1)
+
+	A := (&gfP2{}).Square(&a.z)
+	A.Sub(A, t1)
+
+	B := (&gfP2{}).Square(&a.x)
+	B.MulXi(B)
+	t1.Mul(&a.y, &a.z)
+	B.Sub(B, t1)
+
+	C := (&gfP2{}).Square(&a.y)
+	t1.Mul(&a.x, &a.z)
+	C.Sub(C, t1)
+
+	F := (&gfP2{}).Mul(C, &a.y)
+	F.MulXi(F)
+	t1.Mul(A, &a.z)
+	F.Add(F, t1)
+	t1.Mul(B, &a.x).MulXi(t1)
+	F.Add(F, t1)
+
+	F.InvertVariableTime(F)
+
+	e.x.Mul(C, F)
+	e.y.Mul(B, F)
+	e.z.Mul(A, F)
+	return e
+}

crypto/bn256/cloudflare/twist.go

@@ -176,6 +176,26 @@ func (c *twistPoint) Mul(a *twistPoint, scalar *big.Int) {
 	c.Set(sum)
 }
 
+func (c *twistPoint) MakeAffineVariableTime() {
+	if c.z.IsOne() {
+		return
+	} else if c.z.IsZero() {
+		c.x.SetZero()
+		c.y.SetOne()
+		c.t.SetZero()
+		return
+	}
+
+	zInv := (&gfP2{}).InvertVariableTime(&c.z)
+	t := (&gfP2{}).Mul(&c.y, zInv)
+	zInv2 := (&gfP2{}).Square(zInv)
+	c.y.Mul(t, zInv2)
+	t.Mul(&c.x, zInv2)
+	c.x.Set(t)
+	c.z.SetOne()
+	c.t.SetOne()
+}
+
 func (c *twistPoint) MakeAffine() {
 	if c.z.IsOne() {
 		return