Skip to content

MSM endomorphisms (GLV/GLS) + CIOS field arithmetic (bn254/bls12-381)#75

Open
OBrezhniev wants to merge 17 commits into
masterfrom
feature/msm-signed-buckets
Open

MSM endomorphisms (GLV/GLS) + CIOS field arithmetic (bn254/bls12-381)#75
OBrezhniev wants to merge 17 commits into
masterfrom
feature/msm-signed-buckets

Conversation

@OBrezhniev

Copy link
Copy Markdown
Member

MSM endomorphisms + CIOS field arithmetic (bn254/bls12-381)

Summary

Rewrites the hot arithmetic kernels and adds an endomorphism-accelerated MSM
module. Every change is bit-identical to the previous implementation and
validated against BigInt references, adversarial inputs, and the full test
suite. Measured on a 20-core x86 box (interleaved A/B):

kernel before after
Montgomery mul (bn254) 71.3 ns 52.0 ns (−27%)
Montgomery square 78.8 ns ~55 ns (−30%)
G1 MSM 4k–32k pts −35–39% (GLV + batch-affine)
G2 MSM 4k–6k pts −45% vs plain (GLS)
bls12-381 G1 MSM 4k pts −35–43% (GLV)

Changes

Field arithmetic (build_f1m.js)

  • _mul: product-scanning → CIOS (operand scanning, all in locals, q-limbs as
    immediates). Generic over limb count; frm gets it too via buildF1.
  • no-carry variant when q < R/2 (spare bit in the top limb): the t < 2q
    invariant makes the overflow limb provably dead, so pass-boundary folds
    collapse. Guarded: moduli with q ≥ R/2 keep full carry tracking.
  • _square: dedicated CIOS square (j ≥ i products; doubled cross products via
    a lo/hi split that keeps every step within u64 for any limb count).

Signed-digit Pippenger (build_multiexp.js)

  • Window digits recoded to [−2^(c−1), 2^(c−1)−1] with a carry-guard window;
    bucket count halved (negation is free on the curve).

Batch-affine MSM module (src/as/msm_batch.ts, AssemblyScript → build/msm_batch.wasm)

  • Buckets in affine coordinates, batched affine additions with one field
    inversion per batch (Montgomery's trick). Conflict-free scheduling via two
    counting sorts (bucket, then rank-within-bucket = layer); layers keep
    ascending point order so bases are read near-sequentially.
  • GLV (bn254 + bls12-381 G1): k = k1 + λk2, sub-scalars ≤ 2^127, φ-bases
    materialized once per call; signs ride the signed-digit machinery. Shared
    decomposition code with per-curve baked constants (bls uses the closed-form
    basis (λ,−1),(1,z²)).
  • GLS (bn254 G2): 4-dim decomposition (LLL-derived basis, |k_i| ≤ 2^64),
    ψ-orbit tables; gated on the materialized working set (≤ 3 MiB) with
    internal fallback.
  • The module has no data segments (it shares the caller's linear memory),
    imports all field/group ops, and serves both curves and both groups from
    one binary.

Validation

  • Suite: 114 passing. Decompositions verified against BigInt mirrors
    (values, signs, congruence; 100k+ scalars offline, 2–3k in-wasm) including
    0, 1, r−1, λ. λ/β/ψ variants verified empirically against λ·P on the curve.
  • MSM gauntlets vs the previous implementation: identical points, P/−P
    cancellation, zero scalars, r−1, zero points, deliberately misaligned
    buffers — bit-exact.
  • Field kernels validated on random pairs + edge cases (0, 1, (q−1)²).

Notes for reviewers

  • ffjavascript consumes the prebuilts via its gen-wasm re-vendoring; the
    companion PR there must land together with this one.
  • The square deliberately does NOT use the no-carry simplification: its
    doubled cross products push intermediate t toward 3q, so the overflow limb
    is genuinely live there.

🤖 Generated with Claude Code

…eed to verify if it's prime or not during curve initialization.
Recode each scalar into signed c-bit windows in [-2^(c-1), 2^(c-1)-1] with carry
propagation, so only 2^(c-1) buckets are needed per window instead of 2^c-1.
Negative digits are handled by subtracting the base (point negation is free).

A small carry buffer (1 byte per window/point) is filled in a recode pass; the
signed digit is reconstructed in the chunk loop from the raw window plus the
carry, then reduced with a running-sum pass. An extra guard window absorbs the
final carry so arbitrary scalarSize-byte scalars work, and the window is clamped
to c>=2 (the signed range is degenerate at c=1).

Bit-exact with the previous unsigned multiexp (G1/G2, affine + jacobian).
Regenerated the prebuilt bn128/bls12381 wasm.
A standalone wasm module implementing the signed-digit Pippenger MSM with
batch-affine bucket accumulation: buckets are kept in affine coordinates and
filled with batched affine additions -- one field inversion per batch via
Montgomery's trick -- so a bucket add costs ~2M+1S plus ~3M of amortized
inversion instead of an 11M jacobian addMixed.

The module contains only control flow: all field/group arithmetic is imported
from the main curve module and both share one linear memory, so the same
binary serves bn128 and bls12381, G1 and G2 (the base-field element size is a
runtime parameter). It has no data segments (nothing to clobber in the shared
memory) and grows the memory when its bump allocations pass the current
capacity.

Scheduling is conflict-free by construction: per window, points are counted
per bucket and stably grouped by their rank within the bucket (two counting
sorts); rank r = layer r. A layer touches each bucket at most once, so a
batch never contains two adds into the same bucket, and layers keep ascending
point order so the bases are read near-sequentially (random access stays in
the small bucket array -- this matters when many workers contend for cache).
Doubling (P == bucket) and cancellation (P == -bucket) are detected at
schedule time; the batch is flushed at layer boundaries.

Built with assemblyscript (devDependency): npm run build_msm_batch ->
build/msm_batch.wasm (~3.6KB).
Replace the product-scanning Montgomery multiplication (dual c0/c1 column
accumulators; ~6 dependent ops per partial product) with coarsely integrated
operand scanning: per y-limb, one multiply-accumulate pass over the x-limbs
and one reduction pass, all in i64 locals with a single running carry, and
the q-limbs emitted as immediates instead of loads. Same 32-bit limbs, same
R = 2^(32*n32), bit-identical results, generic over the limb count (bn128
f1m/frm and bls12381 alike -- frm gets it too via buildF1, so FFTs benefit).

Bound check: each step t_j + x_j*y_i + c <= 2^64-1 exactly; result < 2q, so
the single conditional subtraction is unchanged from the previous code.

Measured 71.3 -> 54.7 ns per mul (with the ffjavascript -O2 vendoring; see
that repo). Validated against a BigInt reference on random pairs and the
full wasmcurves/ffjavascript/snarkjs suites; groth16 proves verify.
The dedicated product-scanning square (doubled cross products, dual-carry
column accumulators) computes 22% fewer partial products than a full
multiplication but measured ~28% SLOWER than mul(x,x) on the CIOS mul
(78.8 vs 61.6 ns) -- the column-accumulator dependency chains cost more than
the saved multiplies. A dedicated CIOS square with the doubling trick would
recover at most ~10% over mul(x,x), which is not worth the extra codegen, so
_square now just calls _mul(x,x,r). Bit-identical results; suites pass.
Replace the mul(x,x) delegation with a true CIOS square: pass i computes only
x_j*x_i for j >= i -- the diagonal once, cross products doubled via a lo/hi
split (s += (p & 2^32-1) << 1; carry += (p >> 32) << 1), which keeps every
step inside u64 (running carry < 2^33+8 for any limb count). The reduction
phase is identical to the CIOS _mul. Saves n32*(n32-1)/2 multiplies vs
mul(x,x).

Measured (bn254): 49.6 vs 57.8 ns in the prototype; vendored 57.2 vs 62.4 ns
(delegation) and 78.8 ns (original product-scanning square) -- ~27% total.
Validated on 500 pairs incl 0, 1, q-1 against a BigInt reference and the
previous implementation; all suites pass; groth16 proves verify.
gnark-style no-carry optimization, adapted to 32-bit limbs: when the modulus
has a spare bit in its top limb (q < 2^(32*n32-1)), the CIOS invariant t < 2q
bounds every intermediate T = t + x*y_i + m*q < q*2^33 < 2^32 * R, so the
accumulator provably never outgrows n32+1 limbs. The overflow limb and the
pass-boundary mask/fold work then collapse to plain assignments, and the
final "top limb set" branch disappears (that limb is always zero).

The codegen guards on the condition: moduli with q >= R/2 (e.g. a 2^255-19
style prime packed into 8 limbs) keep the previous full carry tracking.
All four production fields (bn254 q/r, bls12-381 q/r) qualify.

The square keeps its full tracking on purpose: its doubled cross products
push intermediate t toward 3q, so its overflow limb is genuinely live.

Measured: 54.7 -> 52.0 ns per mul. Validated on 206 edge cases (0, 1,
(q-1)^2, ...) and random pairs against a BigInt reference; all suites pass;
groth16 proves verify (authV3 672 ms, sha256 6.92 s medians).
k*P = k1*P + k2*phi(P) with phi(x,y) = (beta*x, y) = lambda*(x,y) and
k = k1 + lambda*k2 (mod r), |k1|,|k2| < 2^127 -- the MSM runs over 2n points
with 128-bit sub-scalars, halving the window count. Decomposition is done
in-wasm (division-free rounding against precomputed floor(|b_i|*2^256/r)
constants; exact congruence regardless of rounding, verified against a
BigInt mirror on 2000+ cases including 0, 1, r-1, lambda). Sub-scalar signs
ride the existing signed-digit machinery (per-item digit flip); the phi
bases are materialized once per call (one field mul per point, beta stored
in Montgomery form). The MSM core is refactored into a shared msmRun() and
perm entries now carry the resolved base pointer.

Fixes a latent alignment bug found during validation: perm entries encoded
base pointers at 8-byte granularity, but the bump allocators only guarantee
4-byte alignment -- a bases buffer at 4 mod 8 silently read every point 4
bytes off. Entries now encode at 4-byte granularity (regression-tested with
a deliberately misaligned allocation stream).

Constants are bn254-specific (hardcoded, no data segments); the export falls
back to the generic path for unexpected scalar/field sizes, and bls12-381
simply does not advertise GLV. Measured single-instance: 35-39% faster than
the batch path at 4k-32k points. Full gauntlet (identical points, P/-P
cancellation, zero scalars, r-1, zero point) bit-exact vs the main module.
k*P = k1*P + k2*psi(P) + k3*psi^2(P) + k4*psi^3(P), where psi is the twist
Frobenius endomorphism acting on G2 as multiplication by lambda = 6x^2 (x =
BN parameter), and k = sum k_i*lambda^(i-1) (mod r) with |k_i| <= 2^64 --
64-bit sub-scalars cut the window count to a quarter.

The 4-dim lattice basis comes from LLL over {z: sum z_i*lambda^i = 0 mod r}
(entries are all small multiples of x, verified congruences); Babai rounding
uses division-free mulhi constants like the G1 GLV path, with the identity
exact regardless of rounding error (bounds verified on 100k+ scalars in the
derivation and 3k in-wasm cases against a BigInt mirror). In coordinates,
psi(x,y) = (conj(x)*Gx, conj(y)*Gy) with Gx = xi^((q-1)/3), Gy = xi^((q-1)/2)
(xi = 9+u; variant verified empirically against lambda*P on the curve);
psi^2 multiplies by the Fq norms and psi^3 by their products, so the three
extra point tables cost ~6 Fq2 muls per point. Requires a new f_conj import
(f2m_conjugate; the G1 instance wires a harmless copy).

The psi-orbit materialization quadruples the bases working set, so the entry
gates on 4n*sG <= 3 MiB (measured: +24% vs batch inside the boundary, -13%
beyond it) and falls back to the plain batch path above it; non-bn254-G2
calls fall through unconditionally.

Measured single-instance G2 MSM at n=4096-6144: 24% faster than batch, 45%
faster than plain. Full gauntlet vs g2m_multiexpAffine bit-exact.
…kets

# Conflicts:
#	build/bn128_wasm_gzip.js
The GLV entry now dispatches on the field size: bn254 (n8f=32) and
bls12-381 (n8f=48) each get a baked constants block in a shared layout
(W1/W2 rounding constants as 7 limbs, the four basis entries as 4 limbs,
beta in the curve's Montgomery form at a fixed offset), and the
decomposition code is shared -- both curves have b1 < 0, b2 > 0, so the
sign structure is identical. bls12-381 uses the closed-form basis
(lambda, -1), (1, z^2) with det exactly r (the generic Euclid degenerates
there because lambda ~ sqrt(r)); |k_i| <= 2^127, verified against a BigInt
mirror on 2000+ cases including r-1, and the beta/lambda pairing verified
empirically against lambda*P. G2 GLS remains bn254-only (bls G2 sizes fall
through to the generic batch path).

Also fixes a scratch-buffer overflow introduced by the generalization: the
W products are 15 limbs (60 bytes) and the pt scratch was still 52 bytes,
corrupting c1 for large scalars (caught by the bn254 regression gauntlet).

Measured (bls12-381 G1, n=4096): 121 ms vs 211 (batch) / 186 (plain) --
~40% faster.
Migrates .eslintrc.js/.eslintignore to eslint.config.mjs. Lint + 114 tests pass.
@OBrezhniev

Copy link
Copy Markdown
Member Author

js-yaml/picomatch bumped in-range; serialize-javascript ^7.0.5 and
diff ^8.0.3 overridden (mocha pins vulnerable ranges; the patched majors
are outside them). npm audit clean; 114 tests pass.
Ported ffiasm's getCriticalNumbers generator (test/fieldasm.js) to hit
32/64-bit limb boundaries and q/2, targeting the CIOS Montgomery
noCarry fast path introduced for the mul/square rewrite.
msm_batch.wasm (signed-digit Pippenger + GLV/GLS endomorphism) was
never instantiated by anything in this repo's own suite -- only
ffjavascript links it at runtime. Reproduces that linking directly
and checks multiexpAffine/multiexpAffineGLV/multiexpAffineGLS against
the trusted g1m/g2m_multiexpAffine at N=0, N=1, all-zero scalars,
duplicate points, an infinity base point, and GLV edge scalars.

Also covers f1m_batchInverse's isZero-skip path (a zero in the batch
must not corrupt other elements' inverses), relevant to batch-affine
addition where equal/inverse-point collisions produce a zero
z-delta.
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

1 participant