perf: emulated-curve scalar multiplication

internal/stats/latest_stats.csv

@@ -111,24 +111,24 @@ pairing_bn254,bn254,groth16,836345,1354248
 pairing_bn254,bn254,plonk,3025774,2907728
 pairing_bw6761,bn254,groth16,3242559,5348585
 pairing_bw6761,bn254,plonk,12295948,11867314
-scalar_mul_G1_bn254,bn254,groth16,167268,262016
-scalar_mul_G1_bn254,bn254,plonk,606389,584504
-scalar_mul_G1_bn254_incomplete,bn254,groth16,76581,124556
-scalar_mul_G1_bn254_incomplete,bn254,plonk,286253,277008
-scalar_mul_G2_bls12381,bn254,groth16,603963,955562
-scalar_mul_G2_bls12381,bn254,plonk,2116766,2029166
-scalar_mul_G2_bn254,bn254,groth16,452657,720753
-scalar_mul_G2_bn254,bn254,plonk,1564257,1499197
+scalar_mul_G1_bn254,bn254,groth16,165192,258804
+scalar_mul_G1_bn254,bn254,plonk,598772,577171
+scalar_mul_G1_bn254_incomplete,bn254,groth16,75656,123090
+scalar_mul_G1_bn254_incomplete,bn254,plonk,282712,273589
+scalar_mul_G2_bls12381,bn254,groth16,598959,947949
+scalar_mul_G2_bls12381,bn254,plonk,2098981,2012192
+scalar_mul_G2_bn254,bn254,groth16,449267,715630
+scalar_mul_G2_bn254,bn254,plonk,1552336,1487843
 scalar_mul_G2_bw6761,bn254,groth16,371608,592187
 scalar_mul_G2_bw6761,bn254,plonk,1330270,1280942
 scalar_mul_P256,bn254,groth16,124998,198925
 scalar_mul_P256,bn254,plonk,459293,443079
 scalar_mul_P256_incomplete,bn254,groth16,92548,148602
 scalar_mul_P256_incomplete,bn254,plonk,341765,329603
-scalar_mul_secp256k1,bn254,groth16,167336,262128
-scalar_mul_secp256k1,bn254,plonk,606674,584781
-scalar_mul_secp256k1_incomplete,bn254,groth16,76649,124668
-scalar_mul_secp256k1_incomplete,bn254,plonk,286535,277279
+scalar_mul_secp256k1,bn254,groth16,165260,258916
+scalar_mul_secp256k1,bn254,plonk,599057,577448
+scalar_mul_secp256k1_incomplete,bn254,groth16,75724,123202
+scalar_mul_secp256k1_incomplete,bn254,plonk,282994,273860
 selector/binaryMux_4,bn254,groth16,5,3
 selector/binaryMux_4,bls12_377,groth16,5,3
 selector/binaryMux_4,bls12_381,groth16,5,3

std/algebra/emulated/sw_bls12381/g1.go

@@ -237,14 +237,23 @@ func (g1 *G1) AssertIsOnG1(P *G1Affine) {
 	// 1- Check P is on the curve
 	g1.AssertIsOnCurve(P)
 
-	// 2- Check P has the right subgroup order
-	// [x²]ϕ(P)
-	phiP := g1.phi(P)
-	_P := g1.scalarMulBySeedSquare(phiP)
-	_P = g1.neg(_P)
-
-	// [r]Q == 0 <==>  P = -[x²]ϕ(P)
-	g1.AssertIsEqual(_P, P)
+	// 2- Check P is in the prime-order subgroup.
+	//
+	// We hint a preimage S and assert [x-1]S == P. With n = (x-1)/3 the
+	// cofactor torsion is Z_n × Z_{3n} (the n-part is rank 2 because the full
+	// n-torsion is rational, the factor 3 is cyclic), so its exponent is
+	// 3n = x-1 and [x-1]E(Fp) = G1 exactly: a point carrying cofactor torsion
+	// has no on-curve preimage. Completeness holds because gcd(x-1, r) = 1
+	// makes [x-1] a bijection on G1. See
+	// [sw_emulated.CurveParams.CofactorClearing].
+	//
+	// This replaces the endomorphism test P = -[x²]ϕ(P) (Bowe, eprint
+	// 2019/814), which needs a 128-bit ladder where this needs a 64-bit one.
+	curve, err := sw_emulated.New[BaseField, ScalarField](g1.api, sw_emulated.GetBLS12381Params())
+	if err != nil {
+		panic(fmt.Sprintf("new emulated curve: %v", err))
+	}
+	curve.AssertIsInSubgroup(P)
 }
 
 // AssertIsEqual asserts that p and q are the same point.

