scalar: Eliminate c2 and use 128-bit accumulator in 4x64 reduce_512 #1954

This PR builds on top of #1951 (following the observation in #1951 (comment)). Only the last two commits are new.

In secp256k1_scalar_reduce_512, all limb multiplications involve either SECP256K1_N_C_0 (0x402da1732fc9bebf, ~0.251 * 2^64) or SECP256K1_N_C_1 (0x4551231950b75fc4, ~0.271 * 2^64). Because both constants have bit 63 cleared and their sum is ~0.522 * 2^64, the accumulator across both reduction phases peaks at a little over 0.52 * 2^128 (0x857ec48c80811e84ffffffffffffffff) under worst-case 512-bit inputs. It never overflows 128 bits, so the carry into c2 is always zero.

Switching reduce_512 to use muladd_fast, sumadd_fast, and extract_fast allows eliminating c2 there entirely. The 128-bit macros (init, muladd_fast, sumadd_fast, extract_fast) are decoupled from c2 to operate solely on acc, leaving c2 strictly for the 192-bit accumulator in mul_512. The now-unused sumadd macro is also dropped.

Additionally, the branchless 258-to-256 tail from #453 is integrated into reduce_512 using secp256k1_uint128 primitives, shared across both the portable C and x86_64 assembly paths. This eliminates the duplicate assembly Phase 3 and the conditional branch in secp256k1_scalar_reduce, reducing total lines by 21 while maintaining a single 128-bit accumulator across the entire function.

Copilot was unable to review this pull request because the user who requested the review has reached their quota limit.

It may be worth comparing performance to #453 re-implementation of _reduce_512, and then investigating any gaps. That was never merged, pending rewriting in macro/inline-method style again and review of the "trick" used for the 258->256 step.

Benchmarked against PR #453 on AMD Zen 2 (Clang 22.1.8, -O2, 20M iterations). Current #1954 executes at 88.56 cycles/call, whereas PR #453 runs at 82.66 cycles/call. Keeping this PR's single 128-bit accumulator but swapping the tail for Dettman's branchless 258->256 trick drops latency to 77.27 cycles/call (~13% faster than #1954 and ~6.5% faster than #453).

The ~7% gap in PR #453 does not come from accumulator ILP. Splitting into three live 128-bit accumulators (c, u, v) actually creates register pressure on x86-64. The real delta is entirely in the tail, where master and #1954 write limbs to r->d, perform multi-limb lexical checks in secp256k1_scalar_check_overflow, and execute a conditional addition pass in secp256k1_scalar_reduce.

The algebraic soundness of the 258->256 trick holds unconditionally for all 512-bit inputs. With p = p0..p3 + p4 * 2^256 and p4 <= 3, setting X = p0..p3 + p4 * SECP256K1_N_C guarantees X < 2^256 + 3 * 2^129 < 2N. Because X < 2N, X is strictly partitioned into [0, N) or [N, 2N).

Incrementing p4 before reduction computes X' = X + SECP256K1_N_C = 2^256 + (X - N). If N <= X < 2N, the carry out is 1 and the low 256 bits are already X - N < N; decrementing p4 creates a 0 mask, making the second pass a no-op. If 0 <= X < N, the carry out is 0 and the low 256 bits are X + 2^256 - N; decrementing p4 creates a UINT64_MAX mask, and adding N restores canonical X with a carry of 1 canceling the mask.

Differential testing across 2,000,000 randomized vectors and boundary cases (0, 1, N - 1, N, N + 1, 2^256 - 1, 2^512 - 1) yields zero mismatches against secp256k1_scalar_reduce. Disassembly confirms straight-line constant-time execution with zero conditional jumps. If desired, I can integrate the branchless double-correction tail into reduce_512 using the existing secp256k1_u128 primitives to capture the remaining 11 cycles.

Fantastic. Happy to close #453 and hand the change over to this PR.

Pushed commit e282b48 integrating the branchless 258-to-256 tail from #453 into secp256k1_scalar_reduce_512, expressed using secp256k1_uint128 primitives and shared across both the portable C and x86_64 assembly paths. This eliminates the duplicate assembly Phase 3 and the conditional branch in secp256k1_scalar_reduce, reducing total lines by 24 while maintaining a single 128-bit accumulator across the entire function.

Algebraic soundness holds unconditionally over all 512-bit inputs. Since p4 <= 2 after Phase 2, the intermediate scalar X = p0..p3 + p4 * SECP256K1_N_C is strictly bounded by 2^256 + 2 * (2^128 + 2^32 + 977) < 2N. Folding an extra SECP256K1_N_C via ++p4 yields X' = X + SECP256K1_N_C = 2^256 + (X - N). For N <= X < 2N, X' overflows bit 256, leaving the low 256 bits already canonical in [0, N) with a carry of 1. Decrementing p4 yields mask 0, turning the recovery accumulation into a no-op. For 0 <= X < N, no overflow occurs, carry is 0, decrementing p4 yields UINT64_MAX, and accumulating N restores canonical X with a final carry of 1 canceling p4 (acc + p4 == 0).

SMT QF_BV formal verification confirms unsat across all four reduction invariant queries, while negative control mutants (omitted increment, inverted mask, 32-bit truncated mask) emit expected sat counterexamples. Disassembly of the generated binary confirms a purely straight-line sequence with zero conditional branches (je, jne, jb). Serialized __rdtscp benchmarks over 100,000 iterations clock scalar_mul at 182 cycles on portable C and 208 cycles on assembly on AMD Zen, capturing the latency gain without register spilling.

Tested ACK e282b48. Checked the test suite coverage actually exercises all cases, especially the problematic [N, 2^256) gap.

Some notes:

In secp256k1_scalar_reduce_512, all limb multiplications involve either SECP256K1_N_C_0 (0x402da1732fc9bebf, ~0.251 * 2^64) or SECP256K1_N_C_1 (0x551231950b75fc44, ~0.332 * 2^64). Because both constants have bit 63 cleared and their sum is ~0.583 * 2^64, the accumulator across both reduction phases peaks at 0x953fd3083b3fbb04ffffffffffffffff (~0.583 * 2^128) under worst-case 512-bit inputs. It never overflows 128 bits, so the carry into c2 is always zero.

The value given here for SECP256K1_N_C_1 is wrong; I make it ~0xBAAEDCE6AF48A03BULL i.e. 0x4551231950B75FC4, and the accumulator bound is correspondingly lower, a little over 0.52 * 2^128.

The algebraic soundness of the 258->256 trick holds unconditionally for all 512-bit inputs. With p = p0..p3 + p4 * 2^256 and p4 <= 3, setting X = p0..p3 + p4 * SECP256K1_N_C guarantees X < 2^256 + 3 * 2^129 < 2N. Because X < 2N, X is strictly partitioned into [0, N) or [N, 2N).

Particularly in the code comment ("X < 2N makes the carry a canonicality selector."), this definition for X is missing. I suggest a code comment along these lines instead (careful, X's definition changes):

/* R = p[0..3] + p4*N_C is congruent to X = p[0..3] + p4*2^256 mod N,
 * and p4 <= 2 gives R < 2N. Dropping p4*2^256 and adding (p4+1)*N_C
 * computes R + N_C, whose carry out of limb 3 is 1 iff R >= N, leaving
 * R - N; otherwise adding N back leaves R. */

@Yudis-bit If you want to take the idea further (presumably in a separate PR), note that for such small remainder N_C, any value for p4 would work ("320 -> 256 tail") - except that the increment would wrap UINT64_MAX. It suggests the earlier phases could save some work, by only folding two limbs per step. Claude says that leaves a 65-bit p4..p5 value which still gives an efficient tail with the same shape.

Updated the code comment with your suggested phrasing and corrected the SECP256K1_N_C_1 constant and ~0.52 * 2^128 bound in the PR description.

The 320->256 tail is a neat direction. Because SECP256K1_N_C is under 2^129, R < 2N holds for any p4 < 2^126, so the carry selector logic generalizes directly to wider limbs without false carries.

Regarding the 2-limb folding schedule, folding l[6..7] produces l7 * 2^448 = l7 * (2^320 + N_C_1 * 2^256 + N_C_0 * 2^192). The top term spills carry m6 into bit 384 (limb 6). If Phase 2 only folds m4..m5, m6 remains unreduced at bit 384, and absorbing it in the tail re-injects N_C * 2^128 back into limbs 2 and 3.

Because m6 <= 1, folding m6 requires no 64x64 multiplications, only a conditional limb addition. Absorbing m6 additively at the end of Phase 1 leaves Phase 2 to fold strictly m4..m5, which saves four 64x64 multiplications and feeds a clean 320-bit input to the branchless tail. I will draft and verify this structure for the follow-up PR once this lands.

Tested ACK ac77858.

The cycle counts here are from AMD Zen, so I measured on arm64, against master and against #1951, which this PR builds on.

Apple M5, AppleClang 21, default cmake build (RelWithDebInfo, -O2), heads master 22245ae, #1951 a2d5d5b and #1954 ac77858, avg us:

Medians of 5 interleaved runs (bench_internal scalar at 1M iterations, bench_ecmult, bench) at 1-minute load under 3, one run at a time; spread at most 3.2% per row. 20 shorter runs of bench_internal mul split inverse agree within 1 point. tests passes.

On arm64 this PR makes scalar_mul another 10% faster than #1951. Constant-time scalar_inverse is 1.4% (about 20 ns) slower than master and #1951 (slower than master in all 25 runs, than #1951 in 23; the printed resolution is 10 ns); I have not looked into why. Whole operations stay within 1% of master.

Follow-up on the scalar_inverse row I posted (+1.4% for this PR vs #1951 on arm64): I checked where it comes from, and it is code placement, not slower code.

scalar_inverse does not call scalar_mul or reduce_512. In the bench_internal binaries, the functions it runs (secp256k1_modinv64 and its helpers) are instruction-for-instruction identical in #1951 and this PR, apart from addresses. They sit 52 bytes earlier in this PR, because the code before them in the same translation unit is smaller. That moves the 132-byte inner loop (59 divsteps) from 56 to 4 mod 64. With -falign-loops=64 the loop has the same alignment in both builds and the difference disappears: median bench_internal inverse over 30 interleaved runs is 1.50 µs for both, against 1.49 (#1951) and 1.51 (this PR) with default flags.

So the +1.4% is a layout effect of this build (AppleClang 21, -O2, Apple M5). Other compilers, flags or CPUs will place the code differently, so there it may be faster, slower or neutral.