std/algebra/emulated/sw_bls12381/g1_test.go

@@ -2,11 +2,14 @@ package sw_bls12381
 
 import (
 	"fmt"
+	"math/big"
 	"testing"
 
 	"github.com/consensys/gnark-crypto/ecc"
 	bls12381 "github.com/consensys/gnark-crypto/ecc/bls12-381"
+	fp_bls12381 "github.com/consensys/gnark-crypto/ecc/bls12-381/fp"
 	"github.com/consensys/gnark/frontend"
+	"github.com/consensys/gnark/std/math/emulated"
 	"github.com/consensys/gnark/test"
 )
 
@@ -58,3 +61,100 @@ func TestTripleG1(t *testing.T) {
 	err = test.IsSolved(&tripleCircuitConsistency{}, &witness2, ecc.BN254.ScalarField())
 	assert.NoError(err)
 }
+
+type assertIsOnG1Circuit struct {
+	P G1Affine
+}
+
+func (c *assertIsOnG1Circuit) Define(api frontend.API) error {
+	g1, err := NewG1(api)
+	if err != nil {
+		return err
+	}
+	g1.AssertIsOnG1(&c.P)
+	return nil
+}
+
+func onG1Witness(p bls12381.G1Affine) *assertIsOnG1Circuit {
+	return &assertIsOnG1Circuit{P: G1Affine{
+		X: emulated.ValueOf[BaseField](p.X),
+		Y: emulated.ValueOf[BaseField](p.Y),
+	}}
+}
+
+// curvePointAtX returns the point of E(Fp): y² = x³ + 4 with the given small
+// integer x-coordinate and the canonical square root as y.
+//
+// The off-subgroup inputs here are built deterministically rather than sampled,
+// so a failure is reproducible and no candidate can be silently skipped for
+// landing in the subgroup. (sw_emulated's randomBLS12381CurvePoint samples
+// instead, because its check is statistical — it asserts [c]P lands in G1 over
+// many random P. The two are not interchangeable, so neither is a copy of the
+// other.)
+func curvePointAtX(t *testing.T, x uint64) bls12381.G1Affine {
+	t.Helper()
+	var xe, y2 fp_bls12381.Element
+	xe.SetUint64(x)
+	y2.Square(&xe).Mul(&y2, &xe).Add(&y2, new(fp_bls12381.Element).SetUint64(4))
+	if y2.Legendre() != 1 {
+		t.Fatalf("x = %d is not the x-coordinate of a curve point", x)
+	}
+	var y fp_bls12381.Element
+	y.Sqrt(&y2)
+	p := bls12381.G1Affine{X: xe, Y: y}
+	if !p.IsOnCurve() {
+		t.Fatalf("constructed point at x = %d is not on the curve", x)
+	}
+	return p
+}
+
+func TestAssertIsOnG1(t *testing.T) {
+	assert := test.NewAssert(t)
+	_, _, g, _ := bls12381.Generators()
+
+	// Completeness: genuine subgroup points, and the (0,0) infinity encoding.
+	for i := 0; i < 3; i++ {
+		var q bls12381.G1Affine
+		q.ScalarMultiplication(&g, big.NewInt(int64(4242421+i*104729)))
+		assert.True(q.IsInSubGroup())
+		assert.CheckCircuit(&assertIsOnG1Circuit{}, test.WithValidAssignment(onG1Witness(q)),
+			test.WithCurves(ecc.BN254), test.NoProverChecks())
+	}
+	var infinity bls12381.G1Affine
+	assert.CheckCircuit(&assertIsOnG1Circuit{}, test.WithValidAssignment(onG1Witness(infinity)),
+		test.WithCurves(ecc.BN254), test.NoProverChecks())
+
+	// Soundness: the rational 3-torsion point (0,2), a subgroup point shifted
+	// by it, and uniformly random off-subgroup curve points.
+	var t0 bls12381.G1Affine
+	t0.X.SetZero()
+	t0.Y.SetUint64(2)
+	assert.True(t0.IsOnCurve())
+	assert.False(t0.IsInSubGroup())
+	assert.CheckCircuit(&assertIsOnG1Circuit{}, test.WithInvalidAssignment(onG1Witness(t0)),
+		test.WithCurves(ecc.BN254), test.NoProverChecks())
+
+	var tainted bls12381.G1Affine
+	tainted.ScalarMultiplication(&g, big.NewInt(424242))
+	tainted.Add(&tainted, &t0)
+	assert.True(tainted.IsOnCurve())
+	assert.False(tainted.IsInSubGroup())
+	assert.CheckCircuit(&assertIsOnG1Circuit{}, test.WithInvalidAssignment(onG1Witness(tainted)),
+		test.WithCurves(ecc.BN254), test.NoProverChecks())
+
+	// Generic off-subgroup curve points, carrying torsion in both the n-part
+	// and the 3-part of the cofactor rather than only the 3-torsion above.
+	for _, x := range []uint64{4, 5, 6} {
+		p := curvePointAtX(t, x)
+		assert.False(p.IsInSubGroup(), "x=%d should be off-subgroup", x)
+		assert.CheckCircuit(&assertIsOnG1Circuit{}, test.WithInvalidAssignment(onG1Witness(p)),
+			test.WithCurves(ecc.BN254), test.NoProverChecks())
+	}
+
+	// A point that is not on the curve at all must also be rejected.
+	notOnCurve := curvePointAtX(t, 4)
+	notOnCurve.Y.Add(&notOnCurve.Y, &fp_bls12381.Element{1})
+	assert.False(notOnCurve.IsOnCurve())
+	assert.CheckCircuit(&assertIsOnG1Circuit{}, test.WithInvalidAssignment(onG1Witness(notOnCurve)),
+		test.WithCurves(ecc.BN254), test.NoProverChecks())
+}

std/algebra/emulated/sw_bls12381/g2.go

@@ -27,7 +27,7 @@ type G2 struct {
 
 	// Precomputed constants for the GLV+FakeGLV scalar mul ([EEMP25] §3.3).
 	g2Gen      *g2AffP // G2 generator
-	g2GenNbits *g2AffP // [2^(nbits-1)]G2 with nbits = (r.BitLen()+3)/4 + 2
+	g2GenNbits *g2AffP // [2^(nbits-1)]G2 with nbits = (r.BitLen()+3)/4 + 1
 }
 
 type g2AffP struct {
@@ -101,15 +101,15 @@ func NewG2(api frontend.API) (*G2, error) {
 			A1: *fp.NewElement("927553665492332455747201965776037880757740193453592970025027978793976877002675564980949289727957565575433344219582"),
 		},
 	}
-	// [2^(nbits-1)]G2 where nbits = (255+3)/4 + 2 = 66, so this is [2^65]G2.
+	// [2^(nbits-1)]G2 where nbits = (255+3)/4 + 1 = 65, so this is [2^64]G2.
 	g2GenNbits := &g2AffP{
 		X: fields_bls12381.E2{
-			A0: *fp.NewElement("1307001654908388153254394944417118155033503188409787277795273489312551176370209873126740711463572657296916966732684"),
-			A1: *fp.NewElement("1066804690119577865989830850277879393407029322116864061755683314318400220056817483617033672656485029228353937929571"),
+			A0: *fp.NewElement("3301848440670724259140706638333507017270073840991444062130369473322504802279997762071965204605602880645098143635804"),
+			A1: *fp.NewElement("1433247473484811575122071252549581543219505413215333519923845874993611449657376371132448509415226221815572147621561"),
 		},
 		Y: fields_bls12381.E2{
-			A0: *fp.NewElement("1233864651366532660795929818904272589705597977637697925481983092108793193162343169655985724823869788077854535468808"),
-			A1: *fp.NewElement("2703972434797875065063829955607449483769333186572810763171217085444622779819503421195150761462859837038921185079043"),
+			A0: *fp.NewElement("1426074667836811492054371797466616700133886925793677035626411405456898211637926213257131071435551354558768304539958"),
+			A1: *fp.NewElement("1551115808906539952939058035985599484510999034415063709273127009249540606507540031231198724146141840052652374305625"),
 		},
 	}
 
@@ -831,8 +831,11 @@ func (g2 *G2) scalarMulGLVAndFakeGLV(Q *G2Affine, s *Scalar, opts ...algopts.Alg
 		panic(err)
 	}
 	var st ScalarField
-	// LLL Hermite bound: u_i, v_i < γ₄·r^(1/4), fits in (BitLen+3)/4 + 2 bits.
-	nbits := (st.Modulus().BitLen()+3)/4 + 2
+	// u_i, v_i < 1.2534·r^(1/4), the LLL δ=0.99 bound on the rank-4 lattice of
+	// determinant r reduced by rationalReconstructExtG2; fits in
+	// (BitLen+3)/4 + 1 bits. Not the Hermite constant — see the derivation on
+	// [sw_emulated.Curve.scalarMulGLVAndFakeGLV].
+	nbits := (st.Modulus().BitLen()+3)/4 + 1
 
 	// handle 0-scalar and (-1)-scalar cases
 	var isScalarZero, isScalarZeroOrMinusOne, isScalarOne, isScalarMinusOne frontend.Variable

std/algebra/emulated/sw_bn254/g2.go

@@ -23,7 +23,7 @@ type G2 struct {
 
 	// Precomputed G2 generator and its multiple for GLV+FakeGLV
 	g2Gen      *g2AffP // G2 generator
-	g2GenNbits *g2AffP // [2^(nbits-1)]G2 where nbits = (r.BitLen()+3)/4 + 2
+	g2GenNbits *g2AffP // [2^(nbits-1)]G2 where nbits = (r.BitLen()+3)/4 + 1
 }
 
 type g2AffP struct {
@@ -87,16 +87,16 @@ func NewG2(api frontend.API) (*G2, error) {
 			A1: *fp.NewElement("4082367875863433681332203403145435568316851327593401208105741076214120093531"),
 		},
 	}
-	// [2^(nbits-1)]G2 where nbits = (254+3)/4 + 2 = 66, so this is [2^65]G2
+	// [2^(nbits-1)]G2 where nbits = (254+3)/4 + 1 = 65, so this is [2^64]G2
 	// The loop does nbits-1 doublings, so the generator accumulates to [2^(nbits-1)]G2
 	g2GenNbits := &g2AffP{
 		X: fields_bn254.E2{
-			A0: *fp.NewElement("6099622139700402640581725571890015148411145321742729577177999911575645303725"),
-			A1: *fp.NewElement("9870328428465937988383794519490899227160817120884239055108452134207619193487"),
+			A0: *fp.NewElement("1217614753444927156396476602882862845417913902931756228448021310943329887746"),
+			A1: *fp.NewElement("17212962051869088512269154357171035171174640071181769173167157021444635233641"),
 		},
 		Y: fields_bn254.E2{
-			A0: *fp.NewElement("16268382111792290652321980382595025991160708296314050973435867558225525677485"),
-			A1: *fp.NewElement("15377126855853471483498618408547895055706247905282062963450025729940352455943"),
+			A0: *fp.NewElement("6256619661487481776330987765107481955241329498344741334923344915358372250312"),
+			A1: *fp.NewElement("6522190482950680652192912716429811254935775044510284086077667530829579214531"),
 		},
 	}
 
@@ -580,8 +580,12 @@ func (g2 *G2) scalarMulGLVAndFakeGLV(Q *G2Affine, s *Scalar, opts ...algopts.Alg
 		panic(err)
 	}
 	var st ScalarField
-	// u1, u2, v1, v2 < c*r^{1/4} where c ≈ 1.25
-	nbits := (st.Modulus().BitLen()+3)/4 + 2
+	// u1, u2, v1, v2 < 1.2534·r^(1/4), the LLL δ=0.99 bound on the rank-4
+	// lattice of determinant r reduced by rationalReconstructExtG2; fits in
+	// (BitLen+3)/4 + 1 bits. The 1.2534 is the LLL approximation factor, not
+	// the Hermite constant — see the derivation on
+	// [sw_emulated.Curve.scalarMulGLVAndFakeGLV].
+	nbits := (st.Modulus().BitLen()+3)/4 + 1
 
 	// handle 0-scalar and (-1)-scalar cases
 	var isScalarZero, isScalarZeroOrMinusOne, isScalarOne, isScalarMinusOne frontend.Variable

std/algebra/emulated/sw_bw6761/g2.go

@@ -507,6 +507,10 @@ var bw6761G2ClassicParams = func() sw_emulated.CurveParams {
 		Eigenvalue:       lambda,
 		ThirdRootOne:     omega,
 		PreferClassicGLV: true,
+		// G2 has a nontrivial cofactor, so PrimeOrder is deliberately left
+		// false and no CofactorClearing is set: classic GLV computes its output
+		// rather than hinting it, so no binding is needed here, and
+		// AssertIsInSubgroup fails loudly instead of silently passing.
 	}
 }()
 

std/algebra/emulated/sw_emulated/fixedbase.go

@@ -48,26 +48,80 @@ import (
 // checks) plus a y materialization (one check). The y-coordinate is
 // materialized once, before the final complete addition.
 
-// combDefaultWindow is the default window width of the fixed-base comb. With
-// 64-bit limb emulation it is supported by all built-in curves (the recoded
-// scalar fits the scalar-field limb capacity) and is close to the
-// constraint-count optimum in R1CS.
-const combDefaultWindow = 8
-
 // combPlonkWindow is the window width used on PLONKish backends. The one-hot
 // selection's wide constant linear combinations are free in R1CS but expand
 // into one addition gate per term in PLONK, making the selector cost scale
 // with 2^w·nbLimbs per window; a smaller window rebalances selector versus
-// chain-addition cost. With a partial top window, w=5 is the current SCS
-// optimum measured on secp256k1.
-const combPlonkWindow = 5
+// chain-addition cost, which is why this sits below the R1CS widths that
+// combOptimalWindow picks.
+//
+// w=6 is the measured SCS optimum, and it is the optimum for every supported
+// curve rather than just one, so a single constant suffices here (SCS counts
+// for w=4..8, see TestScalarMulBaseCombPlonkWindow):
+//
+//	secp256k1   101776  92472  89890  97852  123576
+//	BN254 G1    101538  91162  89652  97614  120402
+//	BLS12-381   148187 132226 129875 141905 176699
+const combPlonkWindow = 6
+
+// combOptimalWindow returns the R1CS-optimal comb window width for the given
+// base-field limb count and scalar bit length. It minimises an estimated
+// constraint count:
+//
+//	cost(w) = (nw−1)·selectCost(w) + selectCost(tw) + (nw−1)·addCost
+//
+// where nw = ⌈n/w⌉, tw = n − w·(nw−1) is the top-window width, and
+// addCost ≈ 100·nbLimbs approximates the R1CS constraints per chain
+// addition (emulated complete addition scales roughly linearly with nbLimbs).
+// selectCost is the minimum-rb combSelect multiplication count at width w
+// (see combSelect for the bilinear one-hot formula).
+//
+// Validated against measured constraint counts:
+//
+//	secp256k1, BN254 G1 (nbLimbs=4, n≈256): optimum w=9 (−2.3 % vs w=8)
+//	BLS12-381 G1        (nbLimbs=6, n≈255): optimum w=10 (−3.2 % vs w=8)
+func combOptimalWindow(nbLimbs, n int) int {
+	oneHotCost := func(k int) int {
+		if k <= 1 {
+			return 0
+		}
+		return (1 << k) - 2
+	}
+	nbOut := 2 * nbLimbs
+	selectCost := func(w int) int {
+		if w <= 0 {
+			return 0
+		}
+		best := oneHotCost(w) // rb=0: full one-hot on all w bits
+		for rb := 1; rb <= w; rb++ {
+			if c := oneHotCost(rb) + oneHotCost(w-rb) + (1<<rb)*nbOut; c < best {
+				best = c
+			}
+		}
+		return best
+	}
+	addCost := 100 * nbLimbs
+	best, bestCost := 8, -1
+	for w := 4; w <= 16; w++ {
+		nw := (n + w - 1) / w
+		tw := n - w*(nw-1)
+		total := (nw-1)*selectCost(w) + selectCost(tw) + (nw-1)*addCost
+		if bestCost < 0 || total < bestCost {
+			bestCost = total
+			best = w
+		}
+	}
+	return best
+}
 
 // combWindow returns the comb window width for the current backend.
 func (c *Curve[B, S]) combWindow() int {
 	if _, ok := c.api.Compiler().(frontend.PlonkAPI); ok {
 		return combPlonkWindow
 	}
-	return combDefaultWindow
+	var fp B
+	var fs S
+	return combOptimalWindow(int(fp.NbLimbs()), fs.Modulus().BitLen())
 }
 
 // combData holds the compile-time data of the comb: the constant window