From 11d8533830e0e0c7cd1f4ecf8b15c6a4ca2c4a1d Mon Sep 17 00:00:00 2001 From: Dave Collins Date: Thu, 9 Jul 2026 01:28:21 -0500 Subject: [PATCH 01/17] secp256k1: Add caveats to documentation. This updates the README.md and doc.go files to include some important caveats about the constant time and memory clearing behavior. Specifically, while the behavior has been manually verified by carefully reviewing the generated assembly as of the most recent version of Go, it is impossible to guarantee the compiler will not change in the future. --- dcrec/secp256k1/README.md | 19 +++++++++++++++++-- dcrec/secp256k1/doc.go | 15 +++++++++++++++ 2 files changed, 32 insertions(+), 2 deletions(-) diff --git a/dcrec/secp256k1/README.md b/dcrec/secp256k1/README.md index 53f3839117..4fcfa18616 100644 --- a/dcrec/secp256k1/README.md +++ b/dcrec/secp256k1/README.md @@ -16,8 +16,8 @@ structures and functions for working with public and private secp256k1 keys. In addition, subpackages are provided to produce, verify, parse, and serialize ECDSA signatures and EC-Schnorr-DCRv0 (a custom Schnorr-based signature scheme -specific to Decred) signatures. See the README.md files in the relevant sub -packages for more details about those aspects. +specific to Decred) signatures. See the README.md files in the relevant +subpackages for more details about those aspects. ## Features @@ -95,6 +95,21 @@ secp256k1 elliptic curve cryptography. The secp256k1 domain parameters are specified in [SEC 2: Recommended Elliptic Curve Domain Parameters](https://www.secg.org/sec2-v2.pdf). +## Caveats + +Constant time implementations are designed to eliminate data dependent timing +variation at the algorithmic level. This property depends on the Go compiler +continuing to compile the relevant operations into constant-time machine +instructions. + +This package also consistently uses temporary stack allocated buffers which it +attempts to clear in order to reduce the likelihood of sensitive data lingering +in memory, but Go does not guarantee such writes are never optimized away. + +Significant effort has been made to ensure correctness including careful review +of the generated assembly, but it is impossible to guarantee exact behavior of +future versions of the Go compiler. + ## secp256k1 use in Decred At the time of this writing, the primary public key cryptography in widespread diff --git a/dcrec/secp256k1/doc.go b/dcrec/secp256k1/doc.go index d946650d73..69121d1164 100644 --- a/dcrec/secp256k1/doc.go +++ b/dcrec/secp256k1/doc.go @@ -50,6 +50,21 @@ use optimized secp256k1 elliptic curve cryptography. Finally, a comprehensive suite of tests is provided to provide a high level of quality assurance. +## Caveats + +Constant time implementations are designed to eliminate data dependent timing +variation at the algorithmic level. This property depends on the Go compiler +continuing to compile the relevant operations into constant time machine +instructions. + +This package also consistently uses temporary stack allocated buffers which it +attempts to clear in order to reduce the likelihood of sensitive data lingering +in memory, but Go does not guarantee such writes are never optimized away. + +Significant effort has been made to ensure correctness, including careful review +of the generated assembly, but it is impossible to guarantee exact behavior of +future versions of the Go compiler. + # Use of secp256k1 in Decred At the time of this writing, the primary public key cryptography in widespread From edc5a2ce64c96dc3186ab0a9e6a592be2160c33c Mon Sep 17 00:00:00 2001 From: Dave Collins Date: Thu, 9 Jul 2026 03:35:47 -0500 Subject: [PATCH 02/17] secp256k1: Add scalar negation algebraic tests. This adds a new test helper for some algebraic properties of scalar negation and calls it from the test that handles specific edge cases as well as the randomized test. --- dcrec/secp256k1/modnscalar_test.go | 26 ++++++++++++++++++++++++++ 1 file changed, 26 insertions(+) diff --git a/dcrec/secp256k1/modnscalar_test.go b/dcrec/secp256k1/modnscalar_test.go index 251708a6ff..b894e9b4ce 100644 --- a/dcrec/secp256k1/modnscalar_test.go +++ b/dcrec/secp256k1/modnscalar_test.go @@ -979,6 +979,26 @@ func TestModNScalarSquareRandom(t *testing.T) { } } +// checkModNScalarNegateProps checks algebraic properties of scalar negation. +func checkModNScalarNegateProps(t *testing.T, val *ModNScalar) { + t.Helper() + + negation := new(ModNScalar).NegateVal(val) + + // Ensure involution produces the original value. That is -(-x) == x. + involution := new(ModNScalar).Set(negation).Negate() + if !involution.Equals(val) { + t.Errorf("mismatched involution -- got: %v, want: %v", involution, + val) + } + + // Ensure x + (-x) == 0. + sum := new(ModNScalar).Add2(val, negation) + if !sum.IsZero() { + t.Errorf("x + (-x) != 0 -- got: %v, want: 0", sum) + } +} + // TestModNScalarNegate ensures that negating scalars works as expected for edge // cases. func TestModNScalarNegate(t *testing.T) { @@ -1051,6 +1071,9 @@ func TestModNScalarNegate(t *testing.T) { result2, expected) continue } + + // Ensure algebraic properties with negation hold. + checkModNScalarNegateProps(t, s) } } @@ -1086,6 +1109,9 @@ func TestModNScalarNegateRandom(t *testing.T) { "big int result: %x\nscalar result %v", bigIntVal, modNVal, bigIntResult, modNValResult) } + + // Ensure algebraic properties with negation hold. + assertModNScalarNegateProperties(t, modNVal) } } From cb8f8c9d2719fc4c1906a37c4d0320e4956019c3 Mon Sep 17 00:00:00 2001 From: Dave Collins Date: Thu, 9 Jul 2026 04:10:05 -0500 Subject: [PATCH 03/17] secp256k1: Add scalar inverse algebraic tests. This adds a new test helper for some algebraic properties of modular multiplicative inverses of scalars and calls it from the test that handles specific edge cases as well as the randomized test. --- dcrec/secp256k1/modnscalar_test.go | 34 +++++++++++++++++++++++++++++- 1 file changed, 33 insertions(+), 1 deletion(-) diff --git a/dcrec/secp256k1/modnscalar_test.go b/dcrec/secp256k1/modnscalar_test.go index b894e9b4ce..e4a9307b63 100644 --- a/dcrec/secp256k1/modnscalar_test.go +++ b/dcrec/secp256k1/modnscalar_test.go @@ -1111,7 +1111,33 @@ func TestModNScalarNegateRandom(t *testing.T) { } // Ensure algebraic properties with negation hold. - assertModNScalarNegateProperties(t, modNVal) + checkModNScalarNegateProps(t, modNVal) + } +} + +// checkModNScalarInverseProps checks algebraic properties of the +// multiplicative modular inverse of scalars. +func checkModNScalarInverseProps(t *testing.T, val *ModNScalar) { + t.Helper() + + // Exception for 0 which does not have a modular multiplicative inverse. + if val.IsZero() { + return + } + + inverse := new(ModNScalar).InverseValNonConst(val) + + // Ensure x * x^-1 ≡ 1 (mod N). + product := new(ModNScalar).Mul2(val, inverse) + one := new(ModNScalar).SetInt(1) + if !product.Equals(one) { + t.Errorf("x * x^-1 != 1 (mod N) -- got: %v, want: %v", product, one) + } + + // Ensure (x^-1)^-1 == x. + doubleInv := new(ModNScalar).InverseValNonConst(inverse) + if !doubleInv.Equals(val) { + t.Errorf("(x^-1)^-1 != x (mod N) -- got: %v, want: %v", doubleInv, val) } } @@ -1189,6 +1215,9 @@ func TestModNScalarInverseNonConst(t *testing.T) { result2, expected) continue } + + // Ensure algebraic properties with inverses hold. + checkModNScalarInverseProps(t, s) } } @@ -1224,6 +1253,9 @@ func TestModNScalarInverseNonConstRandom(t *testing.T) { "big int result: %x\nscalar result %v", bigIntVal, modNVal, bigIntResult, modNValResult) } + + // Ensure algebraic properties with inverses hold. + checkModNScalarInverseProps(t, modNVal) } } From a073c263369e60dd93c09e940061a4dc99e0eac8 Mon Sep 17 00:00:00 2001 From: Dave Collins Date: Thu, 9 Jul 2026 04:30:46 -0500 Subject: [PATCH 04/17] secp256k1: Add scalar half order algebraic tests. This adds a new test helper for some algebraic properties of the scalar half order check and calls it from the test that handles specific edge cases as well as the randomized test. --- dcrec/secp256k1/modnscalar_test.go | 30 +++++++++++++++++++++++++++++- 1 file changed, 29 insertions(+), 1 deletion(-) diff --git a/dcrec/secp256k1/modnscalar_test.go b/dcrec/secp256k1/modnscalar_test.go index e4a9307b63..448378e0cf 100644 --- a/dcrec/secp256k1/modnscalar_test.go +++ b/dcrec/secp256k1/modnscalar_test.go @@ -1259,6 +1259,27 @@ func TestModNScalarInverseNonConstRandom(t *testing.T) { } } +// checkModNScalarHalfOrderProps checks algebraic properties of the half order +// comparison. +func checkModNScalarHalfOrderProps(t *testing.T, val *ModNScalar) { + t.Helper() + + // Exception for 0 which is the only scalar where negation does not produce + // a distinct value. + if val.IsZero() { + return + } + + // Ensure negating a scalar flips whether it is over half order. + // + // x > N/2 implies N-x < N/2 and vice versa. + negation := new(ModNScalar).NegateVal(val) + if val.IsOverHalfOrder() == negation.IsOverHalfOrder() { + t.Errorf("negation did not flip half order -- val: %v, negation: %v", + val, negation) + } +} + // TestModNScalarIsOverHalfOrder ensures that scalars report whether or not they // exceeed the half order works as expected for edge cases. func TestModNScalarIsOverHalfOrder(t *testing.T) { @@ -1305,12 +1326,16 @@ func TestModNScalarIsOverHalfOrder(t *testing.T) { }} for _, test := range tests { - result := mustModNScalarWithOverflow(test.in).IsOverHalfOrder() + s := mustModNScalarWithOverflow(test.in) + result := s.IsOverHalfOrder() if result != test.expected { t.Errorf("%s: unexpected result -- got: %v, want: %v", test.name, result, test.expected) continue } + + // Ensure algebraic properties with the half order hold. + checkModNScalarHalfOrderProps(t, s) } } @@ -1344,5 +1369,8 @@ func TestModNScalarIsOverHalfOrderRandom(t *testing.T) { "in: %v\nbig int result: %v\nscalar result %v", bigIntVal, modNVal, bigIntResult, modNValResult) } + + // Ensure algebraic properties with the half order hold. + checkModNScalarHalfOrderProps(t, modNVal) } } From f857a5ea14527b1323afa0c0a149de40deeb1336 Mon Sep 17 00:00:00 2001 From: Dave Collins Date: Thu, 9 Jul 2026 01:28:29 -0500 Subject: [PATCH 05/17] secp256k1: Implement 4x64 ModNScalar. MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit This significantly optimizes the ModNScalar implementation by rewriting it to use saturated 4x64 arithmetic similar to the new 4x64 field arithmetic. The new implementation retains all of the important properties such as being constant time and only using temporary buffers on the stack that are cleared. Of particular interest is the new modular reduction implementation. It moves away from the 96-bit accumulate design and instead manually tracks and discards carries inline which provides a major speed advantage. Due to the complexity and potential mistakes in the arithmetic, all of the intermediate bounds and carry assumptions documented that the code relies upon have been formally verified to be correct. The formal proof is in a separate commit. The implementation is primarily targeted at 64-bit since it is by far the most common architecture in modern use, however, as the included benchmarks below show, it is faster for almost all operations on 32-bit architectures as well. The few operations where it is a bit slower on 32-bit architectures are not used anywhere near as frequently as the ones that get 10-30% better performance, so the end result is the overall implementation is notably faster on 32 bit too. For 64-bit architectures, every operation is significantly faster. For example, all of the primary operations that get heavy use get ~73-94% better performance. Finally, the tests have been updated to use values specifically crafted to exercise the new limbs versus the old values that were crafted for the old limbs. Comparison for 32-bit platforms: $ benchstat beforescalar_x86.txt afterscalar_x86.txt name old time/op new time/op delta ----------------------------------------------------------------------------------- ModNScalar 22.9ns ± 4% 23.2ns ± 5% ~ (p=0.184 n=10+10) ModNScalarZero 1.36ns ± 3% 6.59ns ± 6% +385.95% (p=0.000 n=9+9) ModNScalarIsZero 0.32ns ± 4% 0.34ns ±14% ~ (p=0.143 n=10+10) ModNScalarEquals 3.69ns ± 4% 0.30ns ± 7% -91.78% (p=0.000 n=10+10) ModNScalarAdd 25.1ns ± 6% 22.5ns ± 5% -10.40% (p=0.000 n=9+10) ModNScalarMul 257ns ± 5% 223ns ± 5% -13.07% (p=0.000 n=10+10) ModNScalarSquare 267ns ± 5% 172ns ± 6% -35.56% (p=0.000 n=10+10) ModNScalarNegate 12.8ns ± 9% 14.5ns ± 2% +13.39% (p=0.000 n=10+10) ModNScalarInverse 686ns ± 7% 623ns ± 3% -9.24% (p=0.000 n=10+9) ModNScalarIsOverHalfOrder 8.92ns ± 6% 6.55ns ± 3% -26.59% (p=0.000 n=10+10) Comparison for 64-bit platforms: $ benchstat beforescalar_x64.txt afterscalar_x64.txt name old time/op new time/op delta ModNScalar 12.1ns ± 1% 3.5ns ± 4% -71.25% (p=0.000 n=7+10) ModNScalarZero 0.64ns ± 9% 0.61ns ± 5% -4.08% (p=0.043 n=10+10) ModNScalarIsZero 0.31ns ± 8% 0.31ns ± 6% ~ (p=0.000 n=10+10) ModNScalarEquals 2.63ns ± 5% 0.33ns ± 9% -87.47% (p=0.000 n=10+10) ModNScalarAdd 15.1ns ± 4% 4.1ns ± 4% -73.25% (p=0.000 n=10+10) ModNScalarMul 214ns ± 6% 29ns ± 6% -86.30% (p=0.000 n=10+10) ModNScalarSquare 215ns ± 4% 23ns ± 6% -89.39% (p=0.000 n=10+10) ModNScalarNegate 5.78ns ± 4% 2.99ns ±10% -48.33% (p=0.000 n=10+10) ModNScalarInverse 452ns ± 4% 408ns ± 6% -9.89% (p=0.000 n=10+10) ModNScalarIsOverHalfOrder 5.52ns ± 6% 0.29ns ± 5% -94.66% (p=0.000 n=10+9) --- dcrec/secp256k1/common.go | 7 +- dcrec/secp256k1/curve.go | 54 +- dcrec/secp256k1/curve_test.go | 20 +- dcrec/secp256k1/field_test.go | 4 +- dcrec/secp256k1/modnscalar.go | 1366 ++++++++++++---------------- dcrec/secp256k1/modnscalar_test.go | 473 +++++----- 6 files changed, 881 insertions(+), 1043 deletions(-) diff --git a/dcrec/secp256k1/common.go b/dcrec/secp256k1/common.go index b45e013bd5..9ff782bf40 100644 --- a/dcrec/secp256k1/common.go +++ b/dcrec/secp256k1/common.go @@ -15,9 +15,10 @@ func constantTimeEq64(a, b uint64) uint32 { return uint32(((t | -t) >> 63) ^ 1) } -// constantTimeNotEq returns 1 if a != b or 0 otherwise in constant time. -func constantTimeNotEq(a, b uint32) uint32 { - return ^uint32((uint64(a^b)-1)>>63) & 1 +// constantTimeNotEq64 returns 1 if a != b or 0 otherwise in constant time. +func constantTimeNotEq64(a, b uint64) uint32 { + t := a ^ b + return uint32(((t | -t) >> 63)) } // constantTimeLess returns 1 if a < b or 0 otherwise in constant time. diff --git a/dcrec/secp256k1/curve.go b/dcrec/secp256k1/curve.go index 84c6d0c414..750b57a327 100644 --- a/dcrec/secp256k1/curve.go +++ b/dcrec/secp256k1/curve.go @@ -720,16 +720,6 @@ func mulAdd64Carry(digit1, digit2, m, c uint64) (hi, lo uint64) { // because it is used for replacing division with multiplication and thus the // intermediate results must be done via a field extension to a larger field. func mul512Rsh320Round(n1, n2 *ModNScalar) ModNScalar { - // Convert n1 and n2 to base 2^64 digits. - n1Digit0 := uint64(n1.n[0]) | uint64(n1.n[1])<<32 - n1Digit1 := uint64(n1.n[2]) | uint64(n1.n[3])<<32 - n1Digit2 := uint64(n1.n[4]) | uint64(n1.n[5])<<32 - n1Digit3 := uint64(n1.n[6]) | uint64(n1.n[7])<<32 - n2Digit0 := uint64(n2.n[0]) | uint64(n2.n[1])<<32 - n2Digit1 := uint64(n2.n[2]) | uint64(n2.n[3])<<32 - n2Digit2 := uint64(n2.n[4]) | uint64(n2.n[5])<<32 - n2Digit3 := uint64(n2.n[6]) | uint64(n2.n[7])<<32 - // Compute the full 512-bit product n1*n2. var r0, r1, r2, r3, r4, r5, r6, r7, c uint64 @@ -738,40 +728,40 @@ func mul512Rsh320Round(n1, n2 *ModNScalar) ModNScalar { // // Note that r0 is ignored because it is not needed to compute the higher // terms and it is shifted out below anyway. - c, _ = bits.Mul64(n2Digit0, n1Digit0) - c, r1 = mulAdd64(n2Digit0, n1Digit1, c) - c, r2 = mulAdd64(n2Digit0, n1Digit2, c) - r4, r3 = mulAdd64(n2Digit0, n1Digit3, c) + c, _ = bits.Mul64(n2.n[0], n1.n[0]) + c, r1 = mulAdd64(n2.n[0], n1.n[1], c) + c, r2 = mulAdd64(n2.n[0], n1.n[2], c) + r4, r3 = mulAdd64(n2.n[0], n1.n[3], c) // Terms resulting from the product of the second digit of the second number // by all digits of the first number. // // Note that r1 is ignored because it is no longer needed to compute the // higher terms and it is shifted out below anyway. - c, _ = mulAdd64(n2Digit1, n1Digit0, r1) - c, r2 = mulAdd64Carry(n2Digit1, n1Digit1, r2, c) - c, r3 = mulAdd64Carry(n2Digit1, n1Digit2, r3, c) - r5, r4 = mulAdd64Carry(n2Digit1, n1Digit3, r4, c) + c, _ = mulAdd64(n2.n[1], n1.n[0], r1) + c, r2 = mulAdd64Carry(n2.n[1], n1.n[1], r2, c) + c, r3 = mulAdd64Carry(n2.n[1], n1.n[2], r3, c) + r5, r4 = mulAdd64Carry(n2.n[1], n1.n[3], r4, c) // Terms resulting from the product of the third digit of the second number // by all digits of the first number. // // Note that r2 is ignored because it is no longer needed to compute the // higher terms and it is shifted out below anyway. - c, _ = mulAdd64(n2Digit2, n1Digit0, r2) - c, r3 = mulAdd64Carry(n2Digit2, n1Digit1, r3, c) - c, r4 = mulAdd64Carry(n2Digit2, n1Digit2, r4, c) - r6, r5 = mulAdd64Carry(n2Digit2, n1Digit3, r5, c) + c, _ = mulAdd64(n2.n[2], n1.n[0], r2) + c, r3 = mulAdd64Carry(n2.n[2], n1.n[1], r3, c) + c, r4 = mulAdd64Carry(n2.n[2], n1.n[2], r4, c) + r6, r5 = mulAdd64Carry(n2.n[2], n1.n[3], r5, c) // Terms resulting from the product of the fourth digit of the second number // by all digits of the first number. // // Note that r3 is ignored because it is no longer needed to compute the // higher terms and it is shifted out below anyway. - c, _ = mulAdd64(n2Digit3, n1Digit0, r3) - c, r4 = mulAdd64Carry(n2Digit3, n1Digit1, r4, c) - c, r5 = mulAdd64Carry(n2Digit3, n1Digit2, r5, c) - r7, r6 = mulAdd64Carry(n2Digit3, n1Digit3, r6, c) + c, _ = mulAdd64(n2.n[3], n1.n[0], r3) + c, r4 = mulAdd64Carry(n2.n[3], n1.n[1], r4, c) + c, r5 = mulAdd64Carry(n2.n[3], n1.n[2], r5, c) + r7, r6 = mulAdd64Carry(n2.n[3], n1.n[3], r6, c) // At this point the upper 256 bits of the full 512-bit product n1*n2 are in // r4..r7 (recall the low order results were discarded as noted above). @@ -796,14 +786,10 @@ func mul512Rsh320Round(n1, n2 *ModNScalar) ModNScalar { // less than the group order given the group order is > 2^255 and the // maximum possible value of the result is 2^192. var result ModNScalar - result.n[0] = uint32(r0) - result.n[1] = uint32(r0 >> 32) - result.n[2] = uint32(r1) - result.n[3] = uint32(r1 >> 32) - result.n[4] = uint32(r2) - result.n[5] = uint32(r2 >> 32) - result.n[6] = uint32(r3) - result.n[7] = uint32(r3 >> 32) + result.n[0] = r0 + result.n[1] = r1 + result.n[2] = r2 + result.n[3] = r3 return result } diff --git a/dcrec/secp256k1/curve_test.go b/dcrec/secp256k1/curve_test.go index b1b4bec009..d3eb497f29 100644 --- a/dcrec/secp256k1/curve_test.go +++ b/dcrec/secp256k1/curve_test.go @@ -877,28 +877,16 @@ func TestScalarBaseMultJacobian(t *testing.T) { // modNBitLen returns the minimum number of bits required to represent the mod n // scalar. The result is 0 when the value is 0. func modNBitLen(s *ModNScalar) uint16 { - if w := s.n[7]; w > 0 { - return uint16(bits.Len32(w)) + 224 - } - if w := s.n[6]; w > 0 { - return uint16(bits.Len32(w)) + 192 - } - if w := s.n[5]; w > 0 { - return uint16(bits.Len32(w)) + 160 - } - if w := s.n[4]; w > 0 { - return uint16(bits.Len32(w)) + 128 - } if w := s.n[3]; w > 0 { - return uint16(bits.Len32(w)) + 96 + return uint16(bits.Len64(w)) + 192 } if w := s.n[2]; w > 0 { - return uint16(bits.Len32(w)) + 64 + return uint16(bits.Len64(w)) + 128 } if w := s.n[1]; w > 0 { - return uint16(bits.Len32(w)) + 32 + return uint16(bits.Len64(w)) + 64 } - return uint16(bits.Len32(s.n[0])) + return uint16(bits.Len64(s.n[0])) } // checkLambdaDecomposition returns an error if the provided decomposed scalars diff --git a/dcrec/secp256k1/field_test.go b/dcrec/secp256k1/field_test.go index a1200b8dd8..41a5a61757 100644 --- a/dcrec/secp256k1/field_test.go +++ b/dcrec/secp256k1/field_test.go @@ -1307,7 +1307,7 @@ func TestFieldMulInt(t *testing.T) { in2: 5, expected: "6cf4eb20f2447c77657fccb172d38c0aa91ea4ac446dc641fa463a6b5091fba7", }, { - name: "random sampling #3", + name: "random sampling #4", in1: "fb4529be3e027a3d1587d8a500b72f2d312e3577340ef5175f96d113be4c2ceb", in2: 8, expected: "da294df1f013d1e8ac3ec52805b979698971abb9a077a8bafcb688a4f261820f", @@ -1369,7 +1369,7 @@ func TestFieldMul(t *testing.T) { in2: "3", expected: "0", }, { - name: "secp256k1 prime * 3", + name: "secp256k1 prime * 8", in1: "fffffffffffffffffffffffffffffffffffffffffffffffffffffffefffffc2e", in2: "8", expected: "fffffffffffffffffffffffffffffffffffffffffffffffffffffffefffffc27", diff --git a/dcrec/secp256k1/modnscalar.go b/dcrec/secp256k1/modnscalar.go index 9405a987f0..6ccaf6306b 100644 --- a/dcrec/secp256k1/modnscalar.go +++ b/dcrec/secp256k1/modnscalar.go @@ -5,8 +5,10 @@ package secp256k1 import ( + "encoding/binary" "encoding/hex" "math/big" + "math/bits" "sync" ) @@ -18,81 +20,76 @@ import ( // https://cacr.uwaterloo.ca/hac/ // Many elliptic curve operations require working with scalars in a finite field -// characterized by the order of the group underlying the secp256k1 curve. -// Given this precision is larger than the biggest available native type, -// obviously some form of bignum math is needed. This code implements -// specialized fixed-precision field arithmetic rather than relying on an -// arbitrary-precision arithmetic package such as math/big for dealing with the -// math modulo the group order since the size is known. As a result, rather -// large performance gains are achieved by taking advantage of many -// optimizations not available to arbitrary-precision arithmetic and generic -// modular arithmetic algorithms. +// defined by the order of the group underlying the secp256k1 curve. Since the +// group order exceeds the width of the largest available native type, some form +// of multi-precision arithmetic is required. This code implements specialized +// fixed-precision field arithmetic rather than relying on an +// arbitrary-precision arithmetic package such as math/big for performing +// arithmetic modulo the fixed group order. As a result, significant +// performance gains are achieved by taking advantage of many optimizations not +// available to arbitrary-precision arithmetic and generic modular arithmetic +// algorithms. // -// There are various ways to internally represent each element. For example, -// the most obvious representation would be to use an array of 4 uint64s (64 -// bits * 4 = 256 bits). However, that representation suffers from the fact -// that there is no native Go type large enough to handle the intermediate -// results while adding or multiplying two 64-bit numbers. +// There are various ways to internally represent each element with their own +// pros and cons. Some common representations are: // -// Given the above, this implementation represents the field elements as 8 -// uint32s with each word (array entry) treated as base 2^32. This was chosen -// because most systems at the current time are 64-bit (or at least have 64-bit -// registers available for specialized purposes such as MMX) so the intermediate -// results can typically be done using a native register (and using uint64s to -// avoid the need for additional half-word arithmetic) - +// - 8 uint32s with base 2^32 limbs (aka 8x32: 32 bits * 8 = 256 bits) +// - 10 uint32s with base 2^26 limbs (aka 10x26: 26 bits * 10 = 260 bits) +// - 4 uint64s with base 2^64 limbs (aka 4x64: 64 bits * 4 = 256 bits) +// - 5 uint64s with base 2^52 limbs (aka 5x52: 52 bits * 5 = 260 bits) +// +// The primary tradeoff is complexity versus performance and the optimal choice +// also largely depends on the available hardware capabilities. +// +// For example, 5x52 and 10x26 allow performing several operations in a row +// before carry propagation or modular reduction becomes necessary since there +// are additional bits available in each limb, but that comes at the cost of +// manual magnitude tracking, multiple non-canonical representations of each +// number, and periodic normalization to propagate the accumulated carries and +// perform the modular reduction all at once. +// +// Conversely, 8x32 and 4x64 involve less complexity, are simpler to use, have +// fewer limbs to manage, and only have a single canonical representation for +// each number, but that comes at the cost of requiring carry propagation and +// modular reduction for every operation and larger intermediate results. +// +// The requirement for larger intermediate results is particularly notable for +// uint64s because there is no native Go type large enough to handle the 128-bit +// intermediate results while adding and multiplying the 64-bit limbs. +// +// This implementation historically used 8x32 for that reason, coupled with the +// fact that Go did not reliably make use of hardware-supported wide arithmetic, +// to ensure all intermediate results fit cleanly into uint64s. +// +// However, most systems are now 64-bit and almost all have instructions capable +// of handling the 128-bit intermediate results without the need for additional +// half-word arithmetic. Further, modern versions of Go automatically use those +// instructions on hardware that supports it. +// +// Given the above, this implementation represents the field elements as 4 +// uint64s with each limb (array entry) treated as base 2^64 (aka 4x64) and +// keeps the value reduced at all times. const ( - // These fields provide convenient access to each of the words of the + // These fields provide convenient access to each of the limbs of the // secp256k1 curve group order N to improve code readability. // // The group order of the curve per [SECG] is: - // 0xffffffff ffffffff ffffffff fffffffe baaedce6 af48a03b bfd25e8c d0364141 - // - // nolint: dupword - orderWordZero uint32 = 0xd0364141 - orderWordOne uint32 = 0xbfd25e8c - orderWordTwo uint32 = 0xaf48a03b - orderWordThree uint32 = 0xbaaedce6 - orderWordFour uint32 = 0xfffffffe - orderWordFive uint32 = 0xffffffff - orderWordSix uint32 = 0xffffffff - orderWordSeven uint32 = 0xffffffff - - // These fields provide convenient access to each of the words of the two's + // 0xfffffffffffffffffffffffffffffffebaaedce6af48a03bbfd25e8cd0364141 + orderLimb0 uint64 = 0xbfd25e8cd0364141 + orderLimb1 uint64 = 0xbaaedce6af48a03b + orderLimb2 uint64 = 0xfffffffffffffffe + orderLimb3 uint64 = 0xffffffffffffffff + + // These fields provide convenient access to each of the limbs of the two's // complement of the secp256k1 curve group order N to improve code // readability. // // The two's complement of the group order is: // 0x00000000 00000000 00000000 00000001 45512319 50b75fc4 402da173 2fc9bebf - orderComplementWordZero uint32 = (^orderWordZero) + 1 - orderComplementWordOne uint32 = ^orderWordOne - orderComplementWordTwo uint32 = ^orderWordTwo - orderComplementWordThree uint32 = ^orderWordThree - // orderComplementWordFour uint32 = ^orderWordFour // unused - // orderComplementWordFive uint32 = ^orderWordFive // unused - // orderComplementWordSix uint32 = ^orderWordSix // unused - // orderComplementWordSeven uint32 = ^orderWordSeven // unused - - // These fields provide convenient access to each of the words of the - // secp256k1 curve group order N / 2 to improve code readability and avoid - // the need to recalculate them. - // - // The half order of the secp256k1 curve group is: - // 0x7fffffff ffffffff ffffffff ffffffff 5d576e73 57a4501d dfe92f46 681b20a0 - // - // nolint: dupword - halfOrderWordZero uint32 = 0x681b20a0 - halfOrderWordOne uint32 = 0xdfe92f46 - halfOrderWordTwo uint32 = 0x57a4501d - halfOrderWordThree uint32 = 0x5d576e73 - halfOrderWordFour uint32 = 0xffffffff - halfOrderWordFive uint32 = 0xffffffff - halfOrderWordSix uint32 = 0xffffffff - halfOrderWordSeven uint32 = 0x7fffffff - - // uint32Mask is simply a mask with all bits set for a uint32 and is used to - // improve the readability of the code. - uint32Mask = 0xffffffff + orderComplementLimb0 uint64 = (^orderLimb0) + 1 + orderComplementLimb1 uint64 = ^orderLimb1 + orderComplementLimb2 uint64 = ^orderLimb2 + // orderComplementLimb3 uint64 = ^orderLimb3 // unused ) var ( @@ -112,39 +109,37 @@ var ( // around if absolutely needed. For example, subtraction can be performed by // adding the negation. // -// Should it be absolutely necessary, conversion to the standard library -// math/big.Int can be accomplished by using the Bytes method, slicing the -// resulting fixed-size array, and feeding it to big.Int.SetBytes. However, -// that should typically be avoided when possible as conversion to big.Ints -// requires allocations, is not constant time, and is slower when working modulo -// the group order. +// Should it be absolutely necessary, conversion to a standard library [big.Int] +// can be accomplished by using [ModNScalar.Bytes], slicing the resulting +// fixed-size array, and feeding it to [big.Int.SetBytes]. However, that should +// typically be avoided when possible as conversion to a [big.Int] requires +// allocations, is not constant time, and is slower when working modulo the +// group order. type ModNScalar struct { - // The scalar is represented as 8 32-bit integers in base 2^32. + // The scalar is represented as 4 64-bit integers in base 2^64. // // The following depicts the internal representation: - // --------------------------------------------------------- - // | n[7] | n[6] | ... | n[0] | - // | 32 bits | 32 bits | ... | 32 bits | - // | Mult: 2^(32*7) | Mult: 2^(32*6) | ... | Mult: 2^(32*0) | - // --------------------------------------------------------- + // -------------------------------------------------------------------- + // | n[3] | n[2] | n[1] | n[0] | + // | 64 bits | 64 bits | 64 bits | 64 bits | + // | Mult: 2^(64*3) | Mult: 2^(64*2) | Mult: 2^(64*1) | Mult: 2^(64*0) | + // -------------------------------------------------------------------- // - // For example, consider the number 2^87 + 2^42 + 1. It would be + // For example, consider the number 2^135 + 2^87 + 1. It would be // represented as: // n[0] = 1 - // n[1] = 2^10 - // n[2] = 2^23 - // n[3..7] = 0 - // - // The full 256-bit value is then calculated by looping i from 7..0 and - // doing sum(n[i] * 2^(32i)) like so: - // n[7] * 2^(32*7) = 0 * 2^224 = 0 - // n[6] * 2^(32*6) = 0 * 2^192 = 0 - // ... - // n[2] * 2^(32*2) = 2^23 * 2^64 = 2^87 - // n[1] * 2^(32*1) = 2^10 * 2^32 = 2^42 - // n[0] * 2^(32*0) = 1 * 2^0 = 1 - // Sum: 0 + 0 + ... + 2^87 + 2^42 + 1 = 2^87 + 2^42 + 1 - n [8]uint32 + // n[1] = 2^23 + // n[2] = 2^7 + // n[3] = 0 + // + // The full 256-bit value is then calculated by looping i from 3..0 and + // doing sum(n[i] * 2^(64i)) like so: + // n[3] * 2^(64*3) = 0 * 2^192 = 0 + // n[2] * 2^(64*2) = 2^7 * 2^128 = 2^135 + // n[1] * 2^(64*1) = 2^23 * 2^64 = 2^87 + // n[0] * 2^(64*0) = 1 * 2^0 = 1 + // Sum: 0 + 2^135 + 2^87 + 1 = 2^135 + 2^87 + 1 + n [4]uint64 } // String returns the scalar as a human-readable hex string. @@ -169,14 +164,7 @@ func (s *ModNScalar) Set(val *ModNScalar) *ModNScalar { // already set to zero. This function can be useful to clear an existing scalar // for reuse. func (s *ModNScalar) Zero() { - s.n[0] = 0 - s.n[1] = 0 - s.n[2] = 0 - s.n[3] = 0 - s.n[4] = 0 - s.n[5] = 0 - s.n[6] = 0 - s.n[7] = 0 + s.n = [4]uint64{} } // IsZeroBit returns 1 when the scalar is equal to zero or 0 otherwise in @@ -187,16 +175,13 @@ func (s *ModNScalar) Zero() { // operations require a numeric value. See [ModNScalar.IsZero] for the version // that returns a bool. func (s *ModNScalar) IsZeroBit() uint32 { - // The scalar can only be zero if no bits are set in any of the words. - bits := s.n[0] | s.n[1] | s.n[2] | s.n[3] | s.n[4] | s.n[5] | s.n[6] | s.n[7] - return constantTimeEq(bits, 0) + return constantTimeEq64(s.n[0]|s.n[1]|s.n[2]|s.n[3], 0) } // IsZero returns whether or not the scalar is equal to zero in constant time. func (s *ModNScalar) IsZero() bool { // The scalar can only be zero if no bits are set in any of the words. - bits := s.n[0] | s.n[1] | s.n[2] | s.n[3] | s.n[4] | s.n[5] | s.n[6] | s.n[7] - return bits == 0 + return (s.n[0] | s.n[1] | s.n[2] | s.n[3]) == 0 } // SetInt sets the scalar to the passed integer in constant time. This is a @@ -206,79 +191,10 @@ func (s *ModNScalar) IsZero() bool { // The scalar is returned to support chaining. This enables syntax like: // s := new(ModNScalar).SetInt(2).Mul(s2) so that s = 2 * s2. func (s *ModNScalar) SetInt(ui uint32) *ModNScalar { - s.Zero() - s.n[0] = ui + s.n = [4]uint64{uint64(ui), 0, 0, 0} return s } -// overflows determines if the current scalar is greater than or equal to the -// group order in constant time and returns 1 if it is or 0 otherwise. -func (s *ModNScalar) overflows() uint32 { - // The intuition here is that the scalar is greater than the group order if - // one of the higher individual words is greater than corresponding word of - // the group order and all higher words in the scalar are equal to their - // corresponding word of the group order. Since this type is modulo the - // group order, being equal is also an overflow back to 0. - // - // Note that the words 5, 6, and 7 are all the max uint32 value, so there is - // no need to test if those individual words of the scalar exceeds them, - // hence, only equality is checked for them. - highWordsEqual := constantTimeEq(s.n[7], orderWordSeven) - highWordsEqual &= constantTimeEq(s.n[6], orderWordSix) - highWordsEqual &= constantTimeEq(s.n[5], orderWordFive) - overflow := highWordsEqual & constantTimeGreater(s.n[4], orderWordFour) - highWordsEqual &= constantTimeEq(s.n[4], orderWordFour) - overflow |= highWordsEqual & constantTimeGreater(s.n[3], orderWordThree) - highWordsEqual &= constantTimeEq(s.n[3], orderWordThree) - overflow |= highWordsEqual & constantTimeGreater(s.n[2], orderWordTwo) - highWordsEqual &= constantTimeEq(s.n[2], orderWordTwo) - overflow |= highWordsEqual & constantTimeGreater(s.n[1], orderWordOne) - highWordsEqual &= constantTimeEq(s.n[1], orderWordOne) - overflow |= highWordsEqual & constantTimeGreaterOrEq(s.n[0], orderWordZero) - - return overflow -} - -// reduce256 reduces the current scalar modulo the group order in accordance -// with the overflows parameter in constant time. The overflows parameter -// specifies whether or not the scalar is known to be greater than the group -// order and MUST either be 1 in the case it is or 0 in the case it is not for a -// correct result. -func (s *ModNScalar) reduce256(overflows uint32) { - // Notice that since s < 2^256 < 2N (where N is the group order), the max - // possible number of reductions required is one. Therefore, in the case a - // reduction is needed, it can be performed with a single subtraction of N. - // Also, recall that subtraction is equivalent to addition by the two's - // complement while ignoring the carry. - // - // When s >= N, the overflows parameter will be 1. Conversely, it will be 0 - // when s < N. Thus multiplying by the overflows parameter will either - // result in 0 or the multiplicand itself. - // - // Combining the above along with the fact that s + 0 = s, the following is - // a constant time implementation that works by either adding 0 or the two's - // complement of N as needed. - // - // The final result will be in the range 0 <= s < N as expected. - overflows64 := uint64(overflows) - c := uint64(s.n[0]) + overflows64*uint64(orderComplementWordZero) - s.n[0] = uint32(c & uint32Mask) - c = (c >> 32) + uint64(s.n[1]) + overflows64*uint64(orderComplementWordOne) - s.n[1] = uint32(c & uint32Mask) - c = (c >> 32) + uint64(s.n[2]) + overflows64*uint64(orderComplementWordTwo) - s.n[2] = uint32(c & uint32Mask) - c = (c >> 32) + uint64(s.n[3]) + overflows64*uint64(orderComplementWordThree) - s.n[3] = uint32(c & uint32Mask) - c = (c >> 32) + uint64(s.n[4]) + overflows64 // * 1 - s.n[4] = uint32(c & uint32Mask) - c = (c >> 32) + uint64(s.n[5]) // + overflows64 * 0 - s.n[5] = uint32(c & uint32Mask) - c = (c >> 32) + uint64(s.n[6]) // + overflows64 * 0 - s.n[6] = uint32(c & uint32Mask) - c = (c >> 32) + uint64(s.n[7]) // + overflows64 * 0 - s.n[7] = uint32(c & uint32Mask) -} - // SetBytes interprets the provided array as a 256-bit big-endian unsigned // integer, reduces it modulo the group order, sets the scalar to the result, // and returns either 1 if it was reduced (aka it overflowed) or 0 otherwise in @@ -288,23 +204,41 @@ func (s *ModNScalar) reduce256(overflows uint32) { // from a bool to numeric value in constant time and many constant-time // operations require a numeric value. func (s *ModNScalar) SetBytes(b *[32]byte) uint32 { - // Pack the 256 total bits across the 8 uint32 words. This could be done + // Pack the 256 total bits across the 4 uint64 words. This could be done // with a for loop, but benchmarks show this unrolled version is about 2 // times faster than the variant that uses a loop. - s.n[0] = uint32(b[31]) | uint32(b[30])<<8 | uint32(b[29])<<16 | uint32(b[28])<<24 - s.n[1] = uint32(b[27]) | uint32(b[26])<<8 | uint32(b[25])<<16 | uint32(b[24])<<24 - s.n[2] = uint32(b[23]) | uint32(b[22])<<8 | uint32(b[21])<<16 | uint32(b[20])<<24 - s.n[3] = uint32(b[19]) | uint32(b[18])<<8 | uint32(b[17])<<16 | uint32(b[16])<<24 - s.n[4] = uint32(b[15]) | uint32(b[14])<<8 | uint32(b[13])<<16 | uint32(b[12])<<24 - s.n[5] = uint32(b[11]) | uint32(b[10])<<8 | uint32(b[9])<<16 | uint32(b[8])<<24 - s.n[6] = uint32(b[7]) | uint32(b[6])<<8 | uint32(b[5])<<16 | uint32(b[4])<<24 - s.n[7] = uint32(b[3]) | uint32(b[2])<<8 | uint32(b[1])<<16 | uint32(b[0])<<24 - - // The value might be >= N, so reduce it as required and return whether or - // not it was reduced. - needsReduce := s.overflows() - s.reduce256(needsReduce) - return needsReduce + s.n[0] = binary.BigEndian.Uint64(b[24:32]) + s.n[1] = binary.BigEndian.Uint64(b[16:24]) + s.n[2] = binary.BigEndian.Uint64(b[8:16]) + s.n[3] = binary.BigEndian.Uint64(b[0:8]) + + // Since s < 2^256 < 2N (where N is the secp256k1 group order), the max + // possible number of reductions required is one. Therefore, in the case a + // reduction is needed, it can be performed with a single subtraction of N. + // + // Since N must only conditionally be subtracted when s ≥ N, the following + // handles it in constant time by always calculating t = s - N and selecting + // the correct case via a constant time select. + + // Subtract N with borrow propagation. borrow is set iff s < N. + // + // In other words, the input overflowed (≥ N) when s - N does NOT borrow. + // + // t = s - N + var t0, t1, t2, t3, borrow uint64 + t0, borrow = bits.Sub64(s.n[0], orderLimb0, 0) + t1, borrow = bits.Sub64(s.n[1], orderLimb1, borrow) + t2, borrow = bits.Sub64(s.n[2], orderLimb2, borrow) + t3, borrow = bits.Sub64(s.n[3], orderLimb3, borrow) + + // Constant-time select. + // + // Set s = s when s < N (aka borrow is set). Otherwise s = t = s - N. + s.n[0] = constantTimeSelect64(borrow, s.n[0], t0) + s.n[1] = constantTimeSelect64(borrow, s.n[1], t1) + s.n[2] = constantTimeSelect64(borrow, s.n[2], t2) + s.n[3] = constantTimeSelect64(borrow, s.n[3], t3) + return uint32(1 - borrow) } // zeroArray32 zeroes the provided 32-byte buffer. @@ -324,6 +258,8 @@ func zeroArray32(b *[32]byte) { // size or it is acceptable to use this function with the described truncation // and overflow behavior. func (s *ModNScalar) SetByteSlice(b []byte) bool { + // Always copy a total of 32 bytes regardless of the input length to avoid + // introducing data-dependent timing. var b32 [32]byte b = b[:constantTimeMin(uint32(len(b)), 32)] copy(b32[:], b32[:32-len(b)]) @@ -337,49 +273,21 @@ func (s *ModNScalar) SetByteSlice(b []byte) bool { // into the passed byte slice in constant time. The target slice must have at // least 32 bytes available or it will panic. // -// There is a similar function, PutBytes, which unpacks the scalar into a -// 32-byte array directly. This version is provided since it can be useful to -// write directly into part of a larger buffer without needing a separate -// allocation. +// There is a similar function, [ModNScalar.PutBytes], which unpacks the scalar +// into a 32-byte array directly. This version is provided since it can be +// useful to write directly into part of a larger buffer without needing a +// separate allocation. // // Preconditions: // - The target slice MUST have at least 32 bytes available func (s *ModNScalar) PutBytesUnchecked(b []byte) { - // Unpack the 256 total bits from the 8 uint32 words. This could be done - // with a for loop, but benchmarks show this unrolled version is about 2 + // Unpack the 256 total bits from the 4 uint64 words. This could be done + // with a for loop, but benchmarks show this unrolled version is about 3 // times faster than the variant which uses a loop. - b[31] = byte(s.n[0]) - b[30] = byte(s.n[0] >> 8) - b[29] = byte(s.n[0] >> 16) - b[28] = byte(s.n[0] >> 24) - b[27] = byte(s.n[1]) - b[26] = byte(s.n[1] >> 8) - b[25] = byte(s.n[1] >> 16) - b[24] = byte(s.n[1] >> 24) - b[23] = byte(s.n[2]) - b[22] = byte(s.n[2] >> 8) - b[21] = byte(s.n[2] >> 16) - b[20] = byte(s.n[2] >> 24) - b[19] = byte(s.n[3]) - b[18] = byte(s.n[3] >> 8) - b[17] = byte(s.n[3] >> 16) - b[16] = byte(s.n[3] >> 24) - b[15] = byte(s.n[4]) - b[14] = byte(s.n[4] >> 8) - b[13] = byte(s.n[4] >> 16) - b[12] = byte(s.n[4] >> 24) - b[11] = byte(s.n[5]) - b[10] = byte(s.n[5] >> 8) - b[9] = byte(s.n[5] >> 16) - b[8] = byte(s.n[5] >> 24) - b[7] = byte(s.n[6]) - b[6] = byte(s.n[6] >> 8) - b[5] = byte(s.n[6] >> 16) - b[4] = byte(s.n[6] >> 24) - b[3] = byte(s.n[7]) - b[2] = byte(s.n[7] >> 8) - b[1] = byte(s.n[7] >> 16) - b[0] = byte(s.n[7] >> 24) + binary.BigEndian.PutUint64(b[0:8], s.n[3]) + binary.BigEndian.PutUint64(b[8:16], s.n[2]) + binary.BigEndian.PutUint64(b[16:24], s.n[1]) + binary.BigEndian.PutUint64(b[24:32], s.n[0]) } // PutBytes unpacks the scalar to a 32-byte big-endian value using the passed @@ -419,11 +327,8 @@ func (s *ModNScalar) IsOdd() bool { func (s *ModNScalar) Equals(val *ModNScalar) bool { // Xor only sets bits when they are different, so the two scalars can only // be the same if no bits are set after xoring each word. - bits := (s.n[0] ^ val.n[0]) | (s.n[1] ^ val.n[1]) | (s.n[2] ^ val.n[2]) | - (s.n[3] ^ val.n[3]) | (s.n[4] ^ val.n[4]) | (s.n[5] ^ val.n[5]) | - (s.n[6] ^ val.n[6]) | (s.n[7] ^ val.n[7]) - - return bits == 0 + return ((s.n[0] ^ val.n[0]) | (s.n[1] ^ val.n[1]) | (s.n[2] ^ val.n[2]) | + (s.n[3] ^ val.n[3])) == 0 } // Add2 adds the passed two scalars together modulo the group order in constant @@ -431,27 +336,42 @@ func (s *ModNScalar) Equals(val *ModNScalar) bool { // // The scalar is returned to support chaining. This enables syntax like: // s3.Add2(s, s2).AddInt(1) so that s3 = s + s2 + 1. -func (s *ModNScalar) Add2(val1, val2 *ModNScalar) *ModNScalar { - c := uint64(val1.n[0]) + uint64(val2.n[0]) - s.n[0] = uint32(c & uint32Mask) - c = (c >> 32) + uint64(val1.n[1]) + uint64(val2.n[1]) - s.n[1] = uint32(c & uint32Mask) - c = (c >> 32) + uint64(val1.n[2]) + uint64(val2.n[2]) - s.n[2] = uint32(c & uint32Mask) - c = (c >> 32) + uint64(val1.n[3]) + uint64(val2.n[3]) - s.n[3] = uint32(c & uint32Mask) - c = (c >> 32) + uint64(val1.n[4]) + uint64(val2.n[4]) - s.n[4] = uint32(c & uint32Mask) - c = (c >> 32) + uint64(val1.n[5]) + uint64(val2.n[5]) - s.n[5] = uint32(c & uint32Mask) - c = (c >> 32) + uint64(val1.n[6]) + uint64(val2.n[6]) - s.n[6] = uint32(c & uint32Mask) - c = (c >> 32) + uint64(val1.n[7]) + uint64(val2.n[7]) - s.n[7] = uint32(c & uint32Mask) - - // The result is now 256 bits, but it might still be >= N, so use the - // existing normal reduce method for 256-bit values. - s.reduce256(uint32(c>>32) + s.overflows()) +func (s *ModNScalar) Add2(a, b *ModNScalar) *ModNScalar { + // Since both values are already in the range 0 ≤ val < N (where N is the + // secp256k1 group order), the maximum possible result is < 2N - 1. So a + // maximum of one subtraction of N is required in the worst case. + // + // Since N must only conditionally be subtracted when a+b ≥ N, the following + // handles it in constant time by calculating both t = a+b and u = a+b - N + // and selecting the correct case via a constant time select. + + // Add with carry propagation. overflow is set iff t = a+b ≥ 2^256. + // + // t = a + b + var t0, t1, t2, t3, overflow, carry uint64 + t0, carry = bits.Add64(a.n[0], b.n[0], 0) + t1, carry = bits.Add64(a.n[1], b.n[1], carry) + t2, carry = bits.Add64(a.n[2], b.n[2], carry) + t3, overflow = bits.Add64(a.n[3], b.n[3], carry) + + // Subtract N with borrow propagation. borrow is set iff t = a+b < N. + // + // u = t - N = a+b - N + var u0, u1, u2, u3, borrow uint64 + u0, borrow = bits.Sub64(t0, orderLimb0, 0) + u1, borrow = bits.Sub64(t1, orderLimb1, borrow) + u2, borrow = bits.Sub64(t2, orderLimb2, borrow) + u3, borrow = bits.Sub64(t3, orderLimb3, borrow) + + // Constant-time select. + // + // Set s = t = a+b only when there was no overflow and t < N (borrow set). + // Otherwise s = u = a+b - N. + cond := (1 - overflow) & borrow + s.n[0] = constantTimeSelect64(cond, t0, u0) + s.n[1] = constantTimeSelect64(cond, t1, u1) + s.n[2] = constantTimeSelect64(cond, t2, u2) + s.n[3] = constantTimeSelect64(cond, t3, u3) return s } @@ -464,315 +384,359 @@ func (s *ModNScalar) Add(val *ModNScalar) *ModNScalar { return s.Add2(s, val) } -// accumulator96 provides a 96-bit accumulator for use in the intermediate -// calculations requiring more than 64-bits. -type accumulator96 struct { - n [3]uint32 -} +// scalar64Reduce512 reduces a 512-bit little-endian limb array modulo the group +// order in constant time and stores the result in r. +func scalar64Reduce512(r *[4]uint64, x *[8]uint64) { + // The overall strategy employed here is: + // 1) Start with the full unreduced 512-bit product of the two scalars. + // 2) Reduce the result modulo the group order via Crandall reduction + // (described below). + // 3) Repeat step 2 noting that each iteration reduces the required number + // of bits by 127 because the two's complement of N has 127 leading zero + // bits. + // 4) Once reduced to less than a max of 2N, perform the final reduction + // with a constant-time conditional subtraction. + // + // Technically the max possible input value here will be (N-1)^2, since the + // only time it is called is with the product of two scalars that are always + // mod N. Nevertheless, it is safer to treat it as the product of two max + // size 256-bit values, (2^256-1)^2 = 2^512 - 2^257 + 1, to ensure + // correctness if that were to ever change. + // + // Per [HAC] section 14.3.4: Reduction method of moduli of special form, + // when the modulus is of the special form m = b^t - c, highly efficient + // reduction can be achieved. While [HAC] only presents the algorithm and + // does not call it out by name or provide the mathematical justification, + // the underlying technique is known as Crandall reduction and is often + // presented as 2^k - c. It is easy to see they are equivalent by setting + // b = 2 and t = k. + // + // The secp256k1 group order is 2^256 - 0x14551231950b75fc4402da1732fc9bebf, + // so it fits this criteria where: + // k = 256 + // c = 432420386565659656852420866394968145599 + // + // Crandall reduction works by taking advantage of the fact that if a prime + // is of the form 2^k - c, then 2^k - c ≡ 0 (mod p), so 2^k ≡ c (mod p). In + // other words, every multiple of 2^k is equivalent to adding c when working + // modulo p. + // + // Since the 512-bit value to reduce is tightly packed into uint64s, the + // upper 4 limbs are all multiples of 2^256. Therefore, reducing modulo the + // group order is equivalent to multiplying those upper limbs by c and + // adding the result to the corresponding lower 4 limbs while propagating + // the carries. + // + // For the specific case of the secp256k1 group order, a max of 4 reductions + // are required because c is 129 bits and so the first round will reduce + // from 512 bits to a max of 385 bits, the second round will reduce to a max + // of 258 bits, and the third round will reduce to within 2N. Then, a + // conditional subtraction of N handles the final reduction. + // + // The reduction is slightly complicated by the fact c needs 3 uint64 limbs + // to represent it, so multiplying the upper limbs by c requires multiplying + // each one by all 3 limbs of c and carefully managing the columns and + // carries. The uppermost limb of c = 1, so multiplication by that limb is + // avoided in the code below. -// Add adds the passed unsigned 64-bit value to the accumulator. -func (a *accumulator96) Add(v uint64) { - low := uint32(v & uint32Mask) - hi := uint32(v >> 32) - a.n[0] += low - hi += constantTimeLess(a.n[0], low) // Carry if overflow in n[0]. - a.n[1] += hi - a.n[2] += constantTimeLess(a.n[1], hi) // Carry if overflow in n[1]. -} + var t0, t1, t2, t3, t4, t5, t6, l0, h0, l1, h1, l2, carry uint64 -// Rsh32 right shifts the accumulator by 32 bits. -func (a *accumulator96) Rsh32() { - a.n[0] = a.n[1] - a.n[1] = a.n[2] - a.n[2] = 0 -} + // ------------------------------------------------------------------------- + // First reduction. + // + // The code in this section is equivalent to the following if uint512 and + // larger numbers and shifts were supported directly: + // + // t0 = uint64(x%(1<<256) + (x>>256)*c) + // t1 = uint64((x%(1<<256) + (x>>256)*c) >> 64) + // t2 = uint64((x%(1<<256) + (x>>256)*c) >> 128) + // t3 = uint64((x%(1<<256) + (x>>256)*c) >> 192) + // t4 = uint64((x%(1<<256) + (x>>256)*c) >> 256) + // t5 = uint64((x%(1<<256) + (x>>256)*c) >> 320) + // t6 = uint64((x%(1<<256) + (x>>256)*c) >> 384) + // + // To visualize the logic, multiply the 3 limbs of c by the upper limbs via + // the standard pencil-and-paper method: + // + // 2^576 2^512 2^448 2^384 2^320 2^256 (place values) + // ----------------------------------- + // x7 x6 x5 x4 + // c2 c1 c0 * + // ----------------------------------- + // x4c2 x4c1 x4c0 + // x5c2 x5c1 x5c0 + + // x6c2 x6c1 x6c0 + + // x7c2 x7c1 x7c0 + + // ----------------------------------- + // x7c2 .... .... .... .... x4c0 + // + // Then add them to the associated lower limbs per the reduction identity: + // + // 2^384 2^320 2^256 2^192 2^128 2^64 2^0 (place values) + // ----------------------------------------- + // x3 x2 x1 x0 + // x4c2 x4c1 x4c0 + + // x5c2 x5c1 x5c0 + + // x6c2 x6c1 x6c0 + + // x7c2 x7c1 x7c0 + + // ----------------------------------------- + // t6 t5 t4 t3 t2 t1 t0 + // + // Where t6 ≤ 1 for the potential carry. + // + // Thus 0 ≤ t < 2^385 after the 1st reduction. + // ------------------------------------------------------------------------- + + // Compute x4*c with result in t0..t2 with carries propagated and the final + // carry in t3. + h0, t0 = bits.Mul64(x[4], orderComplementLimb0) + h1, t1 = bits.Mul64(x[4], orderComplementLimb1) + t1, carry = bits.Add64(t1, h0, 0) + t2, carry = bits.Add64(x[4], h1, carry) + t3 = carry + + // Compute x5*c with result added to t1..t3 with carries propagated and the + // final carry in t4. + // + // Note that since h1 is the upper 64 bits of the product of a uint64 with + // c1 (limb 1 of c) and c1 < 2^63 - 1: + // h1 ≤ floor((2^64-1)(2^63-1) / 2^64) = 2^63 - 2 + // + // Then, because current t3 ≤ 1 and carry ≤ 1, a loose bound for the first + // t3 below is: + // t3 ≤ h1 + 2 = 2^63 < 2^64 + // + // Therefore, it is safe to discard the carry. + h0, l0 = bits.Mul64(x[5], orderComplementLimb0) + h1, l1 = bits.Mul64(x[5], orderComplementLimb1) + t1, carry = bits.Add64(t1, l0, 0) + t2, carry = bits.Add64(t2, h0, carry) + t3, _ = bits.Add64(t3, h1, carry) + t2, carry = bits.Add64(t2, l1, 0) + t3, carry = bits.Add64(t3, x[5], carry) + t4 = carry + + // Compute x6*c with result added to t2..t4 with carries propagated and the + // final carry in t5. + // + // It is safe to discard the carry on the first t4 for the same reason as + // above. + h0, l0 = bits.Mul64(x[6], orderComplementLimb0) + h1, l1 = bits.Mul64(x[6], orderComplementLimb1) + t2, carry = bits.Add64(t2, l0, 0) + t3, carry = bits.Add64(t3, h0, carry) + t4, _ = bits.Add64(t4, h1, carry) + t3, carry = bits.Add64(t3, l1, 0) + t4, carry = bits.Add64(t4, x[6], carry) + t5 = carry + + // Compute x7*c with result added to t3..t5 with carries propagated and the + // final carry in t6. + // + // It is safe to discard the carry on the first t5 for the same reason as + // above. + h0, l0 = bits.Mul64(x[7], orderComplementLimb0) + h1, l1 = bits.Mul64(x[7], orderComplementLimb1) + t3, carry = bits.Add64(t3, l0, 0) + t4, carry = bits.Add64(t4, h0, carry) + t5, _ = bits.Add64(t5, h1, carry) + t4, carry = bits.Add64(t4, l1, 0) + t5, carry = bits.Add64(t5, x[7], carry) + t6 = carry + + // Add result to lower limbs and propagate carries. + // + // It is safe to discard the carry on t6 since the resulting t6 ≤ 1 as + // previously described. + t0, carry = bits.Add64(t0, x[0], 0) + t1, carry = bits.Add64(t1, x[1], carry) + t2, carry = bits.Add64(t2, x[2], carry) + t3, carry = bits.Add64(t3, x[3], carry) + t4, carry = bits.Add64(t4, 0, carry) + t5, carry = bits.Add64(t5, 0, carry) + t6, _ = bits.Add64(t6, 0, carry) + + // ------------------------------------------------------------------------- + // Second reduction. + // + // The value now fits in 385 bits, so reduce it again. Only t4, t5, and t6 + // need to be considered since the higher limb, t7, is ≥ 448 bits and thus + // guaranteed to be 0. Further, as previously noted, t6 ≤ 1. + // + // The code in this section is equivalent to following if uint512 and larger + // numbers and shifts were supported: + // + // t0 = uint64(t%2^256 + (t>>256)*c) + // t1 = uint64((t%2^256 + (t>>256)*c) >> 64) + // t2 = uint64((t%2^256 + (t>>256)*c) >> 128) + // t3 = uint64((t%2^256 + (t>>256)*c) >> 192) + // t4 = uint64((t%2^256 + (t>>256)*c) >> 256) + // + // To visualize the logic, multiply the 3 limbs of c by the new upper limbs + // via the standard pencil-and-paper method: + // + // 2^512 2^448 2^384 2^320 2^256 (place values) + // ----------------------------- + // t6 t5 t4 + // c2 c1 c0 * + // ----------------------------- + // t4c2 t4c1 t4c0 + // t5c2 t5c1 t5c0 + + // t6c2 t6c1 t6c0 + + // ---------------------------- + // t6c2 .... .... .... t4c0 + // + // Then add them to the associated lower limbs per the reduction identity + // noting that the code reuses t for the sums: + // + // 2^256 2^192 2^128 2^64 2^0 (place values) + // ----------------------------- + // t3 t2 t1 t0 + // t4c2 t4c1 t4c0 + + // t5c2 t5c1 t5c0 + + // t6c2 t6c1 t6c0 + + // ----------------------------- + // t4 t3 t2 t1 t0 + // + // With a loose bound for t4 ≤ 3. + // + // Proof: + // + // Known: c2=1, t6 ≤ 1, t3 ≤ 2^64-1, t5 ≤ 2^64-1, and c1 ≤ 2^63-1. + // + // Let s = t3 + t5c2 + t6c1. Then: + // s ≤ 2^64-1 + (2^64-1)*1 + 1*(2^63-1) = 2^65 + 2^63 - 3 + // + // So, the loose bound on the carry in to t4 is carryIn ≤ floor(s/2^64) = 2. + // + // Finally, t4 ≤ t6*c2 + carryIn ≤ 1*1 + 2 ≤ 3. + // + // Thus 0 ≤ t < 2^258 after the 2nd reduction. + // ------------------------------------------------------------------------- -// reduce385 reduces the 385-bit intermediate result in the passed terms modulo -// the group order in constant time and stores the result in s. -func (s *ModNScalar) reduce385(t0, t1, t2, t3, t4, t5, t6, t7, t8, t9, t10, t11, t12 uint64) { - // At this point, the intermediate result in the passed terms has been - // reduced to fit within 385 bits, so reduce it again using the same method - // described in reduce512. As before, the intermediate result will end up - // being reduced by another 127 bits to 258 bits, thus 9 32-bit terms are - // needed for this iteration. The reduced terms are assigned back to t0 - // through t8. - // - // Note that several of the intermediate calculations require adding 64-bit - // products together which would overflow a uint64, so a 96-bit accumulator - // is used instead until the value is reduced enough to use native uint64s. - - // Terms for 2^(32*0). - var acc accumulator96 - acc.n[0] = uint32(t0) // == acc.Add(t0) because acc is guaranteed to be 0. - acc.Add(t8 * uint64(orderComplementWordZero)) - t0 = uint64(acc.n[0]) - acc.Rsh32() - - // Terms for 2^(32*1). - acc.Add(t1) - acc.Add(t8 * uint64(orderComplementWordOne)) - acc.Add(t9 * uint64(orderComplementWordZero)) - t1 = uint64(acc.n[0]) - acc.Rsh32() - - // Terms for 2^(32*2). - acc.Add(t2) - acc.Add(t8 * uint64(orderComplementWordTwo)) - acc.Add(t9 * uint64(orderComplementWordOne)) - acc.Add(t10 * uint64(orderComplementWordZero)) - t2 = uint64(acc.n[0]) - acc.Rsh32() - - // Terms for 2^(32*3). - acc.Add(t3) - acc.Add(t8 * uint64(orderComplementWordThree)) - acc.Add(t9 * uint64(orderComplementWordTwo)) - acc.Add(t10 * uint64(orderComplementWordOne)) - acc.Add(t11 * uint64(orderComplementWordZero)) - t3 = uint64(acc.n[0]) - acc.Rsh32() - - // Terms for 2^(32*4). - acc.Add(t4) - acc.Add(t8) // * uint64(orderComplementWordFour) // * 1 - acc.Add(t9 * uint64(orderComplementWordThree)) - acc.Add(t10 * uint64(orderComplementWordTwo)) - acc.Add(t11 * uint64(orderComplementWordOne)) - acc.Add(t12 * uint64(orderComplementWordZero)) - t4 = uint64(acc.n[0]) - acc.Rsh32() - - // Terms for 2^(32*5). - acc.Add(t5) - // acc.Add(t8 * uint64(orderComplementWordFive)) // 0 - acc.Add(t9) // * uint64(orderComplementWordFour) // * 1 - acc.Add(t10 * uint64(orderComplementWordThree)) - acc.Add(t11 * uint64(orderComplementWordTwo)) - acc.Add(t12 * uint64(orderComplementWordOne)) - t5 = uint64(acc.n[0]) - acc.Rsh32() - - // Terms for 2^(32*6). - acc.Add(t6) - // acc.Add(t8 * uint64(orderComplementWordSix)) // 0 - // acc.Add(t9 * uint64(orderComplementWordFive)) // 0 - acc.Add(t10) // * uint64(orderComplementWordFour) // * 1 - acc.Add(t11 * uint64(orderComplementWordThree)) - acc.Add(t12 * uint64(orderComplementWordTwo)) - t6 = uint64(acc.n[0]) - acc.Rsh32() - - // Terms for 2^(32*7). - acc.Add(t7) - // acc.Add(t8 * uint64(orderComplementWordSeven)) // 0 - // acc.Add(t9 * uint64(orderComplementWordSix)) // 0 - // acc.Add(t10 * uint64(orderComplementWordFive)) // 0 - acc.Add(t11) // * uint64(orderComplementWordFour) // * 1 - acc.Add(t12 * uint64(orderComplementWordThree)) - t7 = uint64(acc.n[0]) - acc.Rsh32() - - // Terms for 2^(32*8). - // acc.Add(t9 * uint64(orderComplementWordSeven)) // 0 - // acc.Add(t10 * uint64(orderComplementWordSix)) // 0 - // acc.Add(t11 * uint64(orderComplementWordFive)) // 0 - acc.Add(t12) // * uint64(orderComplementWordFour) // * 1 - t8 = uint64(acc.n[0]) - // acc.Rsh32() // No need since not used after this. Guaranteed to be 0. - - // NOTE: All of the remaining multiplications for this iteration result in 0 - // as they all involve multiplying by combinations of the fifth, sixth, and - // seventh words of the two's complement of N, which are 0, so skip them. - - // At this point, the result is reduced to fit within 258 bits, so reduce it - // again using a slightly modified version of the same method. The maximum - // value in t8 is 2 at this point and therefore multiplying it by each word - // of the two's complement of N and adding it to a 32-bit term will result - // in a maximum requirement of 33 bits, so it is safe to use native uint64s - // here for the intermediate term carry propagation. - // - // Also, since the maximum value in t8 is 2, this ends up reducing by - // another 2 bits to 256 bits. - c := t0 + t8*uint64(orderComplementWordZero) - s.n[0] = uint32(c & uint32Mask) - c = (c >> 32) + t1 + t8*uint64(orderComplementWordOne) - s.n[1] = uint32(c & uint32Mask) - c = (c >> 32) + t2 + t8*uint64(orderComplementWordTwo) - s.n[2] = uint32(c & uint32Mask) - c = (c >> 32) + t3 + t8*uint64(orderComplementWordThree) - s.n[3] = uint32(c & uint32Mask) - c = (c >> 32) + t4 + t8 // * uint64(orderComplementWordFour) == * 1 - s.n[4] = uint32(c & uint32Mask) - c = (c >> 32) + t5 // + t8*uint64(orderComplementWordFive) == 0 - s.n[5] = uint32(c & uint32Mask) - c = (c >> 32) + t6 // + t8*uint64(orderComplementWordSix) == 0 - s.n[6] = uint32(c & uint32Mask) - c = (c >> 32) + t7 // + t8*uint64(orderComplementWordSeven) == 0 - s.n[7] = uint32(c & uint32Mask) - - // The result is now 256 bits, but it might still be >= N, so use the - // existing normal reduce method for 256-bit values. - s.reduce256(uint32(c>>32) + s.overflows()) -} + // Compute t4*c with result added to t0..t2 with carries propagated and the + // final carry in t4. + // + // It is safe to discard the carry on the final t4 given the carries ≤ 1. + h0, l0 = bits.Mul64(t4, orderComplementLimb0) + h1, l1 = bits.Mul64(t4, orderComplementLimb1) + l2 = t4 // h2, l2 = bits.Mul64(t4, orderComplementLimb2) => h2=0, l2=t4 + t0, carry = bits.Add64(t0, l0, 0) + t1, carry = bits.Add64(t1, h0, carry) + t2, carry = bits.Add64(t2, h1, carry) + t3, carry = bits.Add64(t3, 0, carry) + t4 = carry + t1, carry = bits.Add64(t1, l1, 0) + t2, carry = bits.Add64(t2, l2, carry) + t3, carry = bits.Add64(t3, 0, carry) + t4, _ = bits.Add64(t4, 0, carry) + + // Compute t5*c with result added to t1..t3 with carries propagated. + // + // It is safe to discard the carry on the final t4 given the carries ≤ 1. + h0, l0 = bits.Mul64(t5, orderComplementLimb0) + h1, l1 = bits.Mul64(t5, orderComplementLimb1) + // h2, l2 = bits.Mul64(t5, orderComplementLimb2) => h2=0, l2=t5 + t1, carry = bits.Add64(t1, l0, 0) + t2, carry = bits.Add64(t2, h0, carry) + t3, carry = bits.Add64(t3, h1, carry) + t4, _ = bits.Add64(t4, 0, carry) + t2, carry = bits.Add64(t2, l1, 0) + t3, carry = bits.Add64(t3, t5, carry) + t4, _ = bits.Add64(t4, 0, carry) + + // Compute t6*c with result added to t2..t4 with carries propagated. + // + // Note that since t6 ≤ 1, the product can't overflow a uint64, so it is + // safe to use normal 64-bit multiplication. + // + // It is safe to discard the carry on t4 since the resulting t4 ≤ 3 as + // previously described. + t2, carry = bits.Add64(t2, t6*orderComplementLimb0, 0) + t3, carry = bits.Add64(t3, t6*orderComplementLimb1, carry) + t4, _ = bits.Add64(t4, t6, carry) + + // ------------------------------------------------------------------------- + // Third reduction. + // + // The value now fits in 258 bits, so reduce it again. Only t4 needs to be + // considered since the higher limbs are ≥ 320 bits and thus guaranteed to + // be 0. + // + // The code in this section is equivalent to following if uint512 and larger + // numbers and shifts were supported: + // + // t0 = uint64(t%2^256 + (t>>256)*c) + // t1 = uint64((t%2^256 + (t>>256)*c) >> 64) + // t2 = uint64((t%2^256 + (t>>256)*c) >> 128) + // t3 = uint64((t%2^256 + (t>>256)*c) >> 192) + // t4 = uint64((t%2^256 + (t>>256)*c) >> 256) + // + // To visualize the logic, multiply the 3 limbs of c by the new upper limb + // via the standard pencil-and-paper method: + // + // 2^512 2^448 2^384 2^320 2^256 (place values) + // ----------------------------- + // t4 + // c2 c1 c0 * + // ----------------------------- + // t4c2 t4c1 t4c0 + // + // Then add them to the associated lower limbs noting that the code reuses t + // for the sums: + // + // 2^256 2^192 2^128 2^64 2^0 (place values) + // ----------------------------- + // t3 t2 t1 t0 + // t4c2 t4c1 t4c0 + + // ----------------------------- + // t4 t3 t2 t1 t0 + // + // Where t4 ≤ 1 for the carry. + // ------------------------------------------------------------------------- -// reduce512 reduces the 512-bit intermediate result in the passed terms modulo -// the group order down to 385 bits in constant time and stores the result in s. -func (s *ModNScalar) reduce512(t0, t1, t2, t3, t4, t5, t6, t7, t8, t9, t10, t11, t12, t13, t14, t15 uint64) { - // At this point, the intermediate result in the passed terms is grouped - // into the respective bases. + // Compute t4*c with result added to t0..t2 with carries propagated and the + // final carry in t4. // - // Per [HAC] section 14.3.4: Reduction method of moduli of special form, - // when the modulus is of the special form m = b^t - c, where log_2(c) < t, - // highly efficient reduction can be achieved per the provided algorithm. - // - // The secp256k1 group order fits this criteria since it is: - // 2^256 - 432420386565659656852420866394968145599 - // - // Technically the max possible value here is (N-1)^2 since the two scalars - // being multiplied are always mod N. Nevertheless, it is safer to consider - // it to be (2^256-1)^2 = 2^512 - 2^257 + 1 since it is the product of two - // 256-bit values. - // - // The algorithm is to reduce the result modulo the prime by subtracting - // multiples of the group order N. However, in order simplify carry - // propagation, this adds with the two's complement of N to achieve the same - // result. - // - // Since the two's complement of N has 127 leading zero bits, this will end - // up reducing the intermediate result from 512 bits to 385 bits, resulting - // in 13 32-bit terms. The reduced terms are assigned back to t0 through - // t12. - // - // Note that several of the intermediate calculations require adding 64-bit - // products together which would overflow a uint64, so a 96-bit accumulator - // is used instead. - - // Terms for 2^(32*0). - var acc accumulator96 - acc.n[0] = uint32(t0) // == acc.Add(t0) because acc is guaranteed to be 0. - acc.Add(t8 * uint64(orderComplementWordZero)) - t0 = uint64(acc.n[0]) - acc.Rsh32() - - // Terms for 2^(32*1). - acc.Add(t1) - acc.Add(t8 * uint64(orderComplementWordOne)) - acc.Add(t9 * uint64(orderComplementWordZero)) - t1 = uint64(acc.n[0]) - acc.Rsh32() - - // Terms for 2^(32*2). - acc.Add(t2) - acc.Add(t8 * uint64(orderComplementWordTwo)) - acc.Add(t9 * uint64(orderComplementWordOne)) - acc.Add(t10 * uint64(orderComplementWordZero)) - t2 = uint64(acc.n[0]) - acc.Rsh32() - - // Terms for 2^(32*3). - acc.Add(t3) - acc.Add(t8 * uint64(orderComplementWordThree)) - acc.Add(t9 * uint64(orderComplementWordTwo)) - acc.Add(t10 * uint64(orderComplementWordOne)) - acc.Add(t11 * uint64(orderComplementWordZero)) - t3 = uint64(acc.n[0]) - acc.Rsh32() - - // Terms for 2^(32*4). - acc.Add(t4) - acc.Add(t8) // * uint64(orderComplementWordFour) // * 1 - acc.Add(t9 * uint64(orderComplementWordThree)) - acc.Add(t10 * uint64(orderComplementWordTwo)) - acc.Add(t11 * uint64(orderComplementWordOne)) - acc.Add(t12 * uint64(orderComplementWordZero)) - t4 = uint64(acc.n[0]) - acc.Rsh32() - - // Terms for 2^(32*5). - acc.Add(t5) - // acc.Add(t8 * uint64(orderComplementWordFive)) // 0 - acc.Add(t9) // * uint64(orderComplementWordFour) // * 1 - acc.Add(t10 * uint64(orderComplementWordThree)) - acc.Add(t11 * uint64(orderComplementWordTwo)) - acc.Add(t12 * uint64(orderComplementWordOne)) - acc.Add(t13 * uint64(orderComplementWordZero)) - t5 = uint64(acc.n[0]) - acc.Rsh32() - - // Terms for 2^(32*6). - acc.Add(t6) - // acc.Add(t8 * uint64(orderComplementWordSix)) // 0 - // acc.Add(t9 * uint64(orderComplementWordFive)) // 0 - acc.Add(t10) // * uint64(orderComplementWordFour)) // * 1 - acc.Add(t11 * uint64(orderComplementWordThree)) - acc.Add(t12 * uint64(orderComplementWordTwo)) - acc.Add(t13 * uint64(orderComplementWordOne)) - acc.Add(t14 * uint64(orderComplementWordZero)) - t6 = uint64(acc.n[0]) - acc.Rsh32() - - // Terms for 2^(32*7). - acc.Add(t7) - // acc.Add(t8 * uint64(orderComplementWordSeven)) // 0 - // acc.Add(t9 * uint64(orderComplementWordSix)) // 0 - // acc.Add(t10 * uint64(orderComplementWordFive)) // 0 - acc.Add(t11) // * uint64(orderComplementWordFour) // * 1 - acc.Add(t12 * uint64(orderComplementWordThree)) - acc.Add(t13 * uint64(orderComplementWordTwo)) - acc.Add(t14 * uint64(orderComplementWordOne)) - acc.Add(t15 * uint64(orderComplementWordZero)) - t7 = uint64(acc.n[0]) - acc.Rsh32() - - // Terms for 2^(32*8). - // acc.Add(t9 * uint64(orderComplementWordSeven)) // 0 - // acc.Add(t10 * uint64(orderComplementWordSix)) // 0 - // acc.Add(t11 * uint64(orderComplementWordFive)) // 0 - acc.Add(t12) // * uint64(orderComplementWordFour) // * 1 - acc.Add(t13 * uint64(orderComplementWordThree)) - acc.Add(t14 * uint64(orderComplementWordTwo)) - acc.Add(t15 * uint64(orderComplementWordOne)) - t8 = uint64(acc.n[0]) - acc.Rsh32() - - // Terms for 2^(32*9). - // acc.Add(t10 * uint64(orderComplementWordSeven)) // 0 - // acc.Add(t11 * uint64(orderComplementWordSix)) // 0 - // acc.Add(t12 * uint64(orderComplementWordFive)) // 0 - acc.Add(t13) // * uint64(orderComplementWordFour) // * 1 - acc.Add(t14 * uint64(orderComplementWordThree)) - acc.Add(t15 * uint64(orderComplementWordTwo)) - t9 = uint64(acc.n[0]) - acc.Rsh32() - - // Terms for 2^(32*10). - // acc.Add(t11 * uint64(orderComplementWordSeven)) // 0 - // acc.Add(t12 * uint64(orderComplementWordSix)) // 0 - // acc.Add(t13 * uint64(orderComplementWordFive)) // 0 - acc.Add(t14) // * uint64(orderComplementWordFour) // * 1 - acc.Add(t15 * uint64(orderComplementWordThree)) - t10 = uint64(acc.n[0]) - acc.Rsh32() - - // Terms for 2^(32*11). - // acc.Add(t12 * uint64(orderComplementWordSeven)) // 0 - // acc.Add(t13 * uint64(orderComplementWordSix)) // 0 - // acc.Add(t14 * uint64(orderComplementWordFive)) // 0 - acc.Add(t15) // * uint64(orderComplementWordFour) // * 1 - t11 = uint64(acc.n[0]) - acc.Rsh32() - - // NOTE: All of the remaining multiplications for this iteration result in 0 - // as they all involve multiplying by combinations of the fifth, sixth, and - // seventh words of the two's complement of N, which are 0, so skip them. - - // Terms for 2^(32*12). - t12 = uint64(acc.n[0]) - // acc.Rsh32() // No need since not used after this. Guaranteed to be 0. - - // At this point, the result is reduced to fit within 385 bits, so reduce it - // again using the same method accordingly. - s.reduce385(t0, t1, t2, t3, t4, t5, t6, t7, t8, t9, t10, t11, t12) + // Note that since current t4 ≤ 3, these bounds for the product hold: + // t4*c0 < 2^64 + // t4*c1 < 2^64 + // + // So, uint64 overflow is impossible and therefore it is safe to use normal + // 64-bit multiplication. + t0, carry = bits.Add64(t0, t4*orderComplementLimb0, 0) + t1, carry = bits.Add64(t1, t4*orderComplementLimb1, carry) + t2, carry = bits.Add64(t2, t4, carry) + t3, carry = bits.Add64(t3, 0, carry) + t4 = carry + + // ------------------------------------------------------------------------- + // Final reduction. + // + // The value is now in the range 0 ≤ t < 2N, so one 5-limb conditional + // subtract of N guarantees it is fully reduced. + // ------------------------------------------------------------------------- + + // Subtract N with borrow propagation. borrow is set iff t < N. + // + // In other words, the value overflowed (≥ N) when t - N does NOT borrow. + // + // s = t - N + var s0, s1, s2, s3, borrow uint64 + s0, borrow = bits.Sub64(t0, orderLimb0, 0) + s1, borrow = bits.Sub64(t1, orderLimb1, borrow) + s2, borrow = bits.Sub64(t2, orderLimb2, borrow) + s3, borrow = bits.Sub64(t3, orderLimb3, borrow) + _, borrow = bits.Sub64(t4, 0, borrow) + + // Constant-time select. + // + // Set r = t only when t < N (borrow set). + // Otherwise r = s = t - N. + r[0] = constantTimeSelect64(borrow, t0, s0) + r[1] = constantTimeSelect64(borrow, t1, s1) + r[2] = constantTimeSelect64(borrow, t2, s2) + r[3] = constantTimeSelect64(borrow, t3, s3) } // Mul2 multiplies the passed two scalars together modulo the group order in @@ -780,159 +744,10 @@ func (s *ModNScalar) reduce512(t0, t1, t2, t3, t4, t5, t6, t7, t8, t9, t10, t11, // // The scalar is returned to support chaining. This enables syntax like: // s3.Mul2(s, s2).AddInt(1) so that s3 = (s * s2) + 1. -func (s *ModNScalar) Mul2(val, val2 *ModNScalar) *ModNScalar { - // This could be done with for loops and an array to store the intermediate - // terms, but this unrolled version is significantly faster. - - // The overall strategy employed here is: - // 1) Calculate the 512-bit product of the two scalars using the standard - // pencil-and-paper method. - // 2) Reduce the result modulo the prime by effectively subtracting - // multiples of the group order N (actually performed by adding multiples - // of the two's complement of N to avoid implementing subtraction). - // 3) Repeat step 2 noting that each iteration reduces the required number - // of bits by 127 because the two's complement of N has 127 leading zero - // bits. - // 4) Once reduced to 256 bits, call the existing reduce method to perform - // a final reduction as needed. - // - // Note that several of the intermediate calculations require adding 64-bit - // products together which would overflow a uint64, so a 96-bit accumulator - // is used instead. - - // Terms for 2^(32*0). - var acc accumulator96 - acc.Add(uint64(val.n[0]) * uint64(val2.n[0])) - t0 := uint64(acc.n[0]) - acc.Rsh32() - - // Terms for 2^(32*1). - acc.Add(uint64(val.n[0]) * uint64(val2.n[1])) - acc.Add(uint64(val.n[1]) * uint64(val2.n[0])) - t1 := uint64(acc.n[0]) - acc.Rsh32() - - // Terms for 2^(32*2). - acc.Add(uint64(val.n[0]) * uint64(val2.n[2])) - acc.Add(uint64(val.n[1]) * uint64(val2.n[1])) - acc.Add(uint64(val.n[2]) * uint64(val2.n[0])) - t2 := uint64(acc.n[0]) - acc.Rsh32() - - // Terms for 2^(32*3). - acc.Add(uint64(val.n[0]) * uint64(val2.n[3])) - acc.Add(uint64(val.n[1]) * uint64(val2.n[2])) - acc.Add(uint64(val.n[2]) * uint64(val2.n[1])) - acc.Add(uint64(val.n[3]) * uint64(val2.n[0])) - t3 := uint64(acc.n[0]) - acc.Rsh32() - - // Terms for 2^(32*4). - acc.Add(uint64(val.n[0]) * uint64(val2.n[4])) - acc.Add(uint64(val.n[1]) * uint64(val2.n[3])) - acc.Add(uint64(val.n[2]) * uint64(val2.n[2])) - acc.Add(uint64(val.n[3]) * uint64(val2.n[1])) - acc.Add(uint64(val.n[4]) * uint64(val2.n[0])) - t4 := uint64(acc.n[0]) - acc.Rsh32() - - // Terms for 2^(32*5). - acc.Add(uint64(val.n[0]) * uint64(val2.n[5])) - acc.Add(uint64(val.n[1]) * uint64(val2.n[4])) - acc.Add(uint64(val.n[2]) * uint64(val2.n[3])) - acc.Add(uint64(val.n[3]) * uint64(val2.n[2])) - acc.Add(uint64(val.n[4]) * uint64(val2.n[1])) - acc.Add(uint64(val.n[5]) * uint64(val2.n[0])) - t5 := uint64(acc.n[0]) - acc.Rsh32() - - // Terms for 2^(32*6). - acc.Add(uint64(val.n[0]) * uint64(val2.n[6])) - acc.Add(uint64(val.n[1]) * uint64(val2.n[5])) - acc.Add(uint64(val.n[2]) * uint64(val2.n[4])) - acc.Add(uint64(val.n[3]) * uint64(val2.n[3])) - acc.Add(uint64(val.n[4]) * uint64(val2.n[2])) - acc.Add(uint64(val.n[5]) * uint64(val2.n[1])) - acc.Add(uint64(val.n[6]) * uint64(val2.n[0])) - t6 := uint64(acc.n[0]) - acc.Rsh32() - - // Terms for 2^(32*7). - acc.Add(uint64(val.n[0]) * uint64(val2.n[7])) - acc.Add(uint64(val.n[1]) * uint64(val2.n[6])) - acc.Add(uint64(val.n[2]) * uint64(val2.n[5])) - acc.Add(uint64(val.n[3]) * uint64(val2.n[4])) - acc.Add(uint64(val.n[4]) * uint64(val2.n[3])) - acc.Add(uint64(val.n[5]) * uint64(val2.n[2])) - acc.Add(uint64(val.n[6]) * uint64(val2.n[1])) - acc.Add(uint64(val.n[7]) * uint64(val2.n[0])) - t7 := uint64(acc.n[0]) - acc.Rsh32() - - // Terms for 2^(32*8). - acc.Add(uint64(val.n[1]) * uint64(val2.n[7])) - acc.Add(uint64(val.n[2]) * uint64(val2.n[6])) - acc.Add(uint64(val.n[3]) * uint64(val2.n[5])) - acc.Add(uint64(val.n[4]) * uint64(val2.n[4])) - acc.Add(uint64(val.n[5]) * uint64(val2.n[3])) - acc.Add(uint64(val.n[6]) * uint64(val2.n[2])) - acc.Add(uint64(val.n[7]) * uint64(val2.n[1])) - t8 := uint64(acc.n[0]) - acc.Rsh32() - - // Terms for 2^(32*9). - acc.Add(uint64(val.n[2]) * uint64(val2.n[7])) - acc.Add(uint64(val.n[3]) * uint64(val2.n[6])) - acc.Add(uint64(val.n[4]) * uint64(val2.n[5])) - acc.Add(uint64(val.n[5]) * uint64(val2.n[4])) - acc.Add(uint64(val.n[6]) * uint64(val2.n[3])) - acc.Add(uint64(val.n[7]) * uint64(val2.n[2])) - t9 := uint64(acc.n[0]) - acc.Rsh32() - - // Terms for 2^(32*10). - acc.Add(uint64(val.n[3]) * uint64(val2.n[7])) - acc.Add(uint64(val.n[4]) * uint64(val2.n[6])) - acc.Add(uint64(val.n[5]) * uint64(val2.n[5])) - acc.Add(uint64(val.n[6]) * uint64(val2.n[4])) - acc.Add(uint64(val.n[7]) * uint64(val2.n[3])) - t10 := uint64(acc.n[0]) - acc.Rsh32() - - // Terms for 2^(32*11). - acc.Add(uint64(val.n[4]) * uint64(val2.n[7])) - acc.Add(uint64(val.n[5]) * uint64(val2.n[6])) - acc.Add(uint64(val.n[6]) * uint64(val2.n[5])) - acc.Add(uint64(val.n[7]) * uint64(val2.n[4])) - t11 := uint64(acc.n[0]) - acc.Rsh32() - - // Terms for 2^(32*12). - acc.Add(uint64(val.n[5]) * uint64(val2.n[7])) - acc.Add(uint64(val.n[6]) * uint64(val2.n[6])) - acc.Add(uint64(val.n[7]) * uint64(val2.n[5])) - t12 := uint64(acc.n[0]) - acc.Rsh32() - - // Terms for 2^(32*13). - acc.Add(uint64(val.n[6]) * uint64(val2.n[7])) - acc.Add(uint64(val.n[7]) * uint64(val2.n[6])) - t13 := uint64(acc.n[0]) - acc.Rsh32() - - // Terms for 2^(32*14). - acc.Add(uint64(val.n[7]) * uint64(val2.n[7])) - t14 := uint64(acc.n[0]) - acc.Rsh32() - - // What's left is for 2^(32*15). - t15 := uint64(acc.n[0]) - // acc.Rsh32() // No need since not used after this. Guaranteed to be 0. - - // At this point, all of the terms are grouped into their respective base - // and occupy up to 512 bits. Reduce the result accordingly. - s.reduce512(t0, t1, t2, t3, t4, t5, t6, t7, t8, t9, t10, t11, t12, t13, t14, - t15) +func (s *ModNScalar) Mul2(a, b *ModNScalar) *ModNScalar { + var product [8]uint64 + field64Mul512(&product, &a.n, &b.n) + scalar64Reduce512(&s.n, &product) return s } @@ -951,15 +766,10 @@ func (s *ModNScalar) Mul(val *ModNScalar) *ModNScalar { // The scalar is returned to support chaining. This enables syntax like: // s3.SquareVal(s).Mul(s) so that s3 = s^2 * s = s^3. func (s *ModNScalar) SquareVal(val *ModNScalar) *ModNScalar { - // This could technically be optimized slightly to take advantage of the - // fact that many of the intermediate calculations in squaring are just - // doubling, however, benchmarking has shown that due to the need to use a - // 96-bit accumulator, any savings are essentially offset by that and - // consequently there is no real difference in performance over just - // multiplying the value by itself to justify the extra code for now. This - // can be revisited in the future if it becomes a bottleneck in practice. - - return s.Mul2(val, val) + var product [8]uint64 + field64Square512(&product, &val.n) + scalar64Reduce512(&s.n, &product) + return s } // Square squares the scalar modulo the group order in constant time. The @@ -982,32 +792,38 @@ func (s *ModNScalar) NegateVal(val *ModNScalar) *ModNScalar { // minus the value. This implies that the result will always be in the // desired range with the sole exception of 0 because N - 0 = N itself. // - // Therefore, in order to avoid the need to reduce the result for every - // other case in order to achieve constant time, this creates a mask that is - // all 0s in the case of the scalar being negated is 0 and all 1s otherwise - // and bitwise ands that mask with each word. - // - // Finally, to simplify the carry propagation, this adds the two's - // complement of the scalar to N in order to achieve the same result. - bits := val.n[0] | val.n[1] | val.n[2] | val.n[3] | val.n[4] | val.n[5] | - val.n[6] | val.n[7] - mask := uint64(uint32Mask * constantTimeNotEq(bits, 0)) - c := uint64(orderWordZero) + (uint64(^val.n[0]) + 1) - s.n[0] = uint32(c & mask) - c = (c >> 32) + uint64(orderWordOne) + uint64(^val.n[1]) - s.n[1] = uint32(c & mask) - c = (c >> 32) + uint64(orderWordTwo) + uint64(^val.n[2]) - s.n[2] = uint32(c & mask) - c = (c >> 32) + uint64(orderWordThree) + uint64(^val.n[3]) - s.n[3] = uint32(c & mask) - c = (c >> 32) + uint64(orderWordFour) + uint64(^val.n[4]) - s.n[4] = uint32(c & mask) - c = (c >> 32) + uint64(orderWordFive) + uint64(^val.n[5]) - s.n[5] = uint32(c & mask) - c = (c >> 32) + uint64(orderWordSix) + uint64(^val.n[6]) - s.n[6] = uint32(c & mask) - c = (c >> 32) + uint64(orderWordSeven) + uint64(^val.n[7]) - s.n[7] = uint32(c & mask) + // The following handles that case in constant time by creating a mask that + // is all 0s in the case the scalar being negated is 0 and all 1s otherwise + // and then bitwise ands that mask with each limb. + // + // This approach was chosen over subtracting from 0 and then conditionally + // adding N because it's over twice as fast on 32-bit hardware while only + // being about 3-4% slower on 64-bit hardware. + // + // Determine mask first to allow aliasing. + mask := -uint64(constantTimeNotEq64(val.n[0]|val.n[1]|val.n[2]|val.n[3], 0)) + + // Unconditionally subtract the scalar from the group order. + // + // To simplify the carry propagation, this adds the two's complement of the + // scalar to N in order to achieve the same result. + // + // s = N - val + var c uint64 + s.n[0], c = bits.Add64(orderLimb0, ^val.n[0], 1) + s.n[1], c = bits.Add64(orderLimb1, ^val.n[1], c) + s.n[2], c = bits.Add64(orderLimb2, ^val.n[2], c) + s.n[3], _ = bits.Add64(orderLimb3, ^val.n[3], c) + + // Either keep the result when val != 0 or clear it when val == 0. The + // result is either: + // + // val == 0: s = 0 + // val != 0: s = N - val + s.n[0] &= mask + s.n[1] &= mask + s.n[2] &= mask + s.n[3] &= mask return s } @@ -1079,26 +895,42 @@ func (s *ModNScalar) InverseNonConst() *ModNScalar { // IsOverHalfOrder returns whether or not the scalar exceeds the group order // divided by 2 in constant time. func (s *ModNScalar) IsOverHalfOrder() bool { - // The intuition here is that the scalar is greater than half of the group - // order if one of the higher individual words is greater than the - // corresponding word of the half group order and all higher words in the - // scalar are equal to their corresponding word of the half group order. - // - // Note that the words 4, 5, and 6 are all the max uint32 value, so there is - // no need to test if those individual words of the scalar exceeds them, - // hence, only equality is checked for them. - result := constantTimeGreater(s.n[7], halfOrderWordSeven) - highWordsEqual := constantTimeEq(s.n[7], halfOrderWordSeven) - highWordsEqual &= constantTimeEq(s.n[6], halfOrderWordSix) - highWordsEqual &= constantTimeEq(s.n[5], halfOrderWordFive) - highWordsEqual &= constantTimeEq(s.n[4], halfOrderWordFour) - result |= highWordsEqual & constantTimeGreater(s.n[3], halfOrderWordThree) - highWordsEqual &= constantTimeEq(s.n[3], halfOrderWordThree) - result |= highWordsEqual & constantTimeGreater(s.n[2], halfOrderWordTwo) - highWordsEqual &= constantTimeEq(s.n[2], halfOrderWordTwo) - result |= highWordsEqual & constantTimeGreater(s.n[1], halfOrderWordOne) - highWordsEqual &= constantTimeEq(s.n[1], halfOrderWordOne) - result |= highWordsEqual & constantTimeGreater(s.n[0], halfOrderWordZero) - - return result != 0 + // These fields provide convenient access to each of the limbs of the + // secp256k1 curve group order N / 2 to improve code readability and avoid + // the need to recalculate them. + // + // The half order of the secp256k1 curve group is: + // 0x7fffffffffffffffffffffffffffffff5d576e7357a4501ddfe92f46681b20a0 + // + // Converting that to field representation (base 2^64) is: + // + // n[0] = 0xdfe92f46681b20a0 + // n[1] = 0x5d576e7357a4501d + // n[2] = 0xffffffffffffffff + // n[3] = 0x7fffffffffffffff + // + // This can be verified with the following test code: + // halfOrder := new(big.Int).Div(curveParams.N, big.NewInt(2)) + // var s ModNScalar + // s.SetByteSlice(halfOrder.Bytes()) + // t.Logf("%x", s.n) + const ( + halfOrderLimb0 uint64 = 0xdfe92f46681b20a0 + halfOrderLimb1 uint64 = 0x5d576e7357a4501d + halfOrderLimb2 uint64 = 0xffffffffffffffff + halfOrderLimb3 uint64 = 0x7fffffffffffffff + ) + + // The goal is to return true when the scalar is greater than half of the + // group order. That is, return true when s > N/2, which is trivially + // rearranged to s - (N/2 + 1) ≥ 0. + // + // In other words, the condition is met iff subtracting (N/2 + 1) from s is + // non-negative (aka there was no borrow). + var borrow uint64 + _, borrow = bits.Sub64(s.n[0], halfOrderLimb0+1, borrow) + _, borrow = bits.Sub64(s.n[1], halfOrderLimb1, borrow) + _, borrow = bits.Sub64(s.n[2], halfOrderLimb2, borrow) + _, borrow = bits.Sub64(s.n[3], halfOrderLimb3, borrow) + return borrow == 0 } diff --git a/dcrec/secp256k1/modnscalar_test.go b/dcrec/secp256k1/modnscalar_test.go index 448378e0cf..94e329f237 100644 --- a/dcrec/secp256k1/modnscalar_test.go +++ b/dcrec/secp256k1/modnscalar_test.go @@ -1,4 +1,4 @@ -// Copyright (c) 2020-2023 The Decred developers +// Copyright (c) 2020-2026 The Decred developers // Use of this source code is governed by an ISC // license that can be found in the LICENSE file. @@ -116,23 +116,27 @@ func TestModNScalarSetInt(t *testing.T) { tests := []struct { name string // test description in uint32 // test value - expected [8]uint32 // expected raw ints + expected [4]uint64 // expected limbs }{{ + name: "zero", + in: 0, + expected: [4]uint64{0, 0, 0, 0}, + }, { name: "five", in: 5, - expected: [8]uint32{5, 0, 0, 0, 0, 0, 0, 0}, + expected: [4]uint64{5, 0, 0, 0}, }, { - name: "group order word zero", - in: orderWordZero, - expected: [8]uint32{orderWordZero, 0, 0, 0, 0, 0, 0, 0}, + name: "max uint8 (2^8 - 1)", + in: 255, + expected: [4]uint64{255, 0, 0, 0}, }, { - name: "group order word zero + 1", - in: orderWordZero + 1, - expected: [8]uint32{orderWordZero + 1, 0, 0, 0, 0, 0, 0, 0}, + name: "max uint16 (2^16 - 1)", + in: 65535, + expected: [4]uint64{65535, 0, 0, 0}, }, { - name: "2^32 - 1", + name: "max uint32 (2^32 - 1)", in: 4294967295, - expected: [8]uint32{4294967295, 0, 0, 0, 0, 0, 0, 0}, + expected: [4]uint64{4294967295, 0, 0, 0}, }} for _, test := range tests { @@ -152,88 +156,136 @@ func TestModNScalarSetBytes(t *testing.T) { tests := []struct { name string // test description in string // hex encoded test value - expected [8]uint32 // expected raw ints + expected [4]uint64 // expected raw ints overflow bool // expected overflow result }{{ name: "zero", in: "00", - expected: [8]uint32{0, 0, 0, 0, 0, 0, 0, 0}, + expected: [4]uint64{0, 0, 0, 0}, + overflow: false, + }, { + name: "one (only limb zero set)", + in: "01", + expected: [4]uint64{1, 0, 0, 0}, + overflow: false, + }, { + name: "2^64 - 1 (limb zero all ones)", + in: "ffffffffffffffff", + expected: [4]uint64{0xffffffffffffffff, 0, 0, 0}, + overflow: false, + }, { + name: "2^64 (limb one low bit)", + in: "010000000000000000", + expected: [4]uint64{0, 1, 0, 0}, + overflow: false, + }, { + name: "2^128 - 1 (limbs zero and one all ones)", + in: "ffffffffffffffffffffffffffffffff", + expected: [4]uint64{0xffffffffffffffff, 0xffffffffffffffff, 0, 0}, + overflow: false, + }, { + name: "2^128 (limb two low bit)", + in: "0100000000000000000000000000000000", + expected: [4]uint64{0, 0, 1, 0}, + overflow: false, + }, { + name: "2^192 - 1 (limbs zero, one, and two all ones)", + in: "ffffffffffffffffffffffffffffffffffffffffffffffff", + expected: [4]uint64{ + 0xffffffffffffffff, 0xffffffffffffffff, 0xffffffffffffffff, 0, + }, + overflow: false, + }, { + name: "2^192 (limb three low bit)", + in: "01000000000000000000000000000000000000000000000000", + expected: [4]uint64{0, 0, 0, 1}, + overflow: false, + }, { + name: "2^255 (limb three high bit)", + in: "8000000000000000000000000000000000000000000000000000000000000000", + expected: [4]uint64{0, 0, 0, 0x8000000000000000}, + overflow: false, + }, { + name: "distinct value in every limb (verifies limb ordering)", + in: "0123456789abcdeffedcba98765432100f1e2d3c4b5a69788796a5b4c3d2e1f0", + expected: [4]uint64{ + 0x08796a5b4c3d2e1f0, 0x0f1e2d3c4b5a6978, 0xfedcba9876543210, + 0x0123456789abcdef, + }, overflow: false, }, { name: "group order (aka 0)", in: "fffffffffffffffffffffffffffffffebaaedce6af48a03bbfd25e8cd0364141", - expected: [8]uint32{0, 0, 0, 0, 0, 0, 0, 0}, + expected: [4]uint64{0, 0, 0, 0}, overflow: true, }, { - name: "group order - 1", + name: "group order - 1 (limb zero one below prime, upper limbs maxed)", in: "fffffffffffffffffffffffffffffffebaaedce6af48a03bbfd25e8cd0364140", - expected: [8]uint32{ - 0xd0364140, 0xbfd25e8c, 0xaf48a03b, 0xbaaedce6, - 0xfffffffe, 0xffffffff, 0xffffffff, 0xffffffff, + expected: [4]uint64{ + 0xbfd25e8cd0364140, 0xbaaedce6af48a03b, 0xfffffffffffffffe, + 0xffffffffffffffff, }, overflow: false, }, { - name: "group order + 1 (aka 1, overflow in word zero)", + name: "group order + 1 (aka 1, overflow in limb zero)", in: "fffffffffffffffffffffffffffffffebaaedce6af48a03bbfd25e8cd0364142", - expected: [8]uint32{1, 0, 0, 0, 0, 0, 0, 0}, + expected: [4]uint64{1, 0, 0, 0}, overflow: true, }, { - name: "group order word zero", - in: "d0364141", - expected: [8]uint32{0xd0364141, 0, 0, 0, 0, 0, 0, 0}, - overflow: false, - }, { - name: "group order word zero and one", + name: "group order limb zero", in: "bfd25e8cd0364141", - expected: [8]uint32{0xd0364141, 0xbfd25e8c, 0, 0, 0, 0, 0, 0}, + expected: [4]uint64{0xbfd25e8cd0364141, 0, 0, 0}, overflow: false, }, { - name: "group order words zero, one, and two", - in: "af48a03bbfd25e8cd0364141", - expected: [8]uint32{0xd0364141, 0xbfd25e8c, 0xaf48a03b, 0, 0, 0, 0, 0}, + name: "group order limbs zero and one", + in: "baaedce6af48a03bbfd25e8cd0364141", + expected: [4]uint64{0xbfd25e8cd0364141, 0xbaaedce6af48a03b, 0, 0}, overflow: false, }, { - name: "overflow in word one", - in: "fffffffffffffffffffffffffffffffebaaedce6af48a03bbfd25e8dd0364141", - expected: [8]uint32{0, 1, 0, 0, 0, 0, 0, 0}, - overflow: true, + name: "group order limbs zero, one, and two", + in: "fffffffffffffffebaaedce6af48a03bbfd25e8cd0364141", + expected: [4]uint64{ + 0xbfd25e8cd0364141, 0xbaaedce6af48a03b, 0xfffffffffffffffe, 0, + }, + overflow: false, }, { - name: "overflow in word two", + name: "group order + 2^64 (overflow in limb one)", in: "fffffffffffffffffffffffffffffffebaaedce6af48a03cbfd25e8cd0364141", - expected: [8]uint32{0, 0, 1, 0, 0, 0, 0, 0}, + expected: [4]uint64{0, 1, 0, 0}, overflow: true, }, { - name: "overflow in word three", - in: "fffffffffffffffffffffffffffffffebaaedce7af48a03bbfd25e8cd0364141", - expected: [8]uint32{0, 0, 0, 1, 0, 0, 0, 0}, + name: "group order + 2^128 (overflow in limb two)", + in: "ffffffffffffffffffffffffffffffffbaaedce6af48a03bbfd25e8cd0364141", + expected: [4]uint64{0, 0, 1, 0}, overflow: true, }, { - name: "overflow in word four", - in: "ffffffffffffffffffffffffffffffffbaaedce6af48a03bbfd25e8cd0364141", - expected: [8]uint32{0, 0, 0, 0, 1, 0, 0, 0}, + name: "2^256 - 1 (max input)", + in: "ffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffff", + expected: [4]uint64{0x402da1732fc9bebe, 0x4551231950b75fc4, 1, 0}, + overflow: true, }, { name: "(group order - 1) * 2 NOT mod N, truncated >32 bytes", in: "01fffffffffffffffffffffffffffffffd755db9cd5e9140777fa4bd19a06c8284", - expected: [8]uint32{ - 0x19a06c82, 0x777fa4bd, 0xcd5e9140, 0xfd755db9, - 0xffffffff, 0xffffffff, 0xffffffff, 0x01ffffff, + expected: [4]uint64{ + 0x777fa4bd19a06c82, 0xfd755db9cd5e9140, 0xffffffffffffffff, + 0x01ffffffffffffff, }, overflow: false, }, { - name: "alternating bits", + name: "alternating bits (spans all limbs)", in: "a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5", - expected: [8]uint32{ - 0xa5a5a5a5, 0xa5a5a5a5, 0xa5a5a5a5, 0xa5a5a5a5, - 0xa5a5a5a5, 0xa5a5a5a5, 0xa5a5a5a5, 0xa5a5a5a5, + expected: [4]uint64{ + 0xa5a5a5a5a5a5a5a5, 0xa5a5a5a5a5a5a5a5, 0xa5a5a5a5a5a5a5a5, + 0xa5a5a5a5a5a5a5a5, }, overflow: false, }, { - name: "alternating bits 2", + name: "alternating bits 2 (spans all limbs)", in: "5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a", - expected: [8]uint32{ - 0x5a5a5a5a, 0x5a5a5a5a, 0x5a5a5a5a, 0x5a5a5a5a, - 0x5a5a5a5a, 0x5a5a5a5a, 0x5a5a5a5a, 0x5a5a5a5a, + expected: [4]uint64{ + 0x5a5a5a5a5a5a5a5a, 0x5a5a5a5a5a5a5a5a, 0x5a5a5a5a5a5a5a5a, + 0x5a5a5a5a5a5a5a5a, }, overflow: false, }} @@ -305,35 +357,35 @@ func TestModNScalarBytes(t *testing.T) { in: "fffffffffffffffffffffffffffffffebaaedce6af48a03bbfd25e8cd0364140", expected: "fffffffffffffffffffffffffffffffebaaedce6af48a03bbfd25e8cd0364140", }, { - name: "group order + 1 (aka 1, overflow in word zero)", + name: "group order + 1 (aka 1, overflow in limb zero)", in: "fffffffffffffffffffffffffffffffebaaedce6af48a03bbfd25e8cd0364142", expected: "0000000000000000000000000000000000000000000000000000000000000001", }, { - name: "group order word zero", + name: "group order limb zero", in: "d0364141", expected: "00000000000000000000000000000000000000000000000000000000d0364141", }, { - name: "group order word zero and one", + name: "group order limb zero and one", in: "bfd25e8cd0364141", expected: "000000000000000000000000000000000000000000000000bfd25e8cd0364141", }, { - name: "group order words zero, one, and two", + name: "group order limbs zero, one, and two", in: "af48a03bbfd25e8cd0364141", expected: "0000000000000000000000000000000000000000af48a03bbfd25e8cd0364141", }, { - name: "overflow in word one", + name: "overflow in limb one", in: "fffffffffffffffffffffffffffffffebaaedce6af48a03bbfd25e8dd0364141", expected: "0000000000000000000000000000000000000000000000000000000100000000", }, { - name: "overflow in word two", + name: "overflow in limb two", in: "fffffffffffffffffffffffffffffffebaaedce6af48a03cbfd25e8cd0364141", expected: "0000000000000000000000000000000000000000000000010000000000000000", }, { - name: "overflow in word three", + name: "overflow in limb three", in: "fffffffffffffffffffffffffffffffebaaedce7af48a03bbfd25e8cd0364141", expected: "0000000000000000000000000000000000000001000000000000000000000000", }, { - name: "overflow in word four", + name: "overflow in limb four", in: "ffffffffffffffffffffffffffffffffbaaedce6af48a03bbfd25e8cd0364141", expected: "0000000000000000000000000000000100000000000000000000000000000000", }, { @@ -501,11 +553,11 @@ func TestModNScalarEqualsRandom(t *testing.T) { t.Fatalf("failed equality check\nscalar in: %v", s) } - // Flip a random bit in a random word and ensure it's no longer equal. - randomWord := rng.Int31n(int32(len(s.n))) - randomBit := uint32(1 << uint32(rng.Int31n(32))) + // Flip a random bit in a random limb and ensure it's no longer equal. + randomLimb := rng.Int31n(int32(len(s.n))) + randomBit := uint64(1 << uint64(rng.Int63n(64))) s2 := new(ModNScalar).Set(s) - s2.n[randomWord] ^= randomBit + s2.n[randomLimb] ^= randomBit if s2.Equals(s) { t.Fatalf("failed inequality check\nscalar in: %v", s2) } @@ -531,55 +583,35 @@ func TestModNScalarAdd(t *testing.T) { in2: "0", expected: "1", }, { - name: "group order (aka 0) + 1 (gets reduced, no overflow)", - in1: "fffffffffffffffffffffffffffffffebaaedce6af48a03bbfd25e8cd0364141", - in2: "1", - expected: "1", - }, { - name: "group order - 1 + 1 (aka 0, overflow to prime)", + name: "group order-1 + 1 (aka 0, overflow to order)", in1: "fffffffffffffffffffffffffffffffebaaedce6af48a03bbfd25e8cd0364140", in2: "1", expected: "0", }, { - name: "group order - 1 + 2 (aka 1, overflow in word zero)", - in1: "fffffffffffffffffffffffffffffffebaaedce6af48a03bbfd25e8cd0364140", - in2: "2", + name: "group order (aka 0) + 1 (aka 1, gets reduced, no overflow)", + in1: "fffffffffffffffffffffffffffffffebaaedce6af48a03bbfd25e8cd0364141", + in2: "1", expected: "1", }, { - name: "overflow in word one", - in1: "fffffffffffffffffffffffffffffffebaaedce6af48a03bbfd25e8bd0364141", - in2: "100000001", + name: "group order-1 + 2 (aka 1, overflow in limb zero)", + in1: "fffffffffffffffffffffffffffffffebaaedce6af48a03bbfd25e8cd0364140", + in2: "2", expected: "1", }, { - name: "overflow in word two", + name: "overflow in limb one", in1: "fffffffffffffffffffffffffffffffebaaedce6af48a03abfd25e8cd0364141", in2: "10000000000000001", expected: "1", }, { - name: "overflow in word three", - in1: "fffffffffffffffffffffffffffffffebaaedce5af48a03bbfd25e8cd0364141", - in2: "1000000000000000000000001", - expected: "1", - }, { - name: "overflow in word four", + name: "overflow in limb two", in1: "fffffffffffffffffffffffffffffffdbaaedce6af48a03bbfd25e8cd0364141", in2: "100000000000000000000000000000001", expected: "1", }, { - name: "overflow in word five", - in1: "fffffffffffffffffffffffefffffffebaaedce6af48a03bbfd25e8cd0364141", - in2: "10000000000000000000000000000000000000001", - expected: "1", - }, { - name: "overflow in word six", + name: "overflow in limb three", in1: "fffffffffffffffefffffffffffffffebaaedce6af48a03bbfd25e8cd0364141", in2: "1000000000000000000000000000000000000000000000001", expected: "1", - }, { - name: "overflow in word seven", - in1: "fffffffefffffffffffffffffffffffebaaedce6af48a03bbfd25e8cd0364141", - in2: "100000000000000000000000000000000000000000000000000000001", - expected: "1", }, { name: "alternating bits", in1: "a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5", @@ -601,7 +633,7 @@ func TestModNScalarAdd(t *testing.T) { // Ensure the result has the expected value. s1.Add(s2) if !s1.Equals(expected) { - t.Errorf("%s: unexpected result\ngot: %x\nwant: %x", test.name, + t.Errorf("%s: unexpected result\ngot: %v\nwant: %v", test.name, s1, expected) continue } @@ -645,53 +677,6 @@ func TestModNScalarAddRandom(t *testing.T) { } } -// TestAccumulator96Add ensures that the internal 96-bit accumulator used by -// multiplication works as expected for overflow edge cases including overflow. -func TestAccumulator96Add(t *testing.T) { - tests := []struct { - name string // test description - start accumulator96 // starting value of accumulator - in uint64 // value to add to accumulator - expected accumulator96 // expected value of accumulator after addition - }{{ - name: "0 + 0 = 0", - start: accumulator96{[3]uint32{0, 0, 0}}, - in: 0, - expected: accumulator96{[3]uint32{0, 0, 0}}, - }, { - name: "overflow in word zero", - start: accumulator96{[3]uint32{0xffffffff, 0, 0}}, - in: 1, - expected: accumulator96{[3]uint32{0, 1, 0}}, - }, { - name: "overflow in word one", - start: accumulator96{[3]uint32{0, 0xffffffff, 0}}, - in: 0x100000000, - expected: accumulator96{[3]uint32{0, 0, 1}}, - }, { - name: "overflow in words one and two", - start: accumulator96{[3]uint32{0xffffffff, 0xffffffff, 0}}, - in: 1, - expected: accumulator96{[3]uint32{0, 0, 1}}, - }, { - // Start accumulator at 129127208455837319175 which is the result of - // 4294967295 * 4294967295 accumulated seven times. - name: "max result from eight adds of max uint32 multiplications", - start: accumulator96{[3]uint32{7, 4294967282, 6}}, - in: 18446744065119617025, - expected: accumulator96{[3]uint32{8, 4294967280, 7}}, - }} - - for _, test := range tests { - acc := test.start - acc.Add(test.in) - if acc.n != test.expected.n { - t.Errorf("%s: wrong result\ngot: %v\nwant: %v", test.name, acc.n, - test.expected.n) - } - } -} - // TestModNScalarMul ensures that multiplying two scalars together works as // expected for edge cases. func TestModNScalarMul(t *testing.T) { @@ -720,6 +705,16 @@ func TestModNScalarMul(t *testing.T) { in1: "fffffffffffffffffffffffffffffffebaaedce6af48a03bbfd25e8cd0364140", in2: "2", expected: "fffffffffffffffffffffffffffffffebaaedce6af48a03bbfd25e8cd036413f", + }, { + name: "(group order-1) * 3", + in1: "fffffffffffffffffffffffffffffffebaaedce6af48a03bbfd25e8cd0364140", + in2: "3", + expected: "fffffffffffffffffffffffffffffffebaaedce6af48a03bbfd25e8cd036413e", + }, { + name: "(group order-1) * 4", + in1: "fffffffffffffffffffffffffffffffebaaedce6af48a03bbfd25e8cd0364140", + in2: "4", + expected: "fffffffffffffffffffffffffffffffebaaedce6af48a03bbfd25e8cd036413d", }, { name: "(group order-1) * (group order-1)", in1: "fffffffffffffffffffffffffffffffebaaedce6af48a03bbfd25e8cd0364140", @@ -730,56 +725,76 @@ func TestModNScalarMul(t *testing.T) { in1: "7fffffffffffffffffffffffffffffff5d576e7357a4501ddfe92f46681b20a1", in2: "2", expected: "1", + }, { + name: "highest bit in limb zero squared", + in1: "8000000000000000", + in2: "8000000000000000", + expected: "40000000000000000000000000000000", + }, { + name: "highest bit in limb one squared", + in1: "80000000000000000000000000000000", + in2: "80000000000000000000000000000000", + expected: "4000000000000000000000000000000000000000000000000000000000000000", + }, { + name: "highest bit in limb two squared", + in1: "800000000000000000000000000000000000000000000000", + in2: "800000000000000000000000000000000000000000000000", + expected: "515448c6542dd7f1100b685ccbf26fafc0000000000000000000000000000000", + }, { + name: "highest bit in limb three squared", + in1: "8000000000000000000000000000000000000000000000000000000000000000", + in2: "8000000000000000000000000000000000000000000000000000000000000000", + expected: "2759c7356071a6f179a5fd7916f341f19d0525b0839f3e1e225b3c8519f5f450", + }, { + name: "cross limb product 64x128", + in1: "10000000000000000", + in2: "100000000000000000000000000000000", + expected: "1000000000000000000000000000000000000000000000000", + }, { + name: "cross limb product 64x192", + in1: "10000000000000000", + in2: "1000000000000000000000000000000000000000000000000", + expected: "14551231950b75fc4402da1732fc9bebf", + }, { + name: "cross limb product 128x192", + in1: "100000000000000000000000000000000", + in2: "1000000000000000000000000000000000000000000000000", + expected: "14551231950b75fc4402da1732fc9bebf0000000000000000", }, { name: "group order (aka 0) * 3", in1: "fffffffffffffffffffffffffffffffebaaedce6af48a03bbfd25e8cd0364141", in2: "3", expected: "0", }, { - name: "overflow in word eight", + name: "overflow in limb four", in1: "100000000000000000000000000000000", in2: "100000000000000000000000000000000", expected: "14551231950b75fc4402da1732fc9bebf", }, { - name: "overflow in word nine", - in1: "1000000000000000000000000000000000000", - in2: "1000000000000000000000000000000000000", - expected: "14551231950b75fc4402da1732fc9bebf00000000", - }, { - name: "overflow in word ten", + name: "overflow in limb five", in1: "10000000000000000000000000000000000000000", in2: "10000000000000000000000000000000000000000", expected: "14551231950b75fc4402da1732fc9bebf0000000000000000", }, { - name: "overflow in word eleven", - in1: "100000000000000000000000000000000000000000000", - in2: "100000000000000000000000000000000000000000000", - expected: "14551231950b75fc4402da1732fc9bebf000000000000000000000000", - }, { - name: "overflow in word twelve", + name: "overflow in limb six", in1: "1000000000000000000000000000000000000000000000000", in2: "1000000000000000000000000000000000000000000000000", expected: "4551231950b75fc4402da1732fc9bec04551231950b75fc4402da1732fc9bebf", }, { - name: "overflow in word thirteen", - in1: "10000000000000000000000000000000000000000000000000000", - in2: "10000000000000000000000000000000000000000000000000000", - expected: "50b75fc4402da1732fc9bec09d671cd51b343a1b66926b57d2a4c1c61536bda7", - }, { - name: "overflow in word fourteen", + name: "overflow in limb seven", in1: "100000000000000000000000000000000000000000000000000000000", in2: "100000000000000000000000000000000000000000000000000000000", expected: "402da1732fc9bec09d671cd581c69bc59509b0b074ec0aea8f564d667ec7eb3c", }, { - name: "overflow in word fifteen", - in1: "1000000000000000000000000000000000000000000000000000000000000", - in2: "1000000000000000000000000000000000000000000000000000000000000", - expected: "2fc9bec09d671cd581c69bc5e697f5e41f12c33a0a7b6f4e3302b92ea029cecd", + name: "max limb * one", + in1: "ffffffffffffffff", + in2: "1", + expected: "ffffffffffffffff", }, { - name: "double overflow in internal accumulator", - in1: "fffffffffffffffffffffffffffffffebaaedce6af48a03bbfd25e8cd0364140", - in2: "55555555555555555555555555555554e8e4f44ce51835693ff0ca2ef01215c2", - expected: "aaaaaaaaaaaaaaaaaaaaaaaaaaaaaaa9d1c9e899ca306ad27fe1945de0242b7f", + name: "max limb * max limb", + in1: "ffffffffffffffff", + in2: "ffffffffffffffff", + expected: "fffffffffffffffe0000000000000001", }, { name: "alternating bits", in1: "a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5", @@ -881,35 +896,35 @@ func TestModNScalarSquare(t *testing.T) { in: "fffffffffffffffffffffffffffffffebaaedce6af48a03bbfd25e8cd0364142", expected: "1", }, { - name: "overflow in word eight", + name: "overflow in limb eight", in: "100000000000000000000000000000000", expected: "14551231950b75fc4402da1732fc9bebf", }, { - name: "overflow in word nine", + name: "overflow in limb nine", in: "1000000000000000000000000000000000000", expected: "14551231950b75fc4402da1732fc9bebf00000000", }, { - name: "overflow in word ten", + name: "overflow in limb ten", in: "10000000000000000000000000000000000000000", expected: "14551231950b75fc4402da1732fc9bebf0000000000000000", }, { - name: "overflow in word eleven", + name: "overflow in limb eleven", in: "100000000000000000000000000000000000000000000", expected: "14551231950b75fc4402da1732fc9bebf000000000000000000000000", }, { - name: "overflow in word twelve", + name: "overflow in limb twelve", in: "1000000000000000000000000000000000000000000000000", expected: "4551231950b75fc4402da1732fc9bec04551231950b75fc4402da1732fc9bebf", }, { - name: "overflow in word thirteen", + name: "overflow in limb thirteen", in: "10000000000000000000000000000000000000000000000000000", expected: "50b75fc4402da1732fc9bec09d671cd51b343a1b66926b57d2a4c1c61536bda7", }, { - name: "overflow in word fourteen", + name: "overflow in limb fourteen", in: "100000000000000000000000000000000000000000000000000000000", expected: "402da1732fc9bec09d671cd581c69bc59509b0b074ec0aea8f564d667ec7eb3c", }, { - name: "overflow in word fifteen", + name: "overflow in limb fifteen", in: "1000000000000000000000000000000000000000000000000000000000000", expected: "2fc9bec09d671cd581c69bc5e697f5e41f12c33a0a7b6f4e3302b92ea029cecd", }, { @@ -1015,33 +1030,41 @@ func TestModNScalarNegate(t *testing.T) { in: "1", expected: "fffffffffffffffffffffffffffffffebaaedce6af48a03bbfd25e8cd0364140", }, { - name: "negation in word one", - in: "0000000000000000000000000000000000000000000000000000000100000000", - expected: "fffffffffffffffffffffffffffffffebaaedce6af48a03bbfd25e8bd0364141", + name: "two", + in: "2", + expected: "fffffffffffffffffffffffffffffffebaaedce6af48a03bbfd25e8cd036413f", + }, { + name: "three", + in: "3", + expected: "fffffffffffffffffffffffffffffffebaaedce6af48a03bbfd25e8cd036413e", }, { - name: "negation in word two", + name: "group order - 2", + in: "fffffffffffffffffffffffffffffffebaaedce6af48a03bbfd25e8cd036413f", + expected: "2", + }, { + name: "group order - 1", + in: "fffffffffffffffffffffffffffffffebaaedce6af48a03bbfd25e8cd0364140", + expected: "1", + }, { + name: "negation in limb one", in: "0000000000000000000000000000000000000000000000010000000000000000", expected: "fffffffffffffffffffffffffffffffebaaedce6af48a03abfd25e8cd0364141", }, { - name: "negation in word three", - in: "0000000000000000000000000000000000000001000000000000000000000000", - expected: "fffffffffffffffffffffffffffffffebaaedce5af48a03bbfd25e8cd0364141", - }, { - name: "negation in word four", + name: "negation in limb two", in: "0000000000000000000000000000000100000000000000000000000000000000", expected: "fffffffffffffffffffffffffffffffdbaaedce6af48a03bbfd25e8cd0364141", }, { - name: "negation in word five", - in: "0000000000000000000000010000000000000000000000000000000000000000", - expected: "fffffffffffffffffffffffefffffffebaaedce6af48a03bbfd25e8cd0364141", - }, { - name: "negation in word six", + name: "negation in limb three", in: "0000000000000001000000000000000000000000000000000000000000000000", expected: "fffffffffffffffefffffffffffffffebaaedce6af48a03bbfd25e8cd0364141", }, { - name: "negation in word seven", - in: "0000000100000000000000000000000000000000000000000000000000000000", - expected: "fffffffefffffffffffffffffffffffebaaedce6af48a03bbfd25e8cd0364141", + name: "half group order (floor)", + in: "7fffffffffffffffffffffffffffffff5d576e7357a4501ddfe92f46681b20a0", + expected: "7fffffffffffffffffffffffffffffff5d576e7357a4501ddfe92f46681b20a1", + }, { + name: "half group order (ceil)", + in: "7fffffffffffffffffffffffffffffff5d576e7357a4501ddfe92f46681b20a1", + expected: "7fffffffffffffffffffffffffffffff5d576e7357a4501ddfe92f46681b20a0", }, { name: "alternating bits", in: "a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5", @@ -1157,33 +1180,29 @@ func TestModNScalarInverseNonConst(t *testing.T) { in: "1", expected: "1", }, { - name: "inverse carry in word one", - in: "0000000000000000000000000000000000000000000000000000000100000000", - expected: "5588b13effffffffffffffffffffffff934e5b00ca8417bf50177f7ba415411a", + name: "two", + in: "2", + expected: "7fffffffffffffffffffffffffffffff5d576e7357a4501ddfe92f46681b20a1", + }, { + name: "group order - 1", + in: "fffffffffffffffffffffffffffffffebaaedce6af48a03bbfd25e8cd0364140", + expected: "fffffffffffffffffffffffffffffffebaaedce6af48a03bbfd25e8cd0364140", }, { - name: "inverse carry in word two", + name: "group order - 2", + in: "fffffffffffffffffffffffffffffffebaaedce6af48a03bbfd25e8cd036413f", + expected: "7fffffffffffffffffffffffffffffff5d576e7357a4501ddfe92f46681b20a0", + }, { + name: "inverse carry in limb one", in: "0000000000000000000000000000000000000000000000010000000000000000", expected: "4b0dff665588b13effffffffffffffffa09f710af01555259d4ad302583de6dc", }, { - name: "inverse carry in word three", - in: "0000000000000000000000000000000000000001000000000000000000000000", - expected: "34b9ec244b0dff665588b13effffffffbcff4127932a971a78274c9d74176b38", - }, { - name: "inverse carry in word four", + name: "inverse carry in limb two", in: "0000000000000000000000000000000100000000000000000000000000000000", expected: "50a51ac834b9ec244b0dff665588b13e9984d5b3cf80ef0fd6a23766a3ee9f22", }, { - name: "inverse carry in word five", - in: "0000000000000000000000010000000000000000000000000000000000000000", - expected: "27cfab5e50a51ac834b9ec244b0dff6622f16e85b683d5a059bcd5a3b29d9dff", - }, { - name: "inverse carry in word six", + name: "inverse carry in limb three", in: "0000000000000001000000000000000000000000000000000000000000000000", expected: "897f30c127cfab5e50a51ac834b9ec239c53f268b4700c14f19b9499ac58d8ad", - }, { - name: "inverse carry in word seven", - in: "0000000100000000000000000000000000000000000000000000000000000000", - expected: "6494ef93897f30c127cfab5e50a51ac7b4e8f713e0cddd182234e907286ae6b3", }, { name: "alternating bits", in: "a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5", @@ -1308,21 +1327,33 @@ func TestModNScalarIsOverHalfOrder(t *testing.T) { in: "7fffffffffffffffffffffffffffffff5d576e7357a4501ddfe92f46681b20a1", expected: true, }, { - name: "over half order word one", - in: "7fffffffffffffffffffffffffffffff5d576e7357a4501ddfe92f47681b20a0", + name: "over half order limb one", + in: "7fffffffffffffffffffffffffffffff5d576e7357a4501edfe92f46681b20a0", expected: true, }, { - name: "over half order word two", - in: "7fffffffffffffffffffffffffffffff5d576e7357a4501edfe92f46681b20a0", + name: "group order - 2", + in: "fffffffffffffffffffffffffffffffebaaedce6af48a03bbfd25e8cd036413f", + expected: true, + }, { + name: "group order - 1", + in: "fffffffffffffffffffffffffffffffebaaedce6af48a03bbfd25e8cd0364140", expected: true, }, { - name: "over half order word three", - in: "7fffffffffffffffffffffffffffffff5d576e7457a4501ddfe92f46681b20a0", + name: "highest limb below half order", + in: "7fffffffffffffffffffffffffffffff00000000000000000000000000000000", + expected: false, + }, { + name: "highest bit set", + in: "8000000000000000000000000000000000000000000000000000000000000000", expected: true, }, { - name: "over half order word seven", - in: "8fffffffffffffffffffffffffffffff5d576e7357a4501ddfe92f46681b20a0", + name: "alternating bits", + in: "a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5", expected: true, + }, { + name: "alternating bits 2", + in: "5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a", + expected: false, }} for _, test := range tests { From 688f0f38a0ea3b2560b43eead68f8cf1a88fcf5b Mon Sep 17 00:00:00 2001 From: Dave Collins Date: Thu, 9 Jul 2026 01:28:34 -0500 Subject: [PATCH 06/17] secp256k1: Add scalar reduction proofs. This adds a formal verification proof of the scalar reduction arithmetic in the new 4x64 scalar implementation and updates the proofs README.md accordingly. It includes a formal proof for the following: - 512-bit modular reduction over the group order with saturated 64-bit limbs --- dcrec/secp256k1/internal/proofs/README.md | 9 +- .../proofs/modnscalar_4x64_reduce_prover.py | 571 ++++++++++++++++++ dcrec/secp256k1/modnscalar.go | 4 + 3 files changed, 582 insertions(+), 2 deletions(-) create mode 100644 dcrec/secp256k1/internal/proofs/modnscalar_4x64_reduce_prover.py diff --git a/dcrec/secp256k1/internal/proofs/README.md b/dcrec/secp256k1/internal/proofs/README.md index ba7e919f6e..079069a8f9 100644 --- a/dcrec/secp256k1/internal/proofs/README.md +++ b/dcrec/secp256k1/internal/proofs/README.md @@ -6,8 +6,9 @@ proofs This consists of formal verification artifacts used to help establish the correctness of the security-critical arithmetic implementations. -The current proofs formally verify properties of the optimized field arithmetic -implementations using [z3](https://github.com/Z3Prover/z3): +The current proofs formally verify properties of the optimized scalar and field +arithmetic implementations using the [Z3 Theorem +Prover](https://github.com/Z3Prover/z3): - Correctness of multi-precision multiplication - Bounds on intermediate values @@ -45,3 +46,7 @@ optimized implementations preserve the required arithmetic invariants. - [Field Modular Reduction](field_4x64_reduce_prover.py) Provides formal verification of the reduction of a 512-bit value represented using saturated 64-bit limbs modulo the secp256k1 prime. + +- [Scalar Modular Reduction](modnscalar_4x64_reduce_prover.py) + Provides formal verification of the reduction of a 512-bit value represented + using saturated 64-bit limbs modulo the secp256k1 group order. diff --git a/dcrec/secp256k1/internal/proofs/modnscalar_4x64_reduce_prover.py b/dcrec/secp256k1/internal/proofs/modnscalar_4x64_reduce_prover.py new file mode 100644 index 0000000000..1414738768 --- /dev/null +++ b/dcrec/secp256k1/internal/proofs/modnscalar_4x64_reduce_prover.py @@ -0,0 +1,571 @@ +# Copyright (c) 2026 The Decred developers +# Copyright (c) 2026 Dave Collins +# Use of this source code is governed by an ISC +# license that can be found in the LICENSE file. + +from z3_proof_helpers import * + +# ------- +# Inputs. +# ------- + +x = BitVecs('x0 x1 x2 x3 x4 x5 x6 x7', 64) + +# NOTE: The code in this section reuses the variable names from the Go code +# being proven for easy cross referencing. However, the comments refer to them +# by shorter mathematical variables for clarity. + +# Define N. +orderLimb0 = BitVecVal(0xbfd25e8cd0364141, 64) +orderLimb1 = BitVecVal(0xbaaedce6af48a03b, 64) +orderLimb2 = BitVecVal(0xfffffffffffffffe, 64) +orderLimb3 = BitVecVal(0xffffffffffffffff, 64) +order = Concat(orderLimb3, orderLimb2, orderLimb1, orderLimb0) + +# Define 2N. +twiceOrderLimb0 = BitVecVal(0x7fa4bd19a06c8282, 64) +twiceOrderLimb1 = BitVecVal(0x755db9cd5e914077, 64) +twiceOrderLimb2 = BitVecVal(0xfffffffffffffffd, 64) +twiceOrderLimb3 = BitVecVal(0xffffffffffffffff, 64) +twiceOrderLimb4 = ONE +twiceOrder = Concat( + twiceOrderLimb4, twiceOrderLimb3, twiceOrderLimb2, + twiceOrderLimb1, twiceOrderLimb0 +) + +# Define 2**256 - N. +orderComplementLimb0 = BitVecVal(0x402da1732fc9bebf, 64) +orderComplementLimb1 = BitVecVal(0x4551231950b75fc4, 64) +orderComplementLimb2 = ONE +orderComplement = Concat(orderComplementLimb2, orderComplementLimb1, orderComplementLimb0) + +# -------------------------------------------------------------------------- +# Prove the constant definitions satisfy the required identities before +# trusting anything that uses them. +# -------------------------------------------------------------------------- + +# Prove twiceOrder = 2N. +order320 = ZeroExt(64, order) +prove(twiceOrder == order320 << 1, "twiceOrder != 2N") + +# -------------------------------------------------------------------------- +# Lemma R-identity: N + c == 2**256 +# +# This lemma justifies replacing multiples of 2**256 with multiples of c when +# working modulo N: +# +# N + c ≡ 2**256 (mod N) +# => 0 + c ≡ 2**256 (mod N) +# => c ≡ 2**256 (mod N) +# -------------------------------------------------------------------------- +orderComplement320 = ZeroExt(128, orderComplement) +prove(order320+orderComplement320 == BitVecVal(1, 320)<<256, + "Lemma R-identity: N+c != 2**256") + +# --------------- +# Model the code. +# --------------- + +discards = [] + +# first reduction: 512 bits -> 385 bits. +h0, t0 = mul64(x[4], orderComplementLimb0) +h1, t1 = mul64(x[4], orderComplementLimb1) +t1, carry = add64(t1, h0, ZERO) +t2, carry = add64(x[4], h1, carry) +t3 = carry + +h0, l0 = mul64(x[5], orderComplementLimb0) +h1, l1 = mul64(x[5], orderComplementLimb1) +t1, carry = add64(t1, l0, ZERO) +t2, carry = add64(t2, h0, carry) +t3, discarded = add64(t3, h1, carry) +discards.append(discarded) +t2, carry = add64(t2, l1, ZERO) +t3, carry = add64(t3, x[5], carry) +t4 = carry + +h0, l0 = mul64(x[6], orderComplementLimb0) +h1, l1 = mul64(x[6], orderComplementLimb1) +t2, carry = add64(t2, l0, ZERO) +t3, carry = add64(t3, h0, carry) +t4, discarded = add64(t4, h1, carry) +discards.append(discarded) +t3, carry = add64(t3, l1, ZERO) +t4, carry = add64(t4, x[6], carry) +t5 = carry + +h0, l0 = mul64(x[7], orderComplementLimb0) +h1, l1 = mul64(x[7], orderComplementLimb1) +t3, carry = add64(t3, l0, ZERO) +t4, carry = add64(t4, h0, carry) +t5, discarded = add64(t5, h1, carry) +discards.append(discarded) +t4, carry = add64(t4, l1, ZERO) +t5, carry = add64(t5, x[7], carry) +t6 = carry + +t0, carry = add64(t0, x[0], ZERO) +t1, carry = add64(t1, x[1], carry) +t2, carry = add64(t2, x[2], carry) +t3, carry = add64(t3, x[3], carry) +t4, carry = add64(t4, ZERO, carry) +t5, carry = add64(t5, ZERO, carry) +t6, discarded = add64(t6, ZERO, carry) +discards.append(discarded) + +# Save for analysis. +t6_saved_1 = t6 +t_full_redux_1 = Concat(ZERO, t6, t5, t4, t3, t2, t1, t0) +first_redux_input_high = Concat(x[7], x[6], x[5], x[4]) + +# second reduction: 385 bits -> 258 bits. +h0, l0 = mul64(t4, orderComplementLimb0) +h1, l1 = mul64(t4, orderComplementLimb1) +l2 = t4 +t0, carry = add64(t0, l0, ZERO) +t1, carry = add64(t1, h0, carry) +t2, carry = add64(t2, h1, carry) +t3, carry = add64(t3, ZERO, carry) +t4 = carry +t1, carry = add64(t1, l1, ZERO) +t2, carry = add64(t2, l2, carry) +t3, carry = add64(t3, ZERO, carry) +t4, discarded = add64(t4, ZERO, carry) +discards.append(discarded) + +h0, l0 = mul64(t5, orderComplementLimb0) +h1, l1 = mul64(t5, orderComplementLimb1) +t1, carry = add64(t1, l0, ZERO) +t2, carry = add64(t2, h0, carry) +t3, carry = add64(t3, h1, carry) +t4, discarded = add64(t4, ZERO, carry) +discards.append(discarded) +t2, carry = add64(t2, l1, ZERO) +t3, carry = add64(t3, t5, carry) +t4, discarded = add64(t4, ZERO, carry) +discards.append(discarded) + +t2, carry = add64(t2, t6*orderComplementLimb0, ZERO) +t3, carry = add64(t3, t6*orderComplementLimb1, carry) +t4, discarded = add64(t4, t6, carry) +discards.append(discarded) + +# Save intermediate results for analysis. +t3_carry_2 = carry +t3_saved_2 = t3 +t4_saved_2 = t4 + +# third reduction: max 258 bits -> t < 2N. +t0, carry = add64(t0, t4*orderComplementLimb0, ZERO) +t1, carry = add64(t1, t4*orderComplementLimb1, carry) +t2, carry = add64(t2, t4, carry) +t3, carry = add64(t3, ZERO, carry) +t4 = carry + +# final reduction: t < 2N -> t < N +s0, borrow = sub64(t0, orderLimb0, ZERO) +s1, borrow = sub64(t1, orderLimb1, borrow) +s2, borrow = sub64(t2, orderLimb2, borrow) +s3, borrow = sub64(t3, orderLimb3, borrow) +_, borrow = sub64(t4, ZERO, borrow) +r0 = If(borrow == ONE, t0, s0) +r1 = If(borrow == ONE, t1, s1) +r2 = If(borrow == ONE, t2, s2) +r3 = If(borrow == ONE, t3, s3) + +# ------- +# Proofs. +# ------- + +# Discarded carries are never set. +prove_no_discarded_carries(discards) + +# Top limb after first reduction is at most 1 (t6 ≤ 1). +# +# This along with the proof there is no discarded carry on t6 after the first +# reduction proves the value does not exceed 385 bits (as expected since the +# complement is 129 bits and there were a max of 256 bits in the upper limbs +# before the first reduction, so 256 + 129 = 385). +prove(ULE(t6_saved_1, ONE), "top limb after 1st reduction > 1") + +# Top limb after second reduction is at most 3 (t4 ≤ 3). +# +# This along with the proof there is no discarded carry on t4 after the second +# reduction proves the value does not exceed 258 bits (as expected because the +# complement is 129 bits and there were 129 max bits left in the upper limbs +# after the first reduction, so 129 + 129 = 258). +prove(ULE(t4_saved_2, BitVecVal(3, 64)), "top limb after 2nd reduction > 3") + +# Top limb after second reduction times complement limb 0 never exceeds uint64. +# ((t4*c0)>>64 == 0) +# +# This proves the regular 64-bit product used in the third reduction does not +# overflow. +t4c0hi, _ = mul64(t4_saved_2, orderComplementLimb0) +prove( + t4c0hi == ZERO, + "high word of top limb from 2nd reduction times complement limb 0 != 0", +) + +# Top limb after second reduction times complement limb 1 never exceeds uint64. +# ((t4*c1)>>64 == 0) +# +# This proves the regular 64-bit product used in the third reduction does not +# overflow. +t4c1hi, _ = mul64(t4_saved_2, orderComplementLimb1) +prove( + t4c1hi == ZERO, + "high word of top limb from 2nd reduction times complement limb 1 != 0", +) + +# Top limb after third reduction is at most 1 (t4 ≤ 1). +# +# This proof is not strictly necessary since the bound is actually tigher per +# the next proof which implicitly proves this fact. However, it is fast to +# prove and including it provides nice symmetry with the previous proofs for the +# bounds of the top limb after each reduction. +prove(ULE(t4, ONE), "top limb after 3rd reduction > 1") + +# Value after third reduction is less than twice the group order (t < 2N). +t = Concat(t4, t3, t2, t1, t0) +prove(ULT(t, twiceOrder), "value after 3rd reduction >= 2N") + +# Fully reduced result is less than the group order (r < N). +r = Concat(r3, r2, r1, r0) +prove(ULT(r, order), "fully reduced result >= N") + +# ----------------------------------------------------------------------------- +# Prove functional congruence: r ≡ x (mod N). +# +# The checks above establish bounds and the safety of discarded carries, +# self-consistency of the constants, and the reduction identity (lemma +# R-identity). +# +# The following lemmas establish the arithmetic identity of each transformation +# stage. Together with the previously proven bounds and reduction identity, +# they imply the overall congruence r ≡ x (mod N). +# +# A direct assertion such as t == x_lo + (x >> 256)*P is intractable because +# it introduces a monolithic 256x256 multiplier that appears nowhere in the +# modeled code, so the SAT solver is forced into bit-level equivalence +# checking between two structurally unrelated multiplier circuits. +# +# Instead, because every accumulation block in the implementation is a fixed +# carry-propagation network over the *outputs* of bits.Mul64, whether or not +# those 128-bit inputs are actually products is irrelevant to the network. Each +# block satisfies its per-block identity for ALL 128-bit values. So each block +# identity is proven here as a universal lemma over unconstrained variables, +# which contains no multipliers at all and reduces to adder-network equivalence +# that the solver dispatches quickly. +# +# The Go function, scalar64Reduce512, computes r in four stages: +# +# 1) First reduction: Fold the high 256 bits of x (x4..x7) into the low 256 bits +# (x0..x3). This produces t_1 = x_lo + x_hi*c which fits in 385 bits (7 +# 64-bit limbs with the 7th limb small). +# 2) Second reduction: Fold t_1's 5th, 6th, and 7th limbs back in the same way. +# This produces t_2 = t_1_lo + t_1_hi*c which fits in 258 bits (5 64-bit +# limbs with the 7th limb small) +# 3) Third reduction: Fold t_2's 5th limb back in the same way. This produces +# t_3 = t_2_lo + t_1_hi*c, where t_3 < 2N. +# 4) Final reduction: One conditional subtraction of N which is valid because +# t_3 < 2N. +# +# The lemmas below chain together to prove r ≡ t3 ≡ t2 ≡ t1 ≡ x (mod N) and 0 ≤ +# r < N which is exactly the definition of x ≡ r (mod N). +def prove_functional_congruence_lemmas(): + """ + Prove r == x (mod N) for scalar64Reduce512. + + x is the 512-bit input (8 limbs x0..x7) and r is the fully-reduced 256-bit + output (4 limbs). + """ + + # Use a wide enough congruence width that nothing in a 512-bit input ever + # wraps around. This is comfortably wider than the largest quantity + # of x_hi*c < 2**385. + CONGRUENCE_WIDTH = 512 + + def cwidth(v): + """Zero-extend v to CONGRUENCE_WIDTH bits.""" + return ZeroExt(CONGRUENCE_WIDTH - v.size(), v) + + def pack(*limbs): + """ + Concatenate the provided big-endian 64-bit limbs into a single + CONGRUENCE_WIDTH-bit value: + + limb0*2**(64*num_limbs) + limb1*2**(64*(num_limbs-1)) + ... + """ + return cwidth(Concat(*limbs)) + + def pack_limbs(limbs): + """ + Concatenate the provided iterable little-endian 64-bit limbs into a + single CONGRUENCE_WIDTH-bit value: + + limbs[0] + limbs[1]*2**64 + limbs[2]*2**128 + ... + """ + return cwidth(Concat(*reversed(limbs))) + + def compose_mul_by_c(p0, p1, xj): + """ + Compose the multi-limb product xj*c from the two 128-bit partial + products produced by multiplying xj with the two complement limbs. + + Since the third complement limb is 1, the full product is: + + p0 + (p1 << 64) + (xj << 128) + + The result is equal to xj*c and is exactly how the implementation + accumulates it. + """ + return cwidth(p0) + (cwidth(p1)<<64) + (cwidth(xj)<<128) + + # Universal lemma variables shared by R-open, R-interior, and + # R-second-interior. + # + # x_j models an arbitrary 64-bit input limb. + # + # p0 and p1 model the 128-bit products returned by bits.Mul64(x[j], c0) and + # bits.Mul64(x[j], c1) (the lower two limbs of c) at some point in the real + # code. They are NOT tied to any specific x[j] or t[k]. + # + # h0 and l0 are the upper and lower halves of p0, respectively. Likewise h1 + # and l1 for p1. + # + # x_j_times_c is the composed value of x[j]*c. + x_j = BitVec("x_j", 64) + p0, p1 = BitVecs("p0 p1", 128) + h0, l0 = split128(p0) + h1, l1 = split128(p1) + x_j_times_c = compose_mul_by_c(p0, p1, x_j) + + # ------------------------------------------------------------------------- + # Lemma R-open: The opening row of the first reduction (the x4 line). + # + # h0, t0 = bits.Mul64(x[4], orderComplementLimb0) + # h1, t1 = bits.Mul64(x[4], orderComplementLimb1) + # t1, carry = bits.Add64(t1, h0, 0) + # t2, carry = bits.Add64(x[4], h1, carry) + # t3 = carry + # + # This only writes fresh limbs and has no discarded carries, so no + # assumptions are needed. + # ------------------------------------------------------------------------- + t0 = l0 + t1, carry = add64(l1, h0, ZERO) + t2, carry = add64(x_j, h1, carry) + t3 = carry + prove(pack(t3,t2,t1,t0) == x_j_times_c, "R-open: state != x4*c") + + # ------------------------------------------------------------------------- + # Lemma R-interior: The shape shared by the x5, x6, and x7 blocks in the + # first reduction. + # + # h0, l0 = bits.Mul64(x[j], orderComplementLimb0) + # h1, l1 = bits.Mul64(x[j], orderComplementLimb1) + # t[j], carry = bits.Add64(t[j], l0, 0) + # t[j+1], carry = bits.Add64(t[j+1], h0, carry) + # t[j+2], _ = bits.Add64(t[j+2], h1, carry) + # t[j+1], carry = bits.Add64(t[j+1], l1, 0) + # t[j+2], carry = bits.Add64(t[j+2], x[j], carry) + # t[j+3] = carry + # + # a0..a3 model the incoming running total in t0..t3. + # + # One carry is discarded and assumed zero. This fact is proven above by the + # discard checks for the real x5/x6/x7 blocks, so this does not smuggle in a + # new assumption. + # ------------------------------------------------------------------------- + a = BitVecs("a0 a1 a2", 64) + w = [None]*4 + w[0], carry = add64(a[0], l0, ZERO) + w[1], carry = add64(a[1], h0, carry) + w[2], discarded = add64(a[2], h1, carry) + w[1], carry = add64(w[1], l1, ZERO) + w[2], carry = add64(w[2], x_j, carry) + w[3] = carry + prove( + pack_limbs(w) == pack_limbs(a) + x_j_times_c, + "R-interior: block != state + xj*c", + assumptions=[discarded == ZERO], + ) + + # ------------------------------------------------------------------------- + # Lemma R-tail: Fold x0..x3 into running total after the blocks above. + # + # t0, carry = bits.Add64(t0, x[0], 0) + # t1, carry = bits.Add64(t1, x[1], carry) + # ... + # t5, carry = bits.Add64(t5, 0, carry) + # t6, _ = bits.Add64(t6, 0, carry) + # + # b0..b6 model the incoming running total in t0..t6. + # + # The final line is modeled here as add64 with its own carry-out assumed + # zero. This fact is proven above by the discard checks for the real total, + # so this does not smuggle in a new assumption. + # ------------------------------------------------------------------------- + b = BitVecs("b0 b1 b2 b3 b4 b5 b6", 64) + w = [None]*7 + w[0], carry = add64(b[0], x[0], ZERO) + w[1], carry = add64(b[1], x[1], carry) + w[2], carry = add64(b[2], x[2], carry) + w[3], carry = add64(b[3], x[3], carry) + w[4], carry = add64(b[4], ZERO, carry) + w[5], carry = add64(b[5], ZERO, carry) + w[6], discarded = add64(b[6], ZERO, carry) + prove( + pack_limbs(w) == pack_limbs(b) + pack_limbs(x[:4]), + "R-tail: post-fold total != state + x_lo", + assumptions=[discarded == ZERO], + ) + + # ------------------------------------------------------------------------- + # Lemma R-second-interior: The carry-propagation shape shared by the t4 and + # t5 blocks of the second reduction. + # + # Unlike R-interior, this block propagates carries through an existing fifth + # accumulator limb rather than terminating with a new carry limb. + # + # h0, l0 = bits.Mul64(t[k], orderComplementLimb0) + # h1, l1 = bits.Mul64(t[k], orderComplementLimb1) + # // h2, l2 = bits.Mul64(t[k], orderComplementLimb2) => h2=0, l2=t[k] + # t[j], carry = bits.Add64(t[j], l0, 0) + # t[j+1], carry = bits.Add64(t[j+1], h0, carry) + # ... + # t[j+4], _ = bits.Add64(t[j+4], 0, carry) + # + # d0..d4 model the incoming running total in t0..t4. + # + # d4 is zero when instantiated for the t4 block. + # + # The t4 additions are modeled as add64 with the carry-out assumed zero. + # This fact is proven above by the discard checks for the real totals, so + # this does not smuggle in a new assumption. + # ------------------------------------------------------------------------- + d = BitVecs("d0 d1 d2 d3 d4", 64) + w = [None]*5 + w[0], carry = add64(d[0], l0, ZERO) + w[1], carry = add64(d[1], h0, carry) + w[2], carry = add64(d[2], h1, carry) + w[3], carry = add64(d[3], ZERO, carry) + w[4], discarded1 = add64(d[4], ZERO, carry) + w[1], carry = add64(w[1], l1, ZERO) + w[2], carry = add64(w[2], x_j, carry) + w[3], carry = add64(w[3], ZERO, carry) + w[4], discarded2 = add64(w[4], ZERO, carry) + prove( + pack_limbs(w) == pack_limbs(d) + x_j_times_c, + "R-second-interior: block != state + xj*c", + assumptions=[discarded1 == ZERO, discarded2 == ZERO], + ) + + # ------------------------------------------------------------------------- + # Lemma R-second-tail: The t6 block of the second reduction. + # + # t2, carry = bits.Add64(t2, t6*orderComplementLimb0, 0) + # t3, carry = bits.Add64(t3, t6*orderComplementLimb1, carry) + # t4, _ = bits.Add64(t4, t6, carry) + # + # f0..f2 model the incoming running total in t2..t4. + # + # t6 models the incoming carry limb produced by the preceding block. + # + # The real implementation uses plain (truncating) 64-bit products here rather + # than bits.Mul64. Since t6 ≤ 1, these products are exact, so modeling them + # as ordinary 64-bit multiplications is equivalent to the implementation. + # + # Assumes t6 ≤ 1 and one discarded carry is 0. These facts are both proven + # above, so this is not smuggling in new assumptions. + # ------------------------------------------------------------------------- + f = BitVecs("f0 f1 f2", 64) + t6 = BitVec("t6", 64) + w = [None]*3 + w[0], carry = add64(f[0], t6*orderComplementLimb0, ZERO) + w[1], carry = add64(f[1], t6*orderComplementLimb1, carry) + w[2], discarded = add64(f[2], t6, carry) + t6_times_c = compose_mul_by_c( + ZeroExt(64, t6) * ZeroExt(64, orderComplementLimb0), + ZeroExt(64, t6) * ZeroExt(64, orderComplementLimb1), + t6, + ) + prove( + pack_limbs(w) == pack_limbs(f) + t6_times_c, + "R-second-tail: block != state + t6*c", + assumptions=[ULE(t6, 1), discarded == ZERO], + ) + + # ------------------------------------------------------------------------- + # Lemma R-third: The third reduction. + # + # t0, carry = bits.Add64(t0, t4*orderComplementLimb0, 0) + # t1, carry = bits.Add64(t1, t4*orderComplementLimb1, carry) + # t2, carry = bits.Add64(t2, t4, carry) + # t3, carry = bits.Add64(t3, 0, carry) + # t4 = carry + # + # g0..g3 model the incoming running total in t0..t3. Although the product + # is added only into t0..t2, the carry chain propagates through t3. + # + # t4 models the incoming carry limb produced by the preceding block. + # + # The real implementation uses plain (truncating) 64-bit products here rather + # than bits.Mul64. Since t4 ≤ 3, these products are exact, so modeling them + # as ordinary 64-bit multiplications is equivalent to the implementation. + # + # Assumes t4 ≤ 3. This fact is proven above, so this is not smuggling in + # a new assumption. + # ------------------------------------------------------------------------- + g = BitVecs("g0 g1 g2 g3", 64) + t4 = BitVec("t4", 64) + w = [None]*5 + w[0], carry = add64(g[0], t4*orderComplementLimb0, ZERO) + w[1], carry = add64(g[1], t4*orderComplementLimb1, carry) + w[2], carry = add64(g[2], t4, carry) + w[3], carry = add64(g[3], ZERO, carry) + w[4] = carry + t4_times_c = compose_mul_by_c( + ZeroExt(64, t4) * ZeroExt(64, orderComplementLimb0), + ZeroExt(64, t4) * ZeroExt(64, orderComplementLimb1), + t4, + ) + prove( + pack_limbs(w) == pack_limbs(g) + t4_times_c, + "R-third: 3rd reduction != state + t4*c", + assumptions=[ULE(t4, 3)], + ) + + # ------------------------------------------------------------------------- + # Lemma R-final: The conditional-subtraction final reduction. + # + # s0, borrow = bits.Sub64(t0, orderLimb0, 0) + # ... + # _, borrow = bits.Sub64(t4, 0, borrow) + # r[i] = constantTimeSelect64(borrow, t[i], s[i]) + # + # z0..z4 model t0..t4 going into this stage. + # + # This is only valid because t < 2N as otherwise t - N could still exceed N + # and the 4-limb result wouldn't represent it. This fact is proven above + # for the real total, so this does not smuggle in a new assumption. + # ------------------------------------------------------------------------- + n = cwidth(order) + z = BitVecs("z0 z1 z2 z3 z4", 64) + s = [None]*4 + s[0], borrow = sub64(z[0], orderLimb0, ZERO) + s[1], borrow = sub64(z[1], orderLimb1, borrow) + s[2], borrow = sub64(z[2], orderLimb2, borrow) + s[3], borrow = sub64(z[3], orderLimb3, borrow) + _, borrow = sub64(z[4], ZERO, borrow) + t = pack_limbs(z) + r = If(borrow == ONE, t, pack_limbs(s)) + prove( + r == If(ULT(t, n), t, t - n), + "R-final: conditional-subtract result != t mod N", + assumptions=[ULT(t, cwidth(twiceOrder))], + ) + +prove_functional_congruence_lemmas() \ No newline at end of file diff --git a/dcrec/secp256k1/modnscalar.go b/dcrec/secp256k1/modnscalar.go index 6ccaf6306b..0eaf792c17 100644 --- a/dcrec/secp256k1/modnscalar.go +++ b/dcrec/secp256k1/modnscalar.go @@ -387,6 +387,10 @@ func (s *ModNScalar) Add(val *ModNScalar) *ModNScalar { // scalar64Reduce512 reduces a 512-bit little-endian limb array modulo the group // order in constant time and stores the result in r. func scalar64Reduce512(r *[4]uint64, x *[8]uint64) { + // This algorithm has been formally verified, including its intermediate + // bounds, carry assumptions, and functional correctness. The verification + // artifacts are available in internal/proofs. + // The overall strategy employed here is: // 1) Start with the full unreduced 512-bit product of the two scalars. // 2) Reduce the result modulo the group order via Crandall reduction From f4b5ba907bf496c8dcef051b632a539fd5437e21 Mon Sep 17 00:00:00 2001 From: Dave Collins Date: Sun, 19 Jul 2026 08:06:40 -0500 Subject: [PATCH 07/17] secp256k1: Add width-5 windowed NAF recoder. With the ultimate goal of optimizing scalar multiplication by reducing the average number of point additions, this commit adds code to convert a balanced scalar representative into width-5 windowed non-adjacent form (wNAF). The current scalar multiplication implementation uses ordinary non-adjacent form (NAF), which corresponds to a window width of 2. This commit introduces only the width-5 wNAF recoder and does not yet integrate it into scalar multiplication. Future commits will add comprehensive tests and integrate the new implementation. --- dcrec/secp256k1/curve.go | 228 +++++++++++++++++++++++++++++++++++++++ 1 file changed, 228 insertions(+) diff --git a/dcrec/secp256k1/curve.go b/dcrec/secp256k1/curve.go index 750b57a327..fe93eb8da6 100644 --- a/dcrec/secp256k1/curve.go +++ b/dcrec/secp256k1/curve.go @@ -943,6 +943,234 @@ func splitK(k *ModNScalar) (ModNScalar, ModNScalar) { return k1, k2 } +// wNAFWidth is the window width of the NAF representation used in scalar +// multiplication. +// +// It is important to note that this constant primarily exists to improve code +// readability, to parameterize parts of the scalar multiplication +// implementation, and to allow tests to assert invariants. +// +// However, the overall implementation is intentionally not fully parameterized +// based on this constant because it is carefully optimized for this specific +// window width and those optimizations involve taking advantage of known bounds +// that an arbitrary window size would violate. +// +// Choosing a new window size involves a variety of tradeoffs that must be +// carefully analyzed. For example, recoding costs, precomputation costs, and +// average density. +const wNAFWidth = 5 + +// wnaf5RecodeCodes maps the low five bits of an integer to the encoded +// representative of the width-5 wNAF digit that should be emitted. +// +// These values are not the signed digits themselves. Rather, each is an +// encoded representative that serves directly as an index into precomputed +// tables. See [wnafScalar] for the encoding details. +// +// For each odd residue, the represented signed digit is the unique value in +// {±1, ±3, ±5, ±7, ±9, ±11, ±13, ±15} that is congruent to the residue modulo +// 32: +// +// residue ≡ digit (mod 32) +// +// Thus, subtracting the represented signed digit always produces an integer +// that is divisible by 32. This allows the optimized recoding algorithm to +// skip five one-bit iterations of the canonical algorithm at once. +// +// Even residues always map to zero because they are already divisible by two +// and therefore correspond only to implicit zero digits. +// +// The corresponding arithmetic adjustment for each encoded representative is +// provided by [wnaf5RecodeAdjustments]. +var wnaf5RecodeCodes = [32]uint8{ + 0, 1, 0, 2, 0, 3, 0, 4, 0, 5, 0, 6, 0, 7, 0, 8, + 0, 16, 0, 15, 0, 14, 0, 13, 0, 12, 0, 11, 0, 10, 0, 9, +} + +// wnaf5RecodeAdjustments maps the low five bits of an integer to the signed +// adjustment that must be added to the working integer before shifting during +// the width-5 wNAF recoding algorithm. +// +// While [wnaf5RecodeCodes] produces the encoded representative that is stored +// in the resulting wNAF encoding, this table provides the corresponding +// arithmetic adjustment. The represented signed digit d satisfies: +// +// residue ≡ d (mod 32) +// +// and the recoding algorithm updates the working integer via: +// +// n = (n - d) / 2 +// +// Rather than subtracting d directly, this table stores -d encoded as an +// unsigned 64-bit two's complement value. Negative adjustments therefore +// appear as values such as: +// +// -1 -> ^uint64(0) +// -3 -> ^uint64(2) +// -5 -> ^uint64(4) +// -7 -> ^uint64(6) +// -9 -> ^uint64(8) +// -11 -> ^uint64(10) +// -13 -> ^uint64(12) +// -15 -> ^uint64(14) +// +// The adjustment is sign-extended across the upper limbs during the multiword +// addition which allows the hot loop to remain branch free. +var wnaf5RecodeAdjustments = [32]uint64{ + 0, ^uint64(0), 0, ^uint64(2), 0, ^uint64(4), 0, ^uint64(6), + 0, ^uint64(8), 0, ^uint64(10), 0, ^uint64(12), 0, ^uint64(14), + 0, 15, 0, 13, 0, 11, 0, 9, 0, 7, 0, 5, 0, 3, 0, 1, +} + +// wnafScalar represents a non-negative integer less than 2^129 encoded in +// width-5 windowed non-adjacent form (wNAF). +// +// wNAF is a signed-digit representation where the valid digits are 0, ±1, ±3, +// ±5, ±7, ±9, ±11, ±13, and ±15. +// +// Each entry corresponds to one bit position and encodes the signed digit. +// +// The encoding from code to digit is: +// 0 -> 0 +// 1 -> +1 +// 2 -> +3 +// 3 -> +5 +// 4 -> +7 +// 5 -> +9 +// 6 -> +11 +// 7 -> +13 +// 8 -> +15 +// 9 -> -1 +// 10 -> -3 +// 11 -> -5 +// 12 -> -7 +// 13 -> -9 +// 14 -> -11 +// 15 -> -13 +// 16 -> -15 +// +// The encoding is intentionally not using a more compact form so that it serves +// directly as an index into precomputed tables which allows simplification of +// hot paths. +// +// Formally, letting c denote the stored code value, the decoded signed digit +// d(c) is given by the piecewise function: +// +// {0, where c = 0 +// d(c) = {2c - 1, where 0 < c ≤ 8 +// {17 - 2c, where 8 < c ≤ 16 +// +// The array contains one extra entry because a recoding may produce a carry +// into an additional most-significant digit. +type wnafScalar struct { + codes [130]uint8 + bits uint8 +} + +// wnaf takes a non-negative integer less than 2^129 and returns its width-5 +// windowed non-adjacent form (wNAF) which is a unique signed-digit +// representation such that non-zero digits are separated by at least 4 zeroes. +// See [wnafScalar] for details on how the representation is encoded and how to +// interpret it. +// +// Width-5 wNAF is useful because, on average, only about one in every six +// digits is non-zero. +// +// This property is particularly beneficial for optimizing elliptic curve point +// multiplication because it greatly reduces the number of point additions at +// the cost of precomputing a few odd multiples of the point and, in the worst +// case, one additional point doubling due to a carry introduced by the +// recoding. This is an excellent tradeoff because subtraction of points has +// the same computational complexity as addition of points and point doubling is +// faster than both. +func wnaf(k *ModNScalar) wnafScalar { + const ( + // This is intentionally not using the package constant [wNAFWidth] for + // the window size for the reasons stated by its documentation. + // + // In particular, this implementation is carefully optimized for this + // specific width and involves assumptions that an arbitrary window size + // might violate. + // + // Using a separate constant helps make it clear that updating the + // window size here requires extra care to assert correctness. + windowSize = 5 + windowMask = (1 << windowSize) - 1 + ) + + // The arithmetic for wNAF recoding is performed over the ordinary integers, + // not modulo the curve order. Moreover, this method is required to be + // called with the scalar's balanced integer representative, which is known + // to fit within 129 bits. + // + // Consequently, the scalar is first converted to a uint192 with three + // 64-bit limbs where the top limb is at most 1. + t0, t1, t2 := k.n[0], k.n[1], k.n[2] + + // This implementation is based on the standard width-w NAF recoding + // algorithm presented as Algorithm 3.35 in [GECC]. However, it has been + // modified to skip runs of zero digits instead of processing one bit at a + // time. + // + // The optimizations exploit the fact that runs of zero digits are implicit + // in the output representation. + // + // Also, note that the last non-zero bit is initialized to the maximum value + // so that adding 1 at the end to account for the exclusive endpoint wraps + // around to 0 when no digits are emitted. + var result wnafScalar + var c uint64 + var bit uint8 + var lastNonZeroBit = ^uint8(0) + for t0|t1|t2 != 0 { + // The next five iterations of the canonical bit-by-bit algorithm would + // emit only zero when the low window is zero, so skip directly to the + // next window boundary in that case. + if residue := uint8(t0 & windowMask); residue != 0 { + // Skip the zero digits that the canonical bit-by-bit algorithm + // would emit before reaching the next odd value. + // + // Since 0 < residue < 32, shift is guaranteed to satisfy the + // following bounds: + // + // 0 ≤ shift < 5 + shift := uint8(bits.TrailingZeros8(residue)) + t0 = t0>>shift | t1<<(64-shift) + t1 = t1>>shift | t2<<(64-shift) + t2 >>= shift + bit += shift + + residue = uint8(t0 & windowMask) + result.codes[bit] = wnaf5RecodeCodes[residue] + lastNonZeroBit = bit + + // The selected digit is congruent to the low five bits modulo 32. + // Therefore, adding the stored adjustment (which represents the + // negation of that digit) makes the value divisible by 32. + // + // Negative adjustments are stored in two's complement. Sign + // extending the value across the upper limbs allows the multiword + // addition to operate without branches. + adjustment := wnaf5RecodeAdjustments[residue] + signExtended := uint64(int64(adjustment) >> 63) + t0, c = bits.Add64(t0, adjustment, 0) + t1, c = bits.Add64(t1, signExtended, c) + t2, _ = bits.Add64(t2, signExtended, c) + } + + // Divide by 32 to prepare for the next window. + // + // The adjustment above guarantees the current value is divisible by 32. + t0 = t0>>windowSize | t1<<(64-windowSize) + t1 = t1>>windowSize | t2<<(64-windowSize) + t2 >>= windowSize + bit += windowSize + } + + result.bits = lastNonZeroBit + 1 + return result +} + // nafScalar represents a positive integer up to a maximum value of 2^256 - 1 // encoded in non-adjacent form. // From b3ad956641e0b52f507eeef87e4443ecf21a0085 Mon Sep 17 00:00:00 2001 From: Dave Collins Date: Sun, 19 Jul 2026 08:06:42 -0500 Subject: [PATCH 08/17] secp256k1: Add wNAF tests. This adds comprehensive tests for the new wNAF implementation to help verify correctness. The tests assert the resulting encoding satisfies the required mathematical properties and reconstructs the original value for both deterministic edge cases and randomized balanced scalar representatives. The deterministic edge cases exercise enough small values to exhaustively test every residue modulo 2^w multiple times, and include additional edge cases such as carries introduced by the recoding, maximum values, and alternating bit patterns. --- dcrec/secp256k1/curve_test.go | 172 ++++++++++++++++++++++++++++++++++ 1 file changed, 172 insertions(+) diff --git a/dcrec/secp256k1/curve_test.go b/dcrec/secp256k1/curve_test.go index d3eb497f29..a92a6d1ca0 100644 --- a/dcrec/secp256k1/curve_test.go +++ b/dcrec/secp256k1/curve_test.go @@ -753,6 +753,178 @@ func TestNAFRandom(t *testing.T) { } } +// decodeWNAFDigit converts the provided wNAF code to the signed digit it +// represents. +func decodeWNAFDigit(code uint8) int8 { + const ( + positiveCodeLimit = 1 << (wNAFWidth - 2) + maxCode = 1 << (wNAFWidth - 1) + ) + switch { + case code == 0: + return 0 + case code <= positiveCodeLimit: + return 2*int8(code) - 1 + } + return (maxCode + 1) - 2*int8(code) +} + +// checkWNAFEncoding returns an error if the provided windowed NAF encoding does +// not adhere to the encoding requirements or does not reconstruct the provided +// original value. +func checkWNAFEncoding(s *wnafScalar, origValue *big.Int) error { + // Ensure the reported number of used bits does not exceed the size of the + // array. + if int(s.bits) > len(s.codes) { + return fmt.Errorf("bits %d > %d", s.bits, len(s.codes)) + } + + // wNAF must not have a leading zero and the individual codes must not + // exceed the max allowed code for the window width. + codes := s.codes[:s.bits] + if len(codes) > 0 && codes[len(codes)-1] == 0 { + return fmt.Errorf("leading zero in encoding: %v", codes) + } + const maxCode = 1 << (wNAFWidth - 1) + for _, code := range codes { + if code > maxCode { + return fmt.Errorf("found code %d > %d", code, maxCode) + } + } + + // Ensure each non-zero digit is separated by at least the width of the + // window. + const minDistance = wNAFWidth + prevNonZeroBit := -minDistance + for bit, code := range codes { + if code == 0 { + continue + } + + if distance := bit - prevNonZeroBit; distance < minDistance { + return fmt.Errorf("non-zero digits at bit pos %d and %d are only "+ + "%d bits apart", prevNonZeroBit, bit, distance) + } + + prevNonZeroBit = bit + } + + // Reconstruct the represented value and ensure it matches the original + // value. + sum := new(big.Int) + for bit, code := range codes { + digit := decodeWNAFDigit(code) + term := big.NewInt(int64(digit)) + term.Lsh(term, uint(bit)) + sum.Add(sum, term) + } + if origValue.Cmp(sum) != 0 { + return fmt.Errorf("failed to reconstruct orig value: got %v, want %v", + sum, origValue) + } + + return nil +} + +// TestWNAF ensures encoding various edge cases and values to width-w windowed +// non-adjacent form, where w is [wNAFWidth], produces valid results. +func TestWNAF(t *testing.T) { + type wnafTest struct { + name string // test description + in string // hex encoded test value + } + extraTests := []wnafTest{{ + name: "leading zeroes", + in: "002f20569b90697ad471c1be6107814f", + }, { + name: "highest allowed bit only", + in: "100000000000000000000000000000000", + }, { + name: "largest 129-bit value", + in: "1ffffffffffffffffffffffffffffffff", + }, { + name: "130 bits when NAF encoded", + in: "1f0000000000000000000000000000001", + }, { + name: "alternating bits (0xa5)", + in: "a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5", + }, { + name: "alternating bits (0x5a)", + in: "5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a", + }, { + name: "all ones", + in: "ffffffffffffffffffffffffffffffff", + }, { + name: "first term of balanced length-two representation #1", + in: "b776e53fb55f6b006a270d42d64ec2b1", + }, { + name: "second term balanced length-two representation #1", + in: "d6cc32c857f1174b604eefc544f0c7f7", + }, { + name: "first term of balanced length-two representation #2", + in: "45c53aa1bb56fcd68c011e2dad6758e4", + }, { + name: "second term of balanced length-two representation #2", + in: "a2e79d200f27f2360fba57619936159b", + }} + + // Add tests for the first small values since they collectively exercise + // every possible residue modulo 2^w, where w is the window width, multiple + // times. Then append the extra deterministic tests for edge cases. + const smallValues = 2*(1< Date: Sun, 19 Jul 2026 08:06:43 -0500 Subject: [PATCH 09/17] secp256k1: Add wNAF benchmark. --- dcrec/secp256k1/curve_bench_test.go | 12 ++++++++++++ 1 file changed, 12 insertions(+) diff --git a/dcrec/secp256k1/curve_bench_test.go b/dcrec/secp256k1/curve_bench_test.go index 5100ed8966..5bcd50e171 100644 --- a/dcrec/secp256k1/curve_bench_test.go +++ b/dcrec/secp256k1/curve_bench_test.go @@ -229,6 +229,18 @@ func BenchmarkNAF(b *testing.B) { } } +// BenchmarkWNAF benchmarks conversion of a non-negative integer into its +// width-w windowed non-adjacent form representation, where w is [wNAFWidth]. +func BenchmarkWNAF(b *testing.B) { + k := mustModNScalar("8447e288d34b360bc885cb8ce7c00575") + + b.ReportAllocs() + b.ResetTimer() + for i := 0; i < b.N; i++ { + wnaf(k) + } +} + // BenchmarkJacobianPointEquivalency benchmarks determining if two Jacobian // points represent the same affine point. func BenchmarkJacobianPointEquivalency(b *testing.B) { From 940dde829c5e8369778d9b72449af7a3181abb6a Mon Sep 17 00:00:00 2001 From: Dave Collins Date: Sun, 19 Jul 2026 08:06:45 -0500 Subject: [PATCH 10/17] secp256k1: Integrate width-5 wNAF in scalar mult. MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit This updates the non-constant-time scalar multiplication implementation to use the new width-5 wNAF recoding. In particular, because wNAF requires precomputed odd multiples that were not needed by the previous implementation, this introduces a new dedicated precomputed table type along with a new helper to construct the required multiples. The previous implementation only needed to swap between two precomputed points depending on the scalar signs. However, since the new wNAF implementation uses larger precomputed tables, it instead tracks whether each decomposed scalar was negated and applies the corresponding sign adjustment to the appropriate half of each table. The following benchmark shows a before and after comparison of scalar multiplication: name old time/op new time/op delta --------------------------------------------------------------------------- ScalarMultNonConst 94.4µs ± 1% 84.3µs ± 1% -10.75% (p=0.000 n=10+10) --- dcrec/secp256k1/curve.go | 243 ++++++++++++++++++++++++++------------- 1 file changed, 163 insertions(+), 80 deletions(-) diff --git a/dcrec/secp256k1/curve.go b/dcrec/secp256k1/curve.go index fe93eb8da6..1ce1f1bbf4 100644 --- a/dcrec/secp256k1/curve.go +++ b/dcrec/secp256k1/curve.go @@ -1307,6 +1307,111 @@ func naf(k []byte) nafScalar { return result } +// wNAFPrecompTableSize is the size of the wNAF precomputed point table. +const wNAFPrecompTableSize = 1 << (wNAFWidth - 1) + +// wNAFPrecompTable houses the precomputed odd point multiples used during wNAF +// scalar multiplication. +// +// The first half of the table (starting at index 1) contains the positive odd +// multiples: +// +// P, 3P, 5P, ... +// +// While the second half of the table contains the corresponding negatives: +// +// -P, -3P, -5P, ... +// +// Entry 0 is intentionally left unused so the encoded wNAF digit can be used +// directly as an array index. +type wNAFPrecompTable [wNAFPrecompTableSize + 1]JacobianPoint + +// wNAFNegateStart returns the index of the half of a precomputed table whose +// points must be negated in order to match the sign adjustment made to the +// corresponding scalar. +// +// When the scalar has been negated, the first half is negated so it represents +// the negative odd multiples while the second half remains positive. +// Otherwise, the opposite arrangement is used. +func wNAFNegateStart(isScalarNegated bool) int { + if isScalarNegated { + return 1 + } + return wNAFPrecompTableSize/2 + 1 +} + +// buildOddPrecomputeTables constructs the odd multiple tables for both P and +// φ(P) and arranges the positive and negative representatives so they can be +// indexed directly by encoded wNAF digits. +func buildOddPrecomputeTables(p *JacobianPoint, pTable, phiTable *wNAFPrecompTable, k1Neg, k2Neg bool) { + // oddMultiplesPerHalf is the total number of odd multiples in each half of + // the table. + const oddMultiplesPerHalf = wNAFPrecompTableSize / 2 + + // Build the positive odd multiples of P. + // + // Width-w wNAF only ever emits odd signed digits, so only odd multiples of + // the point need to be precomputed. + // + // Both positive and negative entries are ultimately needed. The second + // half is initially a copy of the positive half and the appropriate half is + // negated below depending on the sign of k1. + // + // The negation is deferred until after the φ(P) table is built because that + // table is constructed by applying φ to the corresponding positive odd + // multiples of P. + var twoP JacobianPoint + pTable[1].Set(p) + DoubleNonConst(&pTable[1], &twoP) + for i := 1; i < oddMultiplesPerHalf; i++ { + AddNonConst(&twoP, &pTable[i], &pTable[i+1]) + } + copy(pTable[oddMultiplesPerHalf+1:], pTable[1:oddMultiplesPerHalf+1]) + + // Build the corresponding odd multiples of φ(P). + // + // As above, both positive and negative entries are ultimately needed. The + // sign of this table is chosen independently from the P table because it + // depends only on the sign of k2. + // + // Note that φ is a group homomorphism, so: + // + // φ(n*P) = n*φ(P) + // + // This allows the second table to be obtained by applying φ to each entry + // of the first table instead of computing each odd multiple independently. + // + // NOTE: φ(x,y) = (βx,y). The Jacobian z coordinates are the same, so this + // math goes through. + for i := 1; i <= oddMultiplesPerHalf; i++ { + phiTable[i].Set(&pTable[i]) + phiTable[i].X.Mul(endoBeta).Normalize() + } + copy(phiTable[oddMultiplesPerHalf+1:], phiTable[1:oddMultiplesPerHalf+1]) + + // Negate whichever half of the table corresponds to the sign adjustment + // that was applied to k1. + // + // This maintains the invariant: + // + // k1*P = -k1*-P + negateStart := wNAFNegateStart(k1Neg) + for i := negateStart; i < negateStart+oddMultiplesPerHalf; i++ { + pTable[i].Y.Negate(1).Normalize() + } + + // Apply the same transformation independently to the φ(P) table based on + // the sign adjustment that was applied to k2. + // + // Similarly, this maintains the invariant: + // + // k2*φ(P) = -k2*-φ(P) + negateStart = wNAFNegateStart(k2Neg) + for i := negateStart; i < negateStart+oddMultiplesPerHalf; i++ { + phiTable[i].Y.Negate(1).Normalize() + } +} + // ScalarMultNonConst multiplies k*P where k is a scalar modulo the curve order // and P is a point in Jacobian projective coordinates and stores the result in // the provided Jacobian point. @@ -1354,25 +1459,6 @@ func ScalarMultNonConst(k *ModNScalar, point, result *JacobianPoint) { // See section 3.5 in [GECC] for a more rigorous treatment. // ------------------------------------------------------------------------- - // Per above, the main equation here to remember is: - // k*P = k1*P + k2*φ(P) - // - // p1 below is P in the equation while p2 is φ(P) in the equation. - // - // NOTE: φ(x,y) = (β*x,y). The Jacobian z coordinates are the same, so this - // math goes through. - // - // Also, calculate -p1 and -p2 for use in the NAF optimization. - p1, p1Neg := new(JacobianPoint), new(JacobianPoint) - p1.Set(point) - p1Neg.Set(p1) - p1Neg.Y.Negate(1).Normalize() - p2, p2Neg := new(JacobianPoint), new(JacobianPoint) - p2.Set(p1) - p2.X.Mul(endoBeta).Normalize() - p2Neg.Set(p2) - p2Neg.Y.Negate(1).Normalize() - // Decompose k into k1 and k2 such that k = k1 + k2*λ (mod n) where k1 and // k2 are around half the bit length of k in order to halve the number of EC // operations. @@ -1384,94 +1470,91 @@ func ScalarMultNonConst(k *ModNScalar, point, result *JacobianPoint) { // which means that when they would otherwise be a small negative magnitude // they will instead be a large positive magnitude. Since the goal is for // the scalars to have a small magnitude to achieve a performance boost, use - // their negation when they are greater than the half order of the group and - // flip the positive and negative values of the corresponding point that - // will be multiplied by to compensate. + // their negation when they are greater than the half order of the group. + // + // In order to compensate, the positive and negative values of the + // corresponding points that will be multiplied by are flipped later. // // In other words, transform the calc when k1 is over the half order to: // k1*P = -k1*-P // // Similarly, transform the calc when k2 is over the half order to: // k2*φ(P) = -k2*-φ(P) + var k1Neg, k2Neg bool k1, k2 := splitK(k) if k1.IsOverHalfOrder() { k1.Negate() - p1, p1Neg = p1Neg, p1 + k1Neg = true } if k2.IsOverHalfOrder() { k2.Negate() - p2, p2Neg = p2Neg, p2 + k2Neg = true } - // Convert k1 and k2 into their NAF representations since NAF has a lot more - // zeros overall on average which minimizes the number of required point - // additions in exchange for a mix of fewer point additions and subtractions - // at the cost of one additional point doubling. + // Per above, the main equation here to remember is: + // k*P = k1*P + k2*φ(P) + // + // The simultaneous multiplication therefore needs two sets of precomputed + // odd multiples: + // + // P, 3P, 5P, ... + // + // and: + // + // φ(P), 3φ(P), 5φ(P), ... + var pPrecomps, phiPrecomps wNAFPrecompTable + buildOddPrecomputeTables(point, &pPrecomps, &phiPrecomps, k1Neg, k2Neg) + + // Convert k1 and k2 into their windowed NAF representations since they have + // a lot more zeroes overall on average which greatly reduces the number of + // point additions at the cost of precomputing a few odd multiples of the + // point and, in the worst case, one additional point doubling due to a + // carry introduced by the recoding. // // This is an excellent tradeoff because subtraction of points has the same // computational complexity as addition of points and point doubling is // faster than both. // // Concretely, on average, 1/2 of all bits will be non-zero with the normal - // binary representation whereas only 1/3rd of the bits will be non-zero - // with NAF. - // - // The Pos version of the bytes contain the +1s and the Neg versions contain - // the -1s. - k1Bytes, k2Bytes := k1.Bytes(), k2.Bytes() - k1NAF, k2NAF := naf(k1Bytes[:]), naf(k2Bytes[:]) - k1PosNAF, k1NegNAF := k1NAF.Pos(), k1NAF.Neg() - k2PosNAF, k2NegNAF := k2NAF.Pos(), k2NAF.Neg() - k1Len, k2Len := len(k1PosNAF), len(k2PosNAF) + // binary representation whereas only 1/6 of the bits will be non-zero with + // width-5 wNAF. + k1NAF, k2NAF := wnaf(&k1), wnaf(&k2) - // Add left-to-right using the NAF optimization. See algorithm 3.77 from - // [GECC]. + // Add left-to-right using the endomorphism and wNAF optimizations. See + // algorithms 3.36 and 3.77 from [GECC]. // // Point Q = ∞ (point at infinity). + // + // Like ordinary binary scalar multiplication, the accumulator is doubled + // once per processed bit position. However, instead of each bit + // contributing either 0 or 1 copies of the point, each non-zero wNAF digit + // contributes a precomputed odd multiple while the preceding doublings + // supply the required power-of-two scaling. + // + // ±P, ±3P, ±5P, ... + // + // Since non-zero digits are sparse and separated by several zero digits, + // relatively few point additions are required compared to the ordinary + // binary representation. var q JacobianPoint - m := k1Len - if m < k2Len { - m = k2Len + maxBits := k1NAF.bits + if maxBits < k2NAF.bits { + maxBits = k2NAF.bits } - for i := 0; i < m; i++ { - // Since k1 and k2 are potentially different lengths and the calculation - // is being done left to right, pad the front of the shorter one with - // 0s. - var k1BytePos, k1ByteNeg, k2BytePos, k2ByteNeg byte - if i >= m-k1Len { - k1BytePos, k1ByteNeg = k1PosNAF[i-(m-k1Len)], k1NegNAF[i-(m-k1Len)] - } - if i >= m-k2Len { - k2BytePos, k2ByteNeg = k2PosNAF[i-(m-k2Len)], k2NegNAF[i-(m-k2Len)] + for bit := maxBits; bit > 0; bit-- { + // Q = 2 * Q + DoubleNonConst(&q, &q) + + // The encoded wNAF digit of k1 serves directly as an index into the + // precomputed odd multiples of P. + if code := k1NAF.codes[bit-1]; code != 0 { + AddNonConst(&q, &pPrecomps[code], &q) } - for mask := uint8(1 << 7); mask > 0; mask >>= 1 { - // Q = 2 * Q - DoubleNonConst(&q, &q) - - // Add or subtract the first point based on the signed digit of the - // NAF representation of k1 at this bit position. - // - // +1: Q = Q + p1 - // -1: Q = Q - p1 - // 0: Q = Q (no change) - if k1BytePos&mask == mask { - AddNonConst(&q, p1, &q) - } else if k1ByteNeg&mask == mask { - AddNonConst(&q, p1Neg, &q) - } - - // Add or subtract the second point based on the signed digit of the - // NAF representation of k2 at this bit position. - // - // +1: Q = Q + p2 - // -1: Q = Q - p2 - // 0: Q = Q (no change) - if k2BytePos&mask == mask { - AddNonConst(&q, p2, &q) - } else if k2ByteNeg&mask == mask { - AddNonConst(&q, p2Neg, &q) - } + // Likewise, the encoded wNAF digit of k2 serves directly as an index + // into the precomputed odd multiples of φ(P). + if code := k2NAF.codes[bit-1]; code != 0 { + AddNonConst(&q, &phiPrecomps[code], &q) } } From 561688b41ba19e376fd7280676d9e9d8118a51b9 Mon Sep 17 00:00:00 2001 From: Dave Collins Date: Sun, 19 Jul 2026 08:49:14 -0500 Subject: [PATCH 11/17] secp256k1: Remove unused NAF code. This removes the old ordinary NAF implementation that is no longer used now that the implementation has been updated to use a new width-5 wNAF recoder instead. It also removes the associated tests and benchmark accordingly. --- dcrec/secp256k1/curve.go | 136 ---------------------------- dcrec/secp256k1/curve_bench_test.go | 13 --- dcrec/secp256k1/curve_test.go | 119 ------------------------ 3 files changed, 268 deletions(-) diff --git a/dcrec/secp256k1/curve.go b/dcrec/secp256k1/curve.go index 1ce1f1bbf4..fb9da3d265 100644 --- a/dcrec/secp256k1/curve.go +++ b/dcrec/secp256k1/curve.go @@ -1171,142 +1171,6 @@ func wnaf(k *ModNScalar) wnafScalar { return result } -// nafScalar represents a positive integer up to a maximum value of 2^256 - 1 -// encoded in non-adjacent form. -// -// NAF is a signed-digit representation where each digit can be +1, 0, or -1. -// -// In order to efficiently encode that information, this type uses two arrays, a -// "positive" array where set bits represent the +1 signed digits and a -// "negative" array where set bits represent the -1 signed digits. 0 is -// represented by neither array having a bit set in that position. -// -// The Pos and Neg methods return the aforementioned positive and negative -// arrays, respectively. -type nafScalar struct { - // pos houses the positive portion of the representation. An additional - // byte is required for the positive portion because the NAF encoding can be - // up to 1 bit longer than the normal binary encoding of the value. - // - // neg houses the negative portion of the representation. Even though the - // additional byte is not required for the negative portion, since it can - // never exceed the length of the normal binary encoding of the value, - // keeping the same length for positive and negative portions simplifies - // working with the representation and allows extra conditional branches to - // be avoided. - // - // start and end specify the starting and ending index to use within the pos - // and neg arrays, respectively. This allows fixed size arrays to be used - // versus needing to dynamically allocate space on the heap. - // - // NOTE: The fields are defined in the order that they are to minimize the - // padding on 32-bit and 64-bit platforms. - pos [33]byte - start, end uint8 - neg [33]byte -} - -// Pos returns the bytes of the encoded value with bits set in the positions -// that represent a signed digit of +1. -func (s *nafScalar) Pos() []byte { - return s.pos[s.start:s.end] -} - -// Neg returns the bytes of the encoded value with bits set in the positions -// that represent a signed digit of -1. -func (s *nafScalar) Neg() []byte { - return s.neg[s.start:s.end] -} - -// naf takes a positive integer up to a maximum value of 2^256 - 1 and returns -// its non-adjacent form (NAF), which is a unique signed-digit representation -// such that no two consecutive digits are nonzero. See the documentation for -// the returned type for details on how the representation is encoded -// efficiently and how to interpret it -// -// NAF is useful in that it has the fewest nonzero digits of any signed digit -// representation, only 1/3rd of its digits are nonzero on average, and at least -// half of the digits will be 0. -// -// The aforementioned properties are particularly beneficial for optimizing -// elliptic curve point multiplication because they effectively minimize the -// number of required point additions in exchange for needing to perform a mix -// of fewer point additions and subtractions and possibly one additional point -// doubling. This is an excellent tradeoff because subtraction of points has -// the same computational complexity as addition of points and point doubling is -// faster than both. -func naf(k []byte) nafScalar { - // Strip leading zero bytes. - for len(k) > 0 && k[0] == 0x00 { - k = k[1:] - } - - // The non-adjacent form (NAF) of a positive integer k is an expression - // k = ∑_(i=0, l-1) k_i * 2^i where k_i ∈ {0,±1}, k_(l-1) != 0, and no two - // consecutive digits k_i are nonzero. - // - // The traditional method of computing the NAF of a positive integer is - // given by algorithm 3.30 in [GECC]. It consists of repeatedly dividing k - // by 2 and choosing the remainder so that the quotient (k−r)/2 is even - // which ensures the next NAF digit is 0. This requires log_2(k) steps. - // - // However, in [BRID], Prodinger notes that a closed form expression for the - // NAF representation is the bitwise difference 3k/2 - k/2. This is more - // efficient as it can be computed in O(1) versus the O(log(n)) of the - // traditional approach. - // - // The following code makes use of that formula to compute the NAF more - // efficiently. - // - // To understand the logic here, observe that the only way the NAF has a - // nonzero digit at a given bit is when either 3k/2 or k/2 has a bit set in - // that position, but not both. In other words, the result of a bitwise - // xor. This can be seen simply by considering that when the bits are the - // same, the subtraction is either 0-0 or 1-1, both of which are 0. - // - // Further, observe that the "+1" digits in the result are contributed by - // 3k/2 while the "-1" digits are from k/2. So, they can be determined by - // taking the bitwise and of each respective value with the result of the - // xor which identifies which bits are nonzero. - // - // Using that information, this loops backwards from the least significant - // byte to the most significant byte while performing the aforementioned - // calculations by propagating the potential carry and high order bit from - // the next word during the right shift. - kLen := len(k) - var result nafScalar - var carry uint8 - for byteNum := kLen - 1; byteNum >= 0; byteNum-- { - // Calculate k/2. Notice the carry from the previous word is added and - // the low order bit from the next word is shifted in accordingly. - kc := uint16(k[byteNum]) + uint16(carry) - var nextWord uint8 - if byteNum > 0 { - nextWord = k[byteNum-1] - } - halfK := kc>>1 | uint16(nextWord<<7) - - // Calculate 3k/2 and determine the non-zero digits in the result. - threeHalfK := kc + halfK - nonZeroResultDigits := threeHalfK ^ halfK - - // Determine the signed digits {0, ±1}. - result.pos[byteNum+1] = uint8(threeHalfK & nonZeroResultDigits) - result.neg[byteNum+1] = uint8(halfK & nonZeroResultDigits) - - // Propagate the potential carry from the 3k/2 calculation. - carry = uint8(threeHalfK >> 8) - } - result.pos[0] = carry - - // Set the starting and ending positions within the fixed size arrays to - // identify the bytes that are actually used. This is important since the - // encoding is big endian and thus trailing zero bytes changes its value. - result.start = 1 - carry - result.end = uint8(kLen + 1) - return result -} - // wNAFPrecompTableSize is the size of the wNAF precomputed point table. const wNAFPrecompTableSize = 1 << (wNAFWidth - 1) diff --git a/dcrec/secp256k1/curve_bench_test.go b/dcrec/secp256k1/curve_bench_test.go index 5bcd50e171..7d8bf3089a 100644 --- a/dcrec/secp256k1/curve_bench_test.go +++ b/dcrec/secp256k1/curve_bench_test.go @@ -216,19 +216,6 @@ func BenchmarkScalarMultNonConst(b *testing.B) { } } -// BenchmarkNAF benchmarks conversion of a positive integer into its -// non-adjacent form representation. -func BenchmarkNAF(b *testing.B) { - k := fromHex("d74bf844b0862475103d96a611cf2d898447e288d34b360bc885cb8ce7c00575") - kBytes := k.Bytes() - - b.ReportAllocs() - b.ResetTimer() - for i := 0; i < b.N; i++ { - naf(kBytes) - } -} - // BenchmarkWNAF benchmarks conversion of a non-negative integer into its // width-w windowed non-adjacent form representation, where w is [wNAFWidth]. func BenchmarkWNAF(b *testing.B) { diff --git a/dcrec/secp256k1/curve_test.go b/dcrec/secp256k1/curve_test.go index a92a6d1ca0..3f32a39a4e 100644 --- a/dcrec/secp256k1/curve_test.go +++ b/dcrec/secp256k1/curve_test.go @@ -634,125 +634,6 @@ func TestDoubleJacobian(t *testing.T) { } } -// checkNAFEncoding returns an error if the provided positive and negative -// portions of an overall NAF encoding do not adhere to the requirements or they -// do not sum back to the provided original value. -func checkNAFEncoding(pos, neg []byte, origValue *big.Int) error { - // NAF must not have a leading zero byte and the number of negative - // bytes must not exceed the positive portion. - if len(pos) > 0 && pos[0] == 0 { - return fmt.Errorf("positive has leading zero -- got %x", pos) - } - if len(neg) > len(pos) { - return fmt.Errorf("negative has len %d > pos len %d", len(neg), - len(pos)) - } - - // Ensure the result doesn't have any adjacent non-zero digits. - gotPos := new(big.Int).SetBytes(pos) - gotNeg := new(big.Int).SetBytes(neg) - posOrNeg := new(big.Int).Or(gotPos, gotNeg) - prevBit := posOrNeg.Bit(0) - for bit := 1; bit < posOrNeg.BitLen(); bit++ { - thisBit := posOrNeg.Bit(bit) - if prevBit == 1 && thisBit == 1 { - return fmt.Errorf("adjacent non-zero digits found at bit pos %d", - bit-1) - } - prevBit = thisBit - } - - // Ensure the resulting positive and negative portions of the overall - // NAF representation sum back to the original value. - gotValue := new(big.Int).Sub(gotPos, gotNeg) - if origValue.Cmp(gotValue) != 0 { - return fmt.Errorf("pos-neg is not original value: got %x, want %x", - gotValue, origValue) - } - - return nil -} - -// TestNAF ensures encoding various edge cases and values to non-adjacent form -// produces valid results. -func TestNAF(t *testing.T) { - tests := []struct { - name string // test description - in string // hex encoded test value - }{{ - name: "empty is zero", - in: "", - }, { - name: "zero", - in: "00", - }, { - name: "just before first carry", - in: "aa", - }, { - name: "first carry", - in: "ab", - }, { - name: "leading zeroes", - in: "002f20569b90697ad471c1be6107814f53f47446be298a3a2a6b686b97d35cf9", - }, { - name: "257 bits when NAF encoded", - in: "c000000000000000000000000000000000000000000000000000000000000001", - }, { - name: "32-byte scalar", - in: "6df2b5d30854069ccdec40ae022f5c948936324a4e9ebed8eb82cfd5a6b6d766", - }, { - name: "first term of balanced length-two representation #1", - in: "b776e53fb55f6b006a270d42d64ec2b1", - }, { - name: "second term balanced length-two representation #1", - in: "d6cc32c857f1174b604eefc544f0c7f7", - }, { - name: "first term of balanced length-two representation #2", - in: "45c53aa1bb56fcd68c011e2dad6758e4", - }, { - name: "second term of balanced length-two representation #2", - in: "a2e79d200f27f2360fba57619936159b", - }} - - for _, test := range tests { - // Ensure the resulting positive and negative portions of the overall - // NAF representation adhere to the requirements of NAF encoding and - // they sum back to the original value. - result := naf(hexToBytes(test.in)) - pos, neg := result.Pos(), result.Neg() - if err := checkNAFEncoding(pos, neg, fromHex(test.in)); err != nil { - t.Errorf("%q: %v", test.name, err) - } - } -} - -// TestNAFRandom ensures that encoding randomly-generated values to non-adjacent -// form produces valid results. -func TestNAFRandom(t *testing.T) { - // Use a unique random seed each test instance and log it if the tests fail. - seed := time.Now().Unix() - rng := mrand.New(mrand.NewSource(seed)) - defer func(t *testing.T, seed int64) { - if t.Failed() { - t.Logf("random seed: %d", seed) - } - }(t, seed) - - for i := 0; i < 100; i++ { - // Ensure the resulting positive and negative portions of the overall - // NAF representation adhere to the requirements of NAF encoding and - // they sum back to the original value. - bigIntVal, modNVal := randIntAndModNScalar(t, rng) - valBytes := modNVal.Bytes() - result := naf(valBytes[:]) - pos, neg := result.Pos(), result.Neg() - if err := checkNAFEncoding(pos, neg, bigIntVal); err != nil { - t.Fatalf("encoding err: %v\nin: %x\npos: %x\nneg: %x", err, - bigIntVal, pos, neg) - } - } -} - // decodeWNAFDigit converts the provided wNAF code to the signed digit it // represents. func decodeWNAFDigit(code uint8) int8 { From 872776ecbbd3a6085dedda3358e3d048e26083f8 Mon Sep 17 00:00:00 2001 From: Dave Collins Date: Sat, 25 Jul 2026 04:01:51 -0500 Subject: [PATCH 12/17] secp256k1: Move constant time helpers to arith. The package currently defines several generic constant time helpers in the main package despite them not being specific to any particular type or algorithm. They are instead arithmetic primitives that are shared by multiple implementations. Notably, the helpers used by the field implementations that will be split into separate packages in future commits to prepare for architecture-specific backend selection. In order to smooth that transition, this moves the helpers into a new internal arith package and updates all call sites to use them. No functional changes are intended. The move is purely organizational and preserves the existing constant time semantics. --- dcrec/secp256k1/common.go | 56 ------------- dcrec/secp256k1/field.go | 82 ++++++++++--------- dcrec/secp256k1/field64.go | 32 ++++---- dcrec/secp256k1/internal/arith/consttime.go | 56 +++++++++++++ .../arith/consttime_test.go} | 4 +- dcrec/secp256k1/modnscalar.go | 33 ++++---- 6 files changed, 135 insertions(+), 128 deletions(-) delete mode 100644 dcrec/secp256k1/common.go create mode 100644 dcrec/secp256k1/internal/arith/consttime.go rename dcrec/secp256k1/{common_test.go => internal/arith/consttime_test.go} (93%) diff --git a/dcrec/secp256k1/common.go b/dcrec/secp256k1/common.go deleted file mode 100644 index 9ff782bf40..0000000000 --- a/dcrec/secp256k1/common.go +++ /dev/null @@ -1,56 +0,0 @@ -// Copyright (c) 2020-2026 The Decred developers -// Use of this source code is governed by an ISC -// license that can be found in the LICENSE file. - -package secp256k1 - -// constantTimeEq returns 1 if a == b or 0 otherwise in constant time. -func constantTimeEq(a, b uint32) uint32 { - return uint32((uint64(a^b) - 1) >> 63) -} - -// constantTimeEq64 returns 1 if a == b or 0 otherwise in constant time. -func constantTimeEq64(a, b uint64) uint32 { - t := a ^ b - return uint32(((t | -t) >> 63) ^ 1) -} - -// constantTimeNotEq64 returns 1 if a != b or 0 otherwise in constant time. -func constantTimeNotEq64(a, b uint64) uint32 { - t := a ^ b - return uint32(((t | -t) >> 63)) -} - -// constantTimeLess returns 1 if a < b or 0 otherwise in constant time. -func constantTimeLess(a, b uint32) uint32 { - return uint32((uint64(a) - uint64(b)) >> 63) -} - -// constantTimeLessOrEq returns 1 if a <= b or 0 otherwise in constant time. -func constantTimeLessOrEq(a, b uint32) uint32 { - return uint32((uint64(a) - uint64(b) - 1) >> 63) -} - -// constantTimeGreater returns 1 if a > b or 0 otherwise in constant time. -func constantTimeGreater(a, b uint32) uint32 { - return constantTimeLess(b, a) -} - -// constantTimeGreaterOrEq returns 1 if a >= b or 0 otherwise in constant time. -func constantTimeGreaterOrEq(a, b uint32) uint32 { - return constantTimeLessOrEq(b, a) -} - -// constantTimeMin returns min(a,b) in constant time. -func constantTimeMin(a, b uint32) uint32 { - return b ^ ((a ^ b) & -constantTimeLess(a, b)) -} - -// constantTimeSelect64 returns a if cond == 1 or b if cond == 0 in constant -// time. -// -// WARNING: The behavior is undefined if cond is anything other than 0 or 1. -func constantTimeSelect64(cond, a, b uint64) uint64 { - mask := -cond - return b ^ (a^b)&mask -} diff --git a/dcrec/secp256k1/field.go b/dcrec/secp256k1/field.go index f633b46dfc..638e01e503 100644 --- a/dcrec/secp256k1/field.go +++ b/dcrec/secp256k1/field.go @@ -6,6 +6,12 @@ package secp256k1 +import ( + "encoding/hex" + + "github.com/decred/dcrd/dcrec/secp256k1/v4/internal/arith" +) + // References: // [HAC]: Handbook of Applied Cryptography Menezes, van Oorschot, Vanstone. // https://cacr.uwaterloo.ca/hac/ @@ -59,10 +65,6 @@ package secp256k1 // performance crypto, this type does not perform any validation where it // ordinarily would. See the documentation for [FieldVal] for more details. -import ( - "encoding/hex" -) - // Constants used to make the code more readable. const ( twoBitsMask = 0x3 @@ -292,17 +294,17 @@ func (f *FieldVal) SetBytes(b *[32]byte) uint32 { // // Thus, there is no need to test if the upper words of the field value // exceeds them, hence, only equality is checked for them. - highWordsEq := constantTimeEq(f.n[9], fieldPrimeWordNine) - highWordsEq &= constantTimeEq(f.n[8], fieldPrimeWordEight) - highWordsEq &= constantTimeEq(f.n[7], fieldPrimeWordSeven) - highWordsEq &= constantTimeEq(f.n[6], fieldPrimeWordSix) - highWordsEq &= constantTimeEq(f.n[5], fieldPrimeWordFive) - highWordsEq &= constantTimeEq(f.n[4], fieldPrimeWordFour) - highWordsEq &= constantTimeEq(f.n[3], fieldPrimeWordThree) - highWordsEq &= constantTimeEq(f.n[2], fieldPrimeWordTwo) - overflow := highWordsEq & constantTimeGreater(f.n[1], fieldPrimeWordOne) - highWordsEq &= constantTimeEq(f.n[1], fieldPrimeWordOne) - overflow |= highWordsEq & constantTimeGreaterOrEq(f.n[0], fieldPrimeWordZero) + highWordsEq := arith.ConstantTimeEq(f.n[9], fieldPrimeWordNine) + highWordsEq &= arith.ConstantTimeEq(f.n[8], fieldPrimeWordEight) + highWordsEq &= arith.ConstantTimeEq(f.n[7], fieldPrimeWordSeven) + highWordsEq &= arith.ConstantTimeEq(f.n[6], fieldPrimeWordSix) + highWordsEq &= arith.ConstantTimeEq(f.n[5], fieldPrimeWordFive) + highWordsEq &= arith.ConstantTimeEq(f.n[4], fieldPrimeWordFour) + highWordsEq &= arith.ConstantTimeEq(f.n[3], fieldPrimeWordThree) + highWordsEq &= arith.ConstantTimeEq(f.n[2], fieldPrimeWordTwo) + overflow := highWordsEq & arith.ConstantTimeGreater(f.n[1], fieldPrimeWordOne) + highWordsEq &= arith.ConstantTimeEq(f.n[1], fieldPrimeWordOne) + overflow |= highWordsEq & arith.ConstantTimeGreaterOrEq(f.n[0], fieldPrimeWordZero) return overflow } @@ -325,7 +327,7 @@ func (f *FieldVal) SetBytes(b *[32]byte) uint32 { // Output Max Magnitude: 1 func (f *FieldVal) SetByteSlice(b []byte) bool { var b32 [32]byte - b = b[:constantTimeMin(uint32(len(b)), 32)] + b = b[:arith.ConstantTimeMin(uint32(len(b)), 32)] copy(b32[:], b32[:32-len(b)]) copy(b32[32-len(b):], b) result := f.SetBytes(&b32) @@ -402,9 +404,9 @@ func (f *FieldVal) Normalize() *FieldVal { // // Also note that 'm' will be zero when neither of the aforementioned // conditions are true and the value will not be changed when 'm' is zero. - m = constantTimeEq(t9, fieldMSBMask) - m &= constantTimeEq(t8&t7&t6&t5&t4&t3&t2, fieldBaseMask) - m &= constantTimeGreater(t1+64+((t0+977)>>fieldBase), fieldBaseMask) + m = arith.ConstantTimeEq(t9, fieldMSBMask) + m &= arith.ConstantTimeEq(t8&t7&t6&t5&t4&t3&t2, fieldBaseMask) + m &= arith.ConstantTimeGreater(t1+64+((t0+977)>>fieldBase), fieldBaseMask) m |= t9 >> fieldMSBBits t0 += m * 977 t1 = (t0 >> fieldBase) + t1 + (m << 6) @@ -541,7 +543,7 @@ func (f *FieldVal) IsZeroBit() uint32 { bits := f.n[0] | f.n[1] | f.n[2] | f.n[3] | f.n[4] | f.n[5] | f.n[6] | f.n[7] | f.n[8] | f.n[9] - return constantTimeEq(bits, 0) + return arith.ConstantTimeEq(bits, 0) } // IsZero returns whether or not the field value is equal to zero in constant @@ -575,7 +577,7 @@ func (f *FieldVal) IsOneBit() uint32 { bits := (f.n[0] ^ 1) | f.n[1] | f.n[2] | f.n[3] | f.n[4] | f.n[5] | f.n[6] | f.n[7] | f.n[8] | f.n[9] - return constantTimeEq(bits, 0) + return arith.ConstantTimeEq(bits, 0) } // IsOne returns whether or not the field value is equal to one in constant @@ -1750,25 +1752,25 @@ func (f *FieldVal) IsGtOrEqPrimeMinusOrder() bool { // The intuition here is that the value is greater than field prime minus // the group order if one of the higher individual words is greater than the // corresponding word and all higher words in the value are equal. - result := constantTimeGreater(f.n[9], pMinusNWordNine) - highWordsEqual := constantTimeEq(f.n[9], pMinusNWordNine) - result |= highWordsEqual & constantTimeGreater(f.n[8], pMinusNWordEight) - highWordsEqual &= constantTimeEq(f.n[8], pMinusNWordEight) - result |= highWordsEqual & constantTimeGreater(f.n[7], pMinusNWordSeven) - highWordsEqual &= constantTimeEq(f.n[7], pMinusNWordSeven) - result |= highWordsEqual & constantTimeGreater(f.n[6], pMinusNWordSix) - highWordsEqual &= constantTimeEq(f.n[6], pMinusNWordSix) - result |= highWordsEqual & constantTimeGreater(f.n[5], pMinusNWordFive) - highWordsEqual &= constantTimeEq(f.n[5], pMinusNWordFive) - result |= highWordsEqual & constantTimeGreater(f.n[4], pMinusNWordFour) - highWordsEqual &= constantTimeEq(f.n[4], pMinusNWordFour) - result |= highWordsEqual & constantTimeGreater(f.n[3], pMinusNWordThree) - highWordsEqual &= constantTimeEq(f.n[3], pMinusNWordThree) - result |= highWordsEqual & constantTimeGreater(f.n[2], pMinusNWordTwo) - highWordsEqual &= constantTimeEq(f.n[2], pMinusNWordTwo) - result |= highWordsEqual & constantTimeGreater(f.n[1], pMinusNWordOne) - highWordsEqual &= constantTimeEq(f.n[1], pMinusNWordOne) - result |= highWordsEqual & constantTimeGreaterOrEq(f.n[0], pMinusNWordZero) + result := arith.ConstantTimeGreater(f.n[9], pMinusNWordNine) + highWordsEqual := arith.ConstantTimeEq(f.n[9], pMinusNWordNine) + result |= highWordsEqual & arith.ConstantTimeGreater(f.n[8], pMinusNWordEight) + highWordsEqual &= arith.ConstantTimeEq(f.n[8], pMinusNWordEight) + result |= highWordsEqual & arith.ConstantTimeGreater(f.n[7], pMinusNWordSeven) + highWordsEqual &= arith.ConstantTimeEq(f.n[7], pMinusNWordSeven) + result |= highWordsEqual & arith.ConstantTimeGreater(f.n[6], pMinusNWordSix) + highWordsEqual &= arith.ConstantTimeEq(f.n[6], pMinusNWordSix) + result |= highWordsEqual & arith.ConstantTimeGreater(f.n[5], pMinusNWordFive) + highWordsEqual &= arith.ConstantTimeEq(f.n[5], pMinusNWordFive) + result |= highWordsEqual & arith.ConstantTimeGreater(f.n[4], pMinusNWordFour) + highWordsEqual &= arith.ConstantTimeEq(f.n[4], pMinusNWordFour) + result |= highWordsEqual & arith.ConstantTimeGreater(f.n[3], pMinusNWordThree) + highWordsEqual &= arith.ConstantTimeEq(f.n[3], pMinusNWordThree) + result |= highWordsEqual & arith.ConstantTimeGreater(f.n[2], pMinusNWordTwo) + highWordsEqual &= arith.ConstantTimeEq(f.n[2], pMinusNWordTwo) + result |= highWordsEqual & arith.ConstantTimeGreater(f.n[1], pMinusNWordOne) + highWordsEqual &= arith.ConstantTimeEq(f.n[1], pMinusNWordOne) + result |= highWordsEqual & arith.ConstantTimeGreaterOrEq(f.n[0], pMinusNWordZero) return result != 0 } diff --git a/dcrec/secp256k1/field64.go b/dcrec/secp256k1/field64.go index 5be951eb0b..30afab1e2b 100644 --- a/dcrec/secp256k1/field64.go +++ b/dcrec/secp256k1/field64.go @@ -8,6 +8,8 @@ import ( "encoding/binary" "encoding/hex" "math/bits" + + "github.com/decred/dcrd/dcrec/secp256k1/v4/internal/arith" ) // References: @@ -137,10 +139,10 @@ func (f *FieldVal64) SetBytes(b *[32]byte) uint32 { // Constant-time select. // // Set f = f when f < p (aka borrow is set). Otherwise f = s = f - p. - f.n[0] = constantTimeSelect64(borrow, f.n[0], s0) - f.n[1] = constantTimeSelect64(borrow, f.n[1], s1) - f.n[2] = constantTimeSelect64(borrow, f.n[2], s2) - f.n[3] = constantTimeSelect64(borrow, f.n[3], s3) + f.n[0] = arith.ConstantTimeSelect64(borrow, f.n[0], s0) + f.n[1] = arith.ConstantTimeSelect64(borrow, f.n[1], s1) + f.n[2] = arith.ConstantTimeSelect64(borrow, f.n[2], s2) + f.n[3] = arith.ConstantTimeSelect64(borrow, f.n[3], s3) return uint32(1 - borrow) } @@ -158,7 +160,7 @@ func (f *FieldVal64) SetBytes(b *[32]byte) uint32 { // overflow behavior. func (f *FieldVal64) SetByteSlice(b []byte) bool { var b32 [32]byte - b = b[:constantTimeMin(uint32(len(b)), 32)] + b = b[:arith.ConstantTimeMin(uint32(len(b)), 32)] copy(b32[:], b32[:32-len(b)]) copy(b32[32-len(b):], b) result := f.SetBytes(&b32) @@ -220,7 +222,7 @@ func (f *FieldVal64) Bytes() *[32]byte { // operations require a numeric value. See IsZero for the version that returns // a bool. func (f *FieldVal64) IsZeroBit() uint32 { - return constantTimeEq64(f.n[0]|f.n[1]|f.n[2]|f.n[3], 0) + return arith.ConstantTimeEq64(f.n[0]|f.n[1]|f.n[2]|f.n[3], 0) } // IsZero returns whether or not the field value is equal to zero in constant @@ -240,7 +242,7 @@ func (f *FieldVal64) IsOneBit() uint32 { // The value can only be one if the single lowest significant bit is set in // the first word and no other bits are set in any of the other words. // This is a constant time implementation. - return constantTimeEq64((f.n[0]^1)|f.n[1]|f.n[2]|f.n[3], 0) + return arith.ConstantTimeEq64((f.n[0]^1)|f.n[1]|f.n[2]|f.n[3], 0) } // IsOne returns whether or not the field value is equal to one in constant @@ -389,10 +391,10 @@ func (f *FieldVal64) Add2(a, b *FieldVal64) *FieldVal64 { // Set f = t = a+b only when there was no overflow and t < p (borrow set). // Otherwise f = s = a+b - p. cond := (1 - overflow) & borrow - f.n[0] = constantTimeSelect64(cond, t0, s0) - f.n[1] = constantTimeSelect64(cond, t1, s1) - f.n[2] = constantTimeSelect64(cond, t2, s2) - f.n[3] = constantTimeSelect64(cond, t3, s3) + f.n[0] = arith.ConstantTimeSelect64(cond, t0, s0) + f.n[1] = arith.ConstantTimeSelect64(cond, t1, s1) + f.n[2] = arith.ConstantTimeSelect64(cond, t2, s2) + f.n[3] = arith.ConstantTimeSelect64(cond, t3, s3) return f } @@ -943,10 +945,10 @@ func field64Reduce512(r *[4]uint64, x *[8]uint64) { s2, borrow = bits.Sub64(t2, field64Prime2, borrow) s3, borrow = bits.Sub64(t3, field64Prime3, borrow) _, borrow = bits.Sub64(t4, 0, borrow) - r[0] = constantTimeSelect64(borrow, t0, s0) - r[1] = constantTimeSelect64(borrow, t1, s1) - r[2] = constantTimeSelect64(borrow, t2, s2) - r[3] = constantTimeSelect64(borrow, t3, s3) + r[0] = arith.ConstantTimeSelect64(borrow, t0, s0) + r[1] = arith.ConstantTimeSelect64(borrow, t1, s1) + r[2] = arith.ConstantTimeSelect64(borrow, t2, s2) + r[3] = arith.ConstantTimeSelect64(borrow, t3, s3) } // Inverse finds the modular multiplicative inverse of the field value in diff --git a/dcrec/secp256k1/internal/arith/consttime.go b/dcrec/secp256k1/internal/arith/consttime.go new file mode 100644 index 0000000000..b798e824d5 --- /dev/null +++ b/dcrec/secp256k1/internal/arith/consttime.go @@ -0,0 +1,56 @@ +// Copyright (c) 2020-2026 The Decred developers +// Use of this source code is governed by an ISC +// license that can be found in the LICENSE file. + +package arith + +// ConstantTimeEq returns 1 if a == b or 0 otherwise in constant time. +func ConstantTimeEq(a, b uint32) uint32 { + return uint32((uint64(a^b) - 1) >> 63) +} + +// ConstantTimeEq64 returns 1 if a == b or 0 otherwise in constant time. +func ConstantTimeEq64(a, b uint64) uint32 { + t := a ^ b + return uint32(((t | -t) >> 63) ^ 1) +} + +// ConstantTimeNotEq64 returns 1 if a != b or 0 otherwise in constant time. +func ConstantTimeNotEq64(a, b uint64) uint32 { + t := a ^ b + return uint32((t | -t) >> 63) +} + +// ConstantTimeLess returns 1 if a < b or 0 otherwise in constant time. +func ConstantTimeLess(a, b uint32) uint32 { + return uint32((uint64(a) - uint64(b)) >> 63) +} + +// ConstantTimeLessOrEq returns 1 if a <= b or 0 otherwise in constant time. +func ConstantTimeLessOrEq(a, b uint32) uint32 { + return uint32((uint64(a) - uint64(b) - 1) >> 63) +} + +// ConstantTimeGreater returns 1 if a > b or 0 otherwise in constant time. +func ConstantTimeGreater(a, b uint32) uint32 { + return ConstantTimeLess(b, a) +} + +// ConstantTimeGreaterOrEq returns 1 if a >= b or 0 otherwise in constant time. +func ConstantTimeGreaterOrEq(a, b uint32) uint32 { + return ConstantTimeLessOrEq(b, a) +} + +// ConstantTimeMin returns min(a,b) in constant time. +func ConstantTimeMin(a, b uint32) uint32 { + return b ^ ((a ^ b) & -ConstantTimeLess(a, b)) +} + +// ConstantTimeSelect64 returns a when cond == 1 or b when cond == 0 in constant +// time. +// +// WARNING: The behavior is undefined if cond is anything other than 0 or 1. +func ConstantTimeSelect64(cond, a, b uint64) uint64 { + mask := -cond + return b ^ (a^b)&mask +} diff --git a/dcrec/secp256k1/common_test.go b/dcrec/secp256k1/internal/arith/consttime_test.go similarity index 93% rename from dcrec/secp256k1/common_test.go rename to dcrec/secp256k1/internal/arith/consttime_test.go index cde06019ff..08a817ea8c 100644 --- a/dcrec/secp256k1/common_test.go +++ b/dcrec/secp256k1/internal/arith/consttime_test.go @@ -3,7 +3,7 @@ // Use of this source code is governed by an ISC // license that can be found in the LICENSE file. -package secp256k1 +package arith import "testing" @@ -28,7 +28,7 @@ func TestConstantTimeSelect64(t *testing.T) { } for _, test := range tests { - got := constantTimeSelect64(test.cond, test.a, test.b) + got := ConstantTimeSelect64(test.cond, test.a, test.b) if got != test.want { t.Errorf("%q: unexpected result -- got %d, want %d", test.name, got, test.want) diff --git a/dcrec/secp256k1/modnscalar.go b/dcrec/secp256k1/modnscalar.go index 0eaf792c17..26d13c5dee 100644 --- a/dcrec/secp256k1/modnscalar.go +++ b/dcrec/secp256k1/modnscalar.go @@ -10,6 +10,8 @@ import ( "math/big" "math/bits" "sync" + + "github.com/decred/dcrd/dcrec/secp256k1/v4/internal/arith" ) // References: @@ -175,7 +177,7 @@ func (s *ModNScalar) Zero() { // operations require a numeric value. See [ModNScalar.IsZero] for the version // that returns a bool. func (s *ModNScalar) IsZeroBit() uint32 { - return constantTimeEq64(s.n[0]|s.n[1]|s.n[2]|s.n[3], 0) + return arith.ConstantTimeEq64(s.n[0]|s.n[1]|s.n[2]|s.n[3], 0) } // IsZero returns whether or not the scalar is equal to zero in constant time. @@ -234,10 +236,10 @@ func (s *ModNScalar) SetBytes(b *[32]byte) uint32 { // Constant-time select. // // Set s = s when s < N (aka borrow is set). Otherwise s = t = s - N. - s.n[0] = constantTimeSelect64(borrow, s.n[0], t0) - s.n[1] = constantTimeSelect64(borrow, s.n[1], t1) - s.n[2] = constantTimeSelect64(borrow, s.n[2], t2) - s.n[3] = constantTimeSelect64(borrow, s.n[3], t3) + s.n[0] = arith.ConstantTimeSelect64(borrow, s.n[0], t0) + s.n[1] = arith.ConstantTimeSelect64(borrow, s.n[1], t1) + s.n[2] = arith.ConstantTimeSelect64(borrow, s.n[2], t2) + s.n[3] = arith.ConstantTimeSelect64(borrow, s.n[3], t3) return uint32(1 - borrow) } @@ -261,7 +263,7 @@ func (s *ModNScalar) SetByteSlice(b []byte) bool { // Always copy a total of 32 bytes regardless of the input length to avoid // introducing data-dependent timing. var b32 [32]byte - b = b[:constantTimeMin(uint32(len(b)), 32)] + b = b[:arith.ConstantTimeMin(uint32(len(b)), 32)] copy(b32[:], b32[:32-len(b)]) copy(b32[32-len(b):], b) result := s.SetBytes(&b32) @@ -368,10 +370,10 @@ func (s *ModNScalar) Add2(a, b *ModNScalar) *ModNScalar { // Set s = t = a+b only when there was no overflow and t < N (borrow set). // Otherwise s = u = a+b - N. cond := (1 - overflow) & borrow - s.n[0] = constantTimeSelect64(cond, t0, u0) - s.n[1] = constantTimeSelect64(cond, t1, u1) - s.n[2] = constantTimeSelect64(cond, t2, u2) - s.n[3] = constantTimeSelect64(cond, t3, u3) + s.n[0] = arith.ConstantTimeSelect64(cond, t0, u0) + s.n[1] = arith.ConstantTimeSelect64(cond, t1, u1) + s.n[2] = arith.ConstantTimeSelect64(cond, t2, u2) + s.n[3] = arith.ConstantTimeSelect64(cond, t3, u3) return s } @@ -737,10 +739,10 @@ func scalar64Reduce512(r *[4]uint64, x *[8]uint64) { // // Set r = t only when t < N (borrow set). // Otherwise r = s = t - N. - r[0] = constantTimeSelect64(borrow, t0, s0) - r[1] = constantTimeSelect64(borrow, t1, s1) - r[2] = constantTimeSelect64(borrow, t2, s2) - r[3] = constantTimeSelect64(borrow, t3, s3) + r[0] = arith.ConstantTimeSelect64(borrow, t0, s0) + r[1] = arith.ConstantTimeSelect64(borrow, t1, s1) + r[2] = arith.ConstantTimeSelect64(borrow, t2, s2) + r[3] = arith.ConstantTimeSelect64(borrow, t3, s3) } // Mul2 multiplies the passed two scalars together modulo the group order in @@ -805,7 +807,8 @@ func (s *ModNScalar) NegateVal(val *ModNScalar) *ModNScalar { // being about 3-4% slower on 64-bit hardware. // // Determine mask first to allow aliasing. - mask := -uint64(constantTimeNotEq64(val.n[0]|val.n[1]|val.n[2]|val.n[3], 0)) + combinedLimbs := val.n[0] | val.n[1] | val.n[2] | val.n[3] + mask := -uint64(arith.ConstantTimeNotEq64(combinedLimbs, 0)) // Unconditionally subtract the scalar from the group order. // From b8fa84b9c986fb86c0bac2dbfd6b0e483dcb2241 Mon Sep 17 00:00:00 2001 From: Dave Collins Date: Sat, 25 Jul 2026 06:04:43 -0500 Subject: [PATCH 13/17] secp256k1: Add constant time helper tests. This adds dedicated tests for the constant time helpers in the new internal arith package. While the helpers are already well exercised by the upper layers that use them, testing the primitives directly helps ensure their correctness in isolation. --- .../internal/arith/consttime_test.go | 280 +++++++++++++++++- 1 file changed, 267 insertions(+), 13 deletions(-) diff --git a/dcrec/secp256k1/internal/arith/consttime_test.go b/dcrec/secp256k1/internal/arith/consttime_test.go index 08a817ea8c..35941ee8bc 100644 --- a/dcrec/secp256k1/internal/arith/consttime_test.go +++ b/dcrec/secp256k1/internal/arith/consttime_test.go @@ -7,31 +7,285 @@ package arith import "testing" -// TestConstantTimeSelect64 ensures that the 64-bit constant time selection -// function works as expected in terms of behavior. +// TestConstantTimeEq ensures [ConstantTimeEq] returns the expected results. +func TestConstantTimeEq(t *testing.T) { + tests := []struct { + a uint32 // first value + b uint32 // second value + want uint32 // expected result (0 or 1) + }{ + {0, 0, 1}, + {1, 1, 1}, + {1<<32 - 1, 1<<32 - 1, 1}, + {0x12345678, 0x12345678, 1}, + {0, 1, 0}, + {1, 0, 0}, + {1<<32 - 1, 0, 0}, + {0x12345678, 0x87654321, 0}, + {0, 1 << 31, 0}, + {1 << 31, 1 << 31, 1}, + {1, 2, 0}, + {1<<32 - 2, 1<<32 - 1, 0}, + } + + for _, test := range tests { + got := ConstantTimeEq(test.a, test.b) + if got != test.want { + t.Errorf("%08x == %08x: got %d, want %d", test.a, test.b, got, + test.want) + } + } +} + +// TestConstantTimeEq64 ensures [ConstantTimeEq64] returns the expected results. +func TestConstantTimeEq64(t *testing.T) { + tests := []struct { + a uint64 // first value + b uint64 // second value + want uint32 // expected result (0 or 1) + }{ + {0, 0, 1}, + {1, 1, 1}, + {1<<64 - 1, 1<<64 - 1, 1}, + {0x123456789abcdef0, 0x123456789abcdef0, 1}, + {0, 1, 0}, + {1, 0, 0}, + {1<<64 - 1, 0, 0}, + {0x123456789abcdef0, 0xfedcba9876543210, 0}, + {0, 1 << 63, 0}, + {1 << 63, 1 << 63, 1}, + {1, 2, 0}, + {1<<64 - 2, 1<<64 - 1, 0}, + } + + for _, test := range tests { + got := ConstantTimeEq64(test.a, test.b) + if got != test.want { + t.Errorf("%016x == %016x: got %d, want %d", test.a, test.b, got, + test.want) + } + } +} + +// TestConstantTimeNotEq64 ensures [ConstantTimeNotEq64] returns the expected +// results. +func TestConstantTimeNotEq64(t *testing.T) { + tests := []struct { + a uint64 // first value + b uint64 // second value + want uint32 // expected result + }{ + {0, 0, 0}, + {1, 1, 0}, + {1<<64 - 1, 1<<64 - 1, 0}, + {0x123456789abcdef0, 0x123456789abcdef0, 0}, + {0, 1, 1}, + {1, 0, 1}, + {1<<64 - 1, 0, 1}, + {0x123456789abcdef0, 0xfedcba9876543210, 1}, + {0, 1 << 63, 1}, + {1 << 63, 1 << 63, 0}, + {1, 2, 1}, + {1<<64 - 2, 1<<64 - 1, 1}, + } + + for _, test := range tests { + got := ConstantTimeNotEq64(test.a, test.b) + if got != test.want { + t.Errorf("%016x != %016x: got %d, want %d", test.a, test.b, got, + test.want) + } + } +} + +// TestConstantTimeLess ensures [ConstantTimeLess] returns the expected results. +func TestConstantTimeLess(t *testing.T) { + tests := []struct { + a uint32 // first value + b uint32 // second value + want uint32 // expected result (0 or 1) + }{ + {0, 1, 1}, + {1, 2, 1}, + {0, 1<<32 - 1, 1}, + {1<<31 - 1, 1 << 31, 1}, + {1<<32 - 2, 1<<32 - 1, 1}, + {0x12345678, 0x87654321, 1}, + {1, 0, 0}, + {1, 1, 0}, + {1<<32 - 1, 0, 0}, + {1<<32 - 1, 1<<32 - 1, 0}, + {0x87654321, 0x12345678, 0}, + {0, 0, 0}, + {1 << 31, 1<<32 - 1, 1}, + } + + for _, test := range tests { + got := ConstantTimeLess(test.a, test.b) + if got != test.want { + t.Errorf("%08x < %08x: got %d, want %d", test.a, test.b, got, + test.want) + } + } +} + +// TestConstantTimeLessOrEq ensures [ConstantTimeLessOrEq] returns the expected +// results. +func TestConstantTimeLessOrEq(t *testing.T) { + tests := []struct { + a uint32 // first value + b uint32 // second value + want uint32 // expected result (0 or 1) + }{ + {0, 0, 1}, + {0, 1, 1}, + {1, 2, 1}, + {1, 1, 1}, + {0, 1<<32 - 1, 1}, + {0x12345678, 0x87654321, 1}, + {1<<32 - 1, 1<<32 - 1, 1}, + {1<<32 - 2, 1<<32 - 1, 1}, + {1, 0, 0}, + {2, 1, 0}, + {1<<32 - 1, 0, 0}, + {0x87654321, 0x12345678, 0}, + {1 << 31, 1<<32 - 1, 1}, + {1 << 31, 1 << 31, 1}, + } + + for _, test := range tests { + got := ConstantTimeLessOrEq(test.a, test.b) + if got != test.want { + t.Errorf("%08x <= %08x: got %d, want %d", test.a, test.b, got, + test.want) + } + } +} + +// TestConstantTimeGreater ensures [ConstantTimeGreater] returns the expected +// results. +func TestConstantTimeGreater(t *testing.T) { + tests := []struct { + a uint32 // first value + b uint32 // second value + want uint32 // expected result (0 or 1) + }{ + {1, 0, 1}, + {2, 1, 1}, + {1<<32 - 1, 0, 1}, + {1<<32 - 1, 1<<32 - 2, 1}, + {0x87654321, 0x12345678, 1}, + {0, 1, 0}, + {1, 1, 0}, + {0, 1<<32 - 1, 0}, + {1<<32 - 1, 1<<32 - 1, 0}, + {0x12345678, 0x87654321, 0}, + {0, 0, 0}, + {1<<32 - 1, 1 << 31, 1}, + } + + for _, test := range tests { + got := ConstantTimeGreater(test.a, test.b) + if got != test.want { + t.Errorf("%08x > %08x: got %d, want %d", test.a, test.b, got, + test.want) + } + } +} + +// TestConstantTimeGreaterOrEq ensures [ConstantTimeGreaterOrEq] returns the +// expected results. +func TestConstantTimeGreaterOrEq(t *testing.T) { + tests := []struct { + a uint32 // first value + b uint32 // second value + want uint32 // expected result (0 or 1) + }{ + {0, 0, 1}, + {1, 0, 1}, + {1, 1, 1}, + {2, 1, 1}, + {1<<32 - 1, 0, 1}, + {1<<32 - 1, 1<<32 - 1, 1}, + {1<<32 - 1, 1<<32 - 2, 1}, + {0x87654321, 0x12345678, 1}, + {0, 1, 0}, + {1, 2, 0}, + {0, 1<<32 - 1, 0}, + {0x12345678, 0x87654321, 0}, + {1 << 31, 1 << 31, 1}, + {1 << 31, 1<<32 - 1, 0}, + } + + for _, test := range tests { + got := ConstantTimeGreaterOrEq(test.a, test.b) + if got != test.want { + t.Errorf("%08x >= %08x: got %d, want %d", test.a, test.b, got, + test.want) + } + } +} + +// TestConstantTimeMin ensures [ConstantTimeMin] returns the expected results. +func TestConstantTimeMin(t *testing.T) { + tests := []struct { + a uint32 // first value + b uint32 // second value + want uint32 // expected minimum value + }{ + {0, 1, 0}, + {1, 2, 1}, + {0, 1<<32 - 1, 0}, + {1<<32 - 2, 1<<32 - 1, 1<<32 - 2}, + {0x12345678, 0x87654321, 0x12345678}, + {1, 0, 0}, + {2, 1, 1}, + {1<<32 - 1, 0, 0}, + {1<<32 - 1, 1<<32 - 2, 1<<32 - 2}, + {0x87654321, 0x12345678, 0x12345678}, + {0, 0, 0}, + {1, 1, 1}, + {1<<32 - 1, 1<<32 - 1, 1<<32 - 1}, + {0x12345678, 0x12345678, 0x12345678}, + {1 << 31, 1<<32 - 1, 1 << 31}, + {1<<32 - 1, 1 << 31, 1 << 31}, + } + + for _, test := range tests { + got := ConstantTimeMin(test.a, test.b) + if got != test.want { + t.Errorf("min(%08x, %08x): got %08x, want %08x", test.a, test.b, + got, test.want) + } + } +} + +// TestConstantTimeSelect64 ensures [ConstantTimeSelect64] returns the expected +// results. func TestConstantTimeSelect64(t *testing.T) { tests := []struct { - name string // test description cond uint64 // condition value a uint64 // first value b uint64 // second value want uint64 // expected selected value }{ - {"sel(0,1,2)==2", 0, 1, 2, 2}, - {"sel(1,1,2)==1", 1, 1, 2, 1}, - {"sel(0,1<<64-1,1)==1", 0, 1<<64 - 1, 1, 1}, - {"sel(1,1<<64-1,1)==1<<64-1", 1, 1<<64 - 1, 1, 1<<64 - 1}, - {"sel(0,1,1<<64-1)==1<<64-1", 0, 1, 1<<64 - 1, 1<<64 - 1}, - {"sel(1,1,1<<64-1)==1", 1, 1, 1<<64 - 1, 1}, - {"sel(0,1<<64-1,1<<64-2)==1<<64-2", 0, 1<<64 - 1, 1<<64 - 2, 1<<64 - 2}, - {"sel(1,1<<64-1,1<<64-2)==1<<64-1", 1, 1<<64 - 1, 1<<64 - 2, 1<<64 - 1}, + {0, 1, 2, 2}, + {1, 1, 2, 1}, + {0, 1<<64 - 1, 1, 1}, + {1, 1<<64 - 1, 1, 1<<64 - 1}, + {0, 1, 1<<64 - 1, 1<<64 - 1}, + {1, 1, 1<<64 - 1, 1}, + {0, 1<<64 - 1, 1<<64 - 2, 1<<64 - 2}, + {1, 1<<64 - 1, 1<<64 - 2, 1<<64 - 1}, + {0, 0x123456789abcdef0, 0x123456789abcdef0, 0x123456789abcdef0}, + {1, 0x123456789abcdef0, 0x123456789abcdef0, 0x123456789abcdef0}, } for _, test := range tests { got := ConstantTimeSelect64(test.cond, test.a, test.b) if got != test.want { - t.Errorf("%q: unexpected result -- got %d, want %d", test.name, - got, test.want) + t.Errorf("sel(%d, %016x, %016x): got %016x, want %016x", test.cond, + test.a, test.b, got, test.want) } } } From 9eb76f8b74e365657751a176c2f72158e20dada9 Mon Sep 17 00:00:00 2001 From: Dave Collins Date: Sat, 25 Jul 2026 04:26:20 -0500 Subject: [PATCH 14/17] secp256k1: Move mul/square primitives to arith. The package currently defines modulus-agnostic 512-bit multiplication and squaring primitives as part of the 64-bit field implementation despite them also being used by the scalar implementation. They are also shared by the field implementations that will be split into separate packages in future commits to prepare for architecture-specific backend selection. In order to smooth that transition, this moves the primitives into the internal arith package and updates all call sites to use them. No functional changes are intended. The move is purely organizational and preserves the existing arithmetic. --- dcrec/secp256k1/field64.go | 203 ----------------------- dcrec/secp256k1/field64_amd64.go | 6 +- dcrec/secp256k1/field64_generic.go | 6 +- dcrec/secp256k1/internal/arith/arith.go | 209 ++++++++++++++++++++++++ dcrec/secp256k1/modnscalar.go | 4 +- 5 files changed, 219 insertions(+), 209 deletions(-) create mode 100644 dcrec/secp256k1/internal/arith/arith.go diff --git a/dcrec/secp256k1/field64.go b/dcrec/secp256k1/field64.go index 30afab1e2b..3968ca6ea6 100644 --- a/dcrec/secp256k1/field64.go +++ b/dcrec/secp256k1/field64.go @@ -653,209 +653,6 @@ func (f *FieldVal64) SquareVal(val *FieldVal64) *FieldVal64 { return f } -// field64Mul512 sets t = x * y as an unreduced 512-bit product via a row-by-row -// schoolbook multiply. -func field64Mul512(t *[8]uint64, x, y *[4]uint64) { - // The intermediate bounds and carry assumptions used by this algorithm have - // been formally verified. The verification artifacts are available in - // internal/proofs. - - a0, a1, a2, a3 := x[0], x[1], x[2], x[3] - b0, b1, b2, b3 := y[0], y[1], y[2], y[3] - - var c uint64 - - // Row 0: p0..p4 = a * b0. - // - // Note that since h3 is the upper 64 bits of the product of two uint64s: - // h3 ≤ floor((2^64-1)^2 / 2^64) = 2^64 - 2 - // - // Without any other considerations, c ≤ 1, so a loose bound is: - // p4 ≤ h3 + 1 = 2^64 - 1 < 2^64 - // - // This already shows that the carryless add in p4 is safe, however, a tight - // upper bound is more useful to prove no overflow is possible in the upper - // words of the subsequent rows. - // - // Claim: p4 ≤ 2^64 - 2 - // - // Consider the row product A*b, where A ≤ 2^256 - 1, b ≤ 2^64 - 1, then: - // A*b ≤ (2^256 - 1)(2^64 - 1) = 2^320 - 2^256 - 2^64 + 1 - // - // Next, expressing the product in base 2^256 gives: - // A*b = p4*2^256 + qlow - // - // Where qlow is the low 256 bits of the product and p4 is the integer - // quotient: - // p4 = floor(A*b / 2^256) - // qlow = A*b (mod 2^256) - // - // Finally, bound the quotient: - // p4 = floor(A*b / 2^256) - // ≤ floor((2^320 - 2^256 - 2^64 + 1) / 2^256) - // = floor(2^64 - 1 - 2^(-192) + 2^(-256)) - // ≤ 2^64 - 2 - // - // So, p4 ≤ 2^64 - 2. - h0, p0 := bits.Mul64(a0, b0) - h1, p1 := bits.Mul64(a1, b0) - h2, p2 := bits.Mul64(a2, b0) - h3, p3 := bits.Mul64(a3, b0) - p1, c = bits.Add64(p1, h0, 0) - p2, c = bits.Add64(p2, h1, c) - p3, c = bits.Add64(p3, h2, c) - p4 := h3 + c - - // Row 1: p1..p5 += a * b1. - // - // Per row 0 above, the tight bound on q4 for this row is: - // q4 ≤ 2^64 - 2 - // - // Since c ≤ 1: - // p5 ≤ q4 + 1 = 2^64 - 1 < 2^64 - // - // So, the carryless add in p5 is safe. - h0, q0 := bits.Mul64(a0, b1) - h1, q1 := bits.Mul64(a1, b1) - h2, q2 := bits.Mul64(a2, b1) - h3, q3 := bits.Mul64(a3, b1) - q1, c = bits.Add64(q1, h0, 0) - q2, c = bits.Add64(q2, h1, c) - q3, c = bits.Add64(q3, h2, c) - q4 := h3 + c - p1, c = bits.Add64(p1, q0, 0) - p2, c = bits.Add64(p2, q1, c) - p3, c = bits.Add64(p3, q2, c) - p4, c = bits.Add64(p4, q3, c) - p5 := q4 + c - - // Row 2: p2..p6 += a * b2. - // - // The same bounds calculation as row 1 applies. - h0, q0 = bits.Mul64(a0, b2) - h1, q1 = bits.Mul64(a1, b2) - h2, q2 = bits.Mul64(a2, b2) - h3, q3 = bits.Mul64(a3, b2) - q1, c = bits.Add64(q1, h0, 0) - q2, c = bits.Add64(q2, h1, c) - q3, c = bits.Add64(q3, h2, c) - q4 = h3 + c - p2, c = bits.Add64(p2, q0, 0) - p3, c = bits.Add64(p3, q1, c) - p4, c = bits.Add64(p4, q2, c) - p5, c = bits.Add64(p5, q3, c) - p6 := q4 + c - - // Row 3: p3..p7 += a * b3. - // - // The same bounds calculation as row 1 applies. - h0, q0 = bits.Mul64(a0, b3) - h1, q1 = bits.Mul64(a1, b3) - h2, q2 = bits.Mul64(a2, b3) - h3, q3 = bits.Mul64(a3, b3) - q1, c = bits.Add64(q1, h0, 0) - q2, c = bits.Add64(q2, h1, c) - q3, c = bits.Add64(q3, h2, c) - q4 = h3 + c - p3, c = bits.Add64(p3, q0, 0) - p4, c = bits.Add64(p4, q1, c) - p5, c = bits.Add64(p5, q2, c) - p6, c = bits.Add64(p6, q3, c) - p7 := q4 + c - - t[0], t[1], t[2], t[3] = p0, p1, p2, p3 - t[4], t[5], t[6], t[7] = p4, p5, p6, p7 -} - -// field64Square512 sets t = a^2 as an unreduced 512-bit product. -func field64Square512(t *[8]uint64, a *[4]uint64) { - // The intermediate bounds and carry assumptions used by this algorithm have - // been formally verified. The verification artifacts are available in - // internal/proofs. - - a0, a1, a2, a3 := a[0], a[1], a[2], a[3] - - var c uint64 - - // Off-diagonal upper-triangle products (not yet doubled). - // - // Note that since h03 is the upper 64 bits of the product of two uint64s: - // h03 ≤ floor((2^64-1)^2 / 2^64) = 2^64 - 2 - // - // Then, because c ≤ 1, a loose bound is: - // p4 ≤ h03 + 1 = 2^64 - 1 < 2^64 - // - // Therefore, it is safe to discard the carry. - p2, p1 := bits.Mul64(a0, a1) - h02, l02 := bits.Mul64(a0, a2) - h03, l03 := bits.Mul64(a0, a3) - p2, c = bits.Add64(p2, l02, 0) - p3, c := bits.Add64(h02, l03, c) - p4, _ := bits.Add64(h03, 0, c) - - h12, l12 := bits.Mul64(a1, a2) - p3, c = bits.Add64(p3, l12, 0) - p4, c = bits.Add64(p4, h12, c) - p5 := c - - // The p5 carry is safe to discard because p5 + h13 + c ≤ 2^64 - 1 (where c - // is the carry from p4 + l13). - // - // A full proof involves case analysis that is omitted here since the - // impossibility of the carry is formally proven in internal/proofs, but the - // key point is that the only way the final add could have a carry is if all - // 3 of the following conditions were simultaneously true: - // - // 1) p5_old = 1 (the carry from the earlier chain, so ≤ 1) - // 2) h13 = 2^64 - 2 (h13 ≤ 2^64 - 2 as proven previously) - // 3) c = 1 (implies p4 + l13 ≥ 2^64) - // - // However, that combination of conditions is impossible because in order - // for condition 2 to be true, a1 = a3 = 2^64 - 1, in which case l13 = 1 - // and so in order for condition 3 to also be true, p4 = 2^64 - 1. But then - // the combination of those conditions forces p5_old = 0. - h13, l13 := bits.Mul64(a1, a3) - p4, c = bits.Add64(p4, l13, 0) - p5, _ = bits.Add64(p5, h13, c) - - // Similarly, the p6 carry is safe to discard because, per above: - // h23 ≤ 2^64 - 2 - // - // Then, again c ≤ 1, so the same loose bound applies: - // p6 ≤ h23 + 1 = 2^64 - 1 < 2^64 - h23, l23 := bits.Mul64(a2, a3) - p5, c = bits.Add64(p5, l23, 0) - p6, _ := bits.Add64(h23, 0, c) - - // Double p1..p6, capturing the top carry into p7. - p1, c = bits.Add64(p1, p1, 0) - p2, c = bits.Add64(p2, p2, c) - p3, c = bits.Add64(p3, p3, c) - p4, c = bits.Add64(p4, p4, c) - p5, c = bits.Add64(p5, p5, c) - p6, c = bits.Add64(p6, p6, c) - p7 := c - - // Add the diagonal squares a[i]^2 at columns 0,2,4,6 in one carry chain. - // - // The carry on the final add is safe to discard because a < p < 2^256, so: - // (2^256 - 1)^2 = 2^512 - 2^257 + 1 < 2^512 - h0, p0 := bits.Mul64(a0, a0) - h1, l1 := bits.Mul64(a1, a1) - h2, l2 := bits.Mul64(a2, a2) - h3, l3 := bits.Mul64(a3, a3) - p1, c = bits.Add64(p1, h0, 0) - p2, c = bits.Add64(p2, l1, c) - p3, c = bits.Add64(p3, h1, c) - p4, c = bits.Add64(p4, l2, c) - p5, c = bits.Add64(p5, h2, c) - p6, c = bits.Add64(p6, l3, c) - p7, _ = bits.Add64(p7, h3, c) - - t[0], t[1], t[2], t[3] = p0, p1, p2, p3 - t[4], t[5], t[6], t[7] = p4, p5, p6, p7 -} - // field64Reduce512 reduces a 512-bit little-endian limb array modulo p in // constant time and stores the result in r. func field64Reduce512(r *[4]uint64, x *[8]uint64) { diff --git a/dcrec/secp256k1/field64_amd64.go b/dcrec/secp256k1/field64_amd64.go index bf068120f4..bef8f12cc2 100644 --- a/dcrec/secp256k1/field64_amd64.go +++ b/dcrec/secp256k1/field64_amd64.go @@ -6,6 +6,8 @@ package secp256k1 +import "github.com/decred/dcrd/dcrec/secp256k1/v4/internal/arith" + // field64UseADX reports whether the CPU supports both BMI2 (MULX) and ADX (ADCX/ADOX). var field64UseADX = func() bool { const ( @@ -59,13 +61,13 @@ func field64Square(r *[4]uint64, a *[4]uint64) { // field64MulGeneric sets r = a * b (mod p) func field64MulGeneric(r *[4]uint64, a, b *[4]uint64) { var product [8]uint64 - field64Mul512(&product, a, b) + arith.Mul512(&product, a, b) field64Reduce512(r, &product) } // field64SquareGeneric sets r = a^2 (mod p) func field64SquareGeneric(r *[4]uint64, a *[4]uint64) { var product [8]uint64 - field64Square512(&product, a) + arith.Square512(&product, a) field64Reduce512(r, &product) } diff --git a/dcrec/secp256k1/field64_generic.go b/dcrec/secp256k1/field64_generic.go index bea665bc69..6ad4be30a7 100644 --- a/dcrec/secp256k1/field64_generic.go +++ b/dcrec/secp256k1/field64_generic.go @@ -6,16 +6,18 @@ package secp256k1 +import "github.com/decred/dcrd/dcrec/secp256k1/v4/internal/arith" + // field64Mul sets r = a * b (mod p). func field64Mul(r *[4]uint64, a, b *[4]uint64) { var product [8]uint64 - field64Mul512(&product, a, b) + arith.Mul512(&product, a, b) field64Reduce512(r, &product) } // field64Square sets r = a^2 (mod p). func field64Square(r *[4]uint64, a *[4]uint64) { var product [8]uint64 - field64Square512(&product, a) + arith.Square512(&product, a) field64Reduce512(r, &product) } diff --git a/dcrec/secp256k1/internal/arith/arith.go b/dcrec/secp256k1/internal/arith/arith.go new file mode 100644 index 0000000000..13b35d6737 --- /dev/null +++ b/dcrec/secp256k1/internal/arith/arith.go @@ -0,0 +1,209 @@ +// Copyright (c) 2026 The Decred developers +// Use of this source code is governed by an ISC +// license that can be found in the LICENSE file. + +package arith + +import "math/bits" + +// Mul512 sets t = x * y as an unreduced 512-bit product. +func Mul512(t *[8]uint64, x, y *[4]uint64) { + // The intermediate bounds and carry assumptions used by this algorithm have + // been formally verified. The verification artifacts are available in + // internal/proofs. + + a0, a1, a2, a3 := x[0], x[1], x[2], x[3] + b0, b1, b2, b3 := y[0], y[1], y[2], y[3] + + var c uint64 + + // Row 0: p0..p4 = a * b0. + // + // Note that since h3 is the upper 64 bits of the product of two uint64s: + // h3 ≤ floor((2^64-1)^2 / 2^64) = 2^64 - 2 + // + // Without any other considerations, c ≤ 1, so a loose bound is: + // p4 ≤ h3 + 1 = 2^64 - 1 < 2^64 + // + // This already shows that the carryless add in p4 is safe, however, a tight + // upper bound is more useful to prove no overflow is possible in the upper + // words of the subsequent rows. + // + // Claim: p4 ≤ 2^64 - 2 + // + // Consider the row product A*b, where A ≤ 2^256 - 1, b ≤ 2^64 - 1, then: + // A*b ≤ (2^256 - 1)(2^64 - 1) = 2^320 - 2^256 - 2^64 + 1 + // + // Next, expressing the product in base 2^256 gives: + // A*b = p4*2^256 + qlow + // + // Where qlow is the low 256 bits of the product and p4 is the integer + // quotient: + // p4 = floor(A*b / 2^256) + // qlow = A*b (mod 2^256) + // + // Finally, bound the quotient: + // p4 = floor(A*b / 2^256) + // ≤ floor((2^320 - 2^256 - 2^64 + 1) / 2^256) + // = floor(2^64 - 1 - 2^(-192) + 2^(-256)) + // ≤ 2^64 - 2 + // + // So, p4 ≤ 2^64 - 2. + h0, p0 := bits.Mul64(a0, b0) + h1, p1 := bits.Mul64(a1, b0) + h2, p2 := bits.Mul64(a2, b0) + h3, p3 := bits.Mul64(a3, b0) + p1, c = bits.Add64(p1, h0, 0) + p2, c = bits.Add64(p2, h1, c) + p3, c = bits.Add64(p3, h2, c) + p4 := h3 + c + + // Row 1: p1..p5 += a * b1. + // + // Per row 0 above, the tight bound on q4 for this row is: + // q4 ≤ 2^64 - 2 + // + // Since c ≤ 1: + // p5 ≤ q4 + 1 = 2^64 - 1 < 2^64 + // + // So, the carryless add in p5 is safe. + h0, q0 := bits.Mul64(a0, b1) + h1, q1 := bits.Mul64(a1, b1) + h2, q2 := bits.Mul64(a2, b1) + h3, q3 := bits.Mul64(a3, b1) + q1, c = bits.Add64(q1, h0, 0) + q2, c = bits.Add64(q2, h1, c) + q3, c = bits.Add64(q3, h2, c) + q4 := h3 + c + p1, c = bits.Add64(p1, q0, 0) + p2, c = bits.Add64(p2, q1, c) + p3, c = bits.Add64(p3, q2, c) + p4, c = bits.Add64(p4, q3, c) + p5 := q4 + c + + // Row 2: p2..p6 += a * b2. + // + // The same bounds calculation as row 1 applies. + h0, q0 = bits.Mul64(a0, b2) + h1, q1 = bits.Mul64(a1, b2) + h2, q2 = bits.Mul64(a2, b2) + h3, q3 = bits.Mul64(a3, b2) + q1, c = bits.Add64(q1, h0, 0) + q2, c = bits.Add64(q2, h1, c) + q3, c = bits.Add64(q3, h2, c) + q4 = h3 + c + p2, c = bits.Add64(p2, q0, 0) + p3, c = bits.Add64(p3, q1, c) + p4, c = bits.Add64(p4, q2, c) + p5, c = bits.Add64(p5, q3, c) + p6 := q4 + c + + // Row 3: p3..p7 += a * b3. + // + // The same bounds calculation as row 1 applies. + h0, q0 = bits.Mul64(a0, b3) + h1, q1 = bits.Mul64(a1, b3) + h2, q2 = bits.Mul64(a2, b3) + h3, q3 = bits.Mul64(a3, b3) + q1, c = bits.Add64(q1, h0, 0) + q2, c = bits.Add64(q2, h1, c) + q3, c = bits.Add64(q3, h2, c) + q4 = h3 + c + p3, c = bits.Add64(p3, q0, 0) + p4, c = bits.Add64(p4, q1, c) + p5, c = bits.Add64(p5, q2, c) + p6, c = bits.Add64(p6, q3, c) + p7 := q4 + c + + t[0], t[1], t[2], t[3] = p0, p1, p2, p3 + t[4], t[5], t[6], t[7] = p4, p5, p6, p7 +} + +// Square512 sets t = a^2 as an unreduced 512-bit product. +func Square512(t *[8]uint64, a *[4]uint64) { + // The intermediate bounds and carry assumptions used by this algorithm have + // been formally verified. The verification artifacts are available in + // internal/proofs. + + a0, a1, a2, a3 := a[0], a[1], a[2], a[3] + + var c uint64 + + // Off-diagonal upper-triangle products (not yet doubled). + // + // Note that since h03 is the upper 64 bits of the product of two uint64s: + // h03 ≤ floor((2^64-1)^2 / 2^64) = 2^64 - 2 + // + // Then, because c ≤ 1, a loose bound is: + // p4 ≤ h03 + 1 = 2^64 - 1 < 2^64 + // + // Therefore, it is safe to discard the carry. + p2, p1 := bits.Mul64(a0, a1) + h02, l02 := bits.Mul64(a0, a2) + h03, l03 := bits.Mul64(a0, a3) + p2, c = bits.Add64(p2, l02, 0) + p3, c := bits.Add64(h02, l03, c) + p4, _ := bits.Add64(h03, 0, c) + + h12, l12 := bits.Mul64(a1, a2) + p3, c = bits.Add64(p3, l12, 0) + p4, c = bits.Add64(p4, h12, c) + p5 := c + + // The p5 carry is safe to discard because p5 + h13 + c ≤ 2^64 - 1 (where c + // is the carry from p4 + l13). + // + // A full proof involves case analysis that is omitted here since the + // impossibility of the carry is formally proven in internal/proofs, but the + // key point is that the only way the final add could have a carry is if all + // 3 of the following conditions were simultaneously true: + // + // 1) p5_old = 1 (the carry from the earlier chain, so ≤ 1) + // 2) h13 = 2^64 - 2 (h13 ≤ 2^64 - 2 as proven previously) + // 3) c = 1 (implies p4 + l13 ≥ 2^64) + // + // However, that combination of conditions is impossible because in order + // for condition 2 to be true, a1 = a3 = 2^64 - 1, in which case l13 = 1 + // and so in order for condition 3 to also be true, p4 = 2^64 - 1. But then + // the combination of those conditions forces p5_old = 0. + h13, l13 := bits.Mul64(a1, a3) + p4, c = bits.Add64(p4, l13, 0) + p5, _ = bits.Add64(p5, h13, c) + + // Similarly, the p6 carry is safe to discard because, per above: + // h23 ≤ 2^64 - 2 + // + // Then, again c ≤ 1, so the same loose bound applies: + // p6 ≤ h23 + 1 = 2^64 - 1 < 2^64 + h23, l23 := bits.Mul64(a2, a3) + p5, c = bits.Add64(p5, l23, 0) + p6, _ := bits.Add64(h23, 0, c) + + // Double p1..p6, capturing the top carry into p7. + p1, c = bits.Add64(p1, p1, 0) + p2, c = bits.Add64(p2, p2, c) + p3, c = bits.Add64(p3, p3, c) + p4, c = bits.Add64(p4, p4, c) + p5, c = bits.Add64(p5, p5, c) + p6, c = bits.Add64(p6, p6, c) + p7 := c + + // Add the diagonal squares a[i]^2 at columns 0,2,4,6 in one carry chain. + // + // The carry on the final add is safe to discard because a < p < 2^256, so: + // (2^256 - 1)^2 = 2^512 - 2^257 + 1 < 2^512 + h0, p0 := bits.Mul64(a0, a0) + h1, l1 := bits.Mul64(a1, a1) + h2, l2 := bits.Mul64(a2, a2) + h3, l3 := bits.Mul64(a3, a3) + p1, c = bits.Add64(p1, h0, 0) + p2, c = bits.Add64(p2, l1, c) + p3, c = bits.Add64(p3, h1, c) + p4, c = bits.Add64(p4, l2, c) + p5, c = bits.Add64(p5, h2, c) + p6, c = bits.Add64(p6, l3, c) + p7, _ = bits.Add64(p7, h3, c) + + t[0], t[1], t[2], t[3] = p0, p1, p2, p3 + t[4], t[5], t[6], t[7] = p4, p5, p6, p7 +} diff --git a/dcrec/secp256k1/modnscalar.go b/dcrec/secp256k1/modnscalar.go index 26d13c5dee..34a6d49a66 100644 --- a/dcrec/secp256k1/modnscalar.go +++ b/dcrec/secp256k1/modnscalar.go @@ -752,7 +752,7 @@ func scalar64Reduce512(r *[4]uint64, x *[8]uint64) { // s3.Mul2(s, s2).AddInt(1) so that s3 = (s * s2) + 1. func (s *ModNScalar) Mul2(a, b *ModNScalar) *ModNScalar { var product [8]uint64 - field64Mul512(&product, &a.n, &b.n) + arith.Mul512(&product, &a.n, &b.n) scalar64Reduce512(&s.n, &product) return s } @@ -773,7 +773,7 @@ func (s *ModNScalar) Mul(val *ModNScalar) *ModNScalar { // s3.SquareVal(s).Mul(s) so that s3 = s^2 * s = s^3. func (s *ModNScalar) SquareVal(val *ModNScalar) *ModNScalar { var product [8]uint64 - field64Square512(&product, &val.n) + arith.Square512(&product, &val.n) scalar64Reduce512(&s.n, &product) return s } From a8bfc6b717dda571162f02ceaaf44a6c3bcf3358 Mon Sep 17 00:00:00 2001 From: Dave Collins Date: Sun, 26 Jul 2026 04:21:21 -0500 Subject: [PATCH 15/17] secp256k1: Add Mul512 tests. This adds dedicated tests for the 512-bit multiplication in the new internal arith package. While it is already already well exercised by the upper layers that use it, testing primitives directly helps ensure correctness in isolation. --- dcrec/secp256k1/internal/arith/arith_test.go | 244 +++++++++++++++++++ 1 file changed, 244 insertions(+) create mode 100644 dcrec/secp256k1/internal/arith/arith_test.go diff --git a/dcrec/secp256k1/internal/arith/arith_test.go b/dcrec/secp256k1/internal/arith/arith_test.go new file mode 100644 index 0000000000..c2bc42abce --- /dev/null +++ b/dcrec/secp256k1/internal/arith/arith_test.go @@ -0,0 +1,244 @@ +// Copyright (c) 2026 The Decred developers +// Use of this source code is governed by an ISC +// license that can be found in the LICENSE file. + +package arith + +import ( + "encoding/binary" + "fmt" + "math/big" + mrand "math/rand" + "testing" + "time" +) + +// mustBig converts the passed hex string into a big integer and will panic if +// there is an error. This is only provided for the hard-coded constants so +// errors in the source code can be detected. It will only (and must only) be +// called with hard-coded values. +func mustBig(s string) *big.Int { + val, ok := new(big.Int).SetString(s, 16) + if !ok { + panic("failed to parse big integer from hex: " + s) + } + return val +} + +// mustBigToUint256 converts a [big.Int] to an array of 4 uint64s representing a +// 256-bit little-endian value. It will panic if the big integer is larger than +// the maximum 256-bit value. This is only provided for the hard-coded +// constants and randomly-generated values that are known to be in range, so +// errors in the source code can be detected. +func mustBigToUint256(v *big.Int) [4]uint64 { + if v.BitLen() > 256 { + panic(fmt.Sprintf("big integer %x is larger than max uint256", v)) + } + + var buf [32]byte + v.FillBytes(buf[:]) + + var result [4]uint64 + for i := 0; i < 4; i++ { + result[i] = binary.BigEndian.Uint64(buf[32-((i+1)*8):]) + } + return result +} + +// uint512ToBig converts an array of 8 uint64s representing a 512-bit +// little-endian value to a [big.Int]. +func uint512ToBig(v [8]uint64) *big.Int { + var buf [64]byte + for i := 0; i < 8; i++ { + binary.BigEndian.PutUint64(buf[i*8:], v[7-i]) + } + + return new(big.Int).SetBytes(buf[:]) +} + +// randBig returns a random 256-bit [big.Int] created from the passed rng. +func randBig(t *testing.T, rng *mrand.Rand) *big.Int { + t.Helper() + + var buf [32]byte + if _, err := rng.Read(buf[:]); err != nil { + t.Fatalf("failed to read random: %v", err) + } + + return new(big.Int).SetBytes(buf[:]) +} + +// TestMul512 ensures [Mul512] returns the expected result by comparing the +// product against the [big.Int] result. It also tests commutativity and full +// buffer replacement. +func TestMul512(t *testing.T) { + tests := []struct { + name string // test description + x, y string // hex encoded multiplicands + }{{ + name: "all zero", + x: "0", + y: "0", + }, { + name: "zero * identity", + x: "0", + y: "1", + }, { + name: "zero * max uint256", + x: "0", + y: "ffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffff", + }, { + name: "identity", + x: "1", + y: "1", + }, { + name: "identity * max uint256", + x: "1", + y: "ffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffff", + }, { + name: "small values", + x: "2", + y: "3", + }, { + name: "unbalanced small values", + x: "1", + y: "ffffffffffffffff", + }, { + name: "small * max uint256", + x: "2", + y: "ffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffff", + }, { + name: "max uint64 * max uint64", + x: "ffffffffffffffff", + y: "ffffffffffffffff", + }, { + name: "max uint64 * 2*max uint64", + x: "ffffffffffffffff", + y: "1fffffffffffffffe", + }, { + name: "2^1 * 2^1", + x: "2", + y: "2", + }, { + name: "2^8 * 2^8", + x: "100", + y: "100", + }, { + name: "2^16 * 2^16", + x: "10000", + y: "10000", + }, { + name: "2^32 * 2^32", + x: "100000000", + y: "100000000", + }, { + name: "2^64 * 2^64", + x: "10000000000000000", + y: "10000000000000000", + }, { + name: "2^128 * 2^128", + x: "100000000000000000000000000000000", + y: "100000000000000000000000000000000", + }, { + name: "2^192 * 2^64", + x: "1000000000000000000000000000000000000000000000000", + y: "10000000000000000", + }, { + name: "2^255 * 2^255", + x: "8000000000000000000000000000000000000000000000000000000000000000", + y: "8000000000000000000000000000000000000000000000000000000000000000", + }, { + name: "max uint256 minus one", + x: "fffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffe", + y: "fffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffe", + }, { + name: "max uint256 * max uint256", + x: "ffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffff", + y: "ffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffff", + }, { + name: "alternating bits", + x: "a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5", + y: "5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a", + }, { + name: "alternating bits 2", + x: "aaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaa", + y: "5555555555555555555555555555555555555555555555555555555555555555", + }, { + name: "carry propagation through zero limbs", + x: "ffffffffffffffffffffffffffffffff00000000000000000000000000000000", + y: "100000000000000000000000000000001", + }} + + for _, test := range tests { + xBig, yBig := mustBig(test.x), mustBig(test.y) + x, y := mustBigToUint256(xBig), mustBigToUint256(yBig) + want := new(big.Int).Mul(xBig, yBig) + + // Fill the array that will be used to store the product with non-zero + // values to ensure the entire array is overwritten. + var product [8]uint64 + for i := 0; i < 8; i++ { + product[i] = 0xffffffffffffffff + } + + // Compute the product with [Mul512] and ensure it matches the same + // result produced by [big.Int]. + Mul512(&product, &x, &y) + if got := uint512ToBig(product); got.Cmp(want) != 0 { + t.Errorf("%s: incorrect product\n x: %064x\n y: %064x\n"+ + "got: %0128x\nwant: %0128x", test.name, xBig, yBig, got, want) + } + + // Ensure commutativity works properly with [Mul512] by computing the + // product again with the operands swapped. Poison the array used to + // store the product again to ensure the opposite order overwrites the + // entire array too. + for i := 0; i < 8; i++ { + product[i] = 0xffffffffffffffff + } + Mul512(&product, &y, &x) + if got := uint512ToBig(product); got.Cmp(want) != 0 { + t.Errorf("%s (swapped): incorrect product\n y: %064x\n x: %064x\n"+ + "got: %0128x\nwant: %0128x", test.name, yBig, xBig, got, want) + } + } +} + +// TestMul512Random ensures that multiplying randomly-generated values via +// [Mul512] returns the expected result by comparing the product against the +// [big.Int] result. It also tests commutativity. +func TestMul512Random(t *testing.T) { + // Use a unique random seed each test instance and log it if the tests fail. + seed := time.Now().UnixNano() + rng := mrand.New(mrand.NewSource(seed)) + defer func(t *testing.T, seed int64) { + if t.Failed() { + t.Logf("random seed: %d", seed) + } + }(t, seed) + + for i := 0; i < 1000; i++ { + // Generate random [big.Int] operands and the expected product. + xBig, yBig := randBig(t, rng), randBig(t, rng) + x, y := mustBigToUint256(xBig), mustBigToUint256(yBig) + want := new(big.Int).Mul(xBig, yBig) + + // Compute the product with [Mul512] and ensure it matches the same + // result produced by [big.Int]. + var product [8]uint64 + Mul512(&product, &x, &y) + if got := uint512ToBig(product); got.Cmp(want) != 0 { + t.Errorf("incorrect product\n x: %064x\n y: %064x\ngot: %0128x\n"+ + "want: %0128x", xBig, yBig, got, want) + } + + // Ensure commutativity works properly with [Mul512] by computing the + // product again with the operands swapped. + var product2 [8]uint64 + Mul512(&product2, &y, &x) + if got := uint512ToBig(product2); got.Cmp(want) != 0 { + t.Errorf("incorrect product\n y: %064x\n x: %064x\ngot: %0128x\n"+ + "want: %0128x", yBig, xBig, got, want) + } + } +} From 4f158fe57974dd4b36252aa69f874c4bf8175e68 Mon Sep 17 00:00:00 2001 From: Dave Collins Date: Sun, 26 Jul 2026 05:07:27 -0500 Subject: [PATCH 16/17] secp256k1: Add Square512 tests. This adds dedicated tests for the 512-bit squaring in the new internal arith package. While it is already already well exercised by the upper layers that use it, testing primitives directly helps ensure correctness in isolation. --- dcrec/secp256k1/internal/arith/arith_test.go | 151 +++++++++++++++++++ 1 file changed, 151 insertions(+) diff --git a/dcrec/secp256k1/internal/arith/arith_test.go b/dcrec/secp256k1/internal/arith/arith_test.go index c2bc42abce..31f9469409 100644 --- a/dcrec/secp256k1/internal/arith/arith_test.go +++ b/dcrec/secp256k1/internal/arith/arith_test.go @@ -242,3 +242,154 @@ func TestMul512Random(t *testing.T) { } } } + +// TestSquare512 ensures [Square512] returns the expected result by comparing it +// against the [big.Int] result. It also tests full buffer replacement and that +// [Mul512] and [Square512] agree on the result. +func TestSquare512(t *testing.T) { + tests := []struct { + name string // test description + a string // hex encoded value to square + }{{ + name: "zero", + a: "0", + }, { + name: "identity", + a: "1", + }, { + name: "small value", + a: "2", + }, { + name: "2^4", + a: "10", + }, { + name: "2^8", + a: "100", + }, { + name: "2^16", + a: "10000", + }, { + name: "2^32", + a: "100000000", + }, { + name: "max uint64", + a: "ffffffffffffffff", + }, { + name: "2^64", + a: "10000000000000000", + }, { + name: "max uint128 minus one", + a: "fffffffffffffffffffffffffffffffe", + }, { + name: "max uint128", + a: "ffffffffffffffffffffffffffffffff", + }, { + name: "2^128", + a: "100000000000000000000000000000000", + }, { + name: "max uint192", + a: "ffffffffffffffffffffffffffffffffffffffffffffffff", + }, { + name: "2^192", + a: "1000000000000000000000000000000000000000000000000", + }, { + name: "2^255", + a: "8000000000000000000000000000000000000000000000000000000000000000", + }, { + name: "near max uint256", + a: "ffffffffffffffffffffffffffffffffffffffffffffffffffffffffffff0001", + }, { + name: "max uint256 minus one", + a: "fffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffe", + }, { + name: "max uint256", + a: "ffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffff", + }, { + name: "sparse bits", + a: "8000000000000000000000000000000100000000000000000000000000000001", + }, { + name: "high and low bits", + a: "8000000000000000000000000000000000000000000000000000000000000001", + }, { + name: "alternating bits", + a: "a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5a5", + }, { + name: "alternating bits 2", + a: "aaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaa", + }, { + name: "alternating bits 3", + a: "5555555555555555555555555555555555555555555555555555555555555555", + }, { + name: "carry propagation through zero limbs", + a: "ffffffffffffffffffffffffffffffff00000000000000000000000000000001", + }} + + for _, test := range tests { + aBig := mustBig(test.a) + a := mustBigToUint256(aBig) + want := new(big.Int).Mul(aBig, aBig) + + // Fill the array that will be used to store the result with non-zero + // values to ensure the entire array is overwritten. + var result [8]uint64 + for i := 0; i < 8; i++ { + result[i] = 0xffffffffffffffff + } + + // Compute the result with [Square512] and ensure it matches the same + // result produced by [big.Int]. + Square512(&result, &a) + if got := uint512ToBig(result); got.Cmp(want) != 0 { + t.Errorf("%s: incorrect square\n a: %064x\n got: %0128x\n"+ + "want: %0128x", test.name, aBig, got, want) + } + + // Ensure [Mul512] and [Square512] agree. + var mulResult [8]uint64 + Mul512(&mulResult, &a, &a) + if result != mulResult { + t.Errorf("%s: square and mul disagree\n a: %064x\n"+ + " square: %0128x\n mul: %0128x", test.name, aBig, + uint512ToBig(result), uint512ToBig(mulResult)) + } + } +} + +// TestSquare512Random ensures that squaring randomly-generated values via +// [Square512] returns the expected result by comparing it against the [big.Int] +// result. It also tests that [Mul512] and [Square512] agree on the result. +func TestSquare512Random(t *testing.T) { + // Use a unique random seed each test instance and log it if the tests fail. + seed := time.Now().UnixNano() + rng := mrand.New(mrand.NewSource(seed)) + defer func(t *testing.T, seed int64) { + if t.Failed() { + t.Logf("random seed: %d", seed) + } + }(t, seed) + + for i := 0; i < 1000; i++ { + // Generate a random [big.Int] operand and the expected result. + aBig := randBig(t, rng) + a := mustBigToUint256(aBig) + want := new(big.Int).Mul(aBig, aBig) + + // Compute the result with [Square512] and ensure it matches the same + // result produced by [big.Int]. + var result [8]uint64 + Square512(&result, &a) + if got := uint512ToBig(result); got.Cmp(want) != 0 { + t.Errorf("incorrect square\n a: %064x\n got: %0128x\n"+ + "want: %0128x", aBig, got, want) + } + + // Ensure [Mul512] and [Square512] agree. + var mulResult [8]uint64 + Mul512(&mulResult, &a, &a) + if result != mulResult { + t.Errorf("square and mul disagree\n a: %064x\n"+ + " square: %0128x\n mul: %0128x", aBig, + uint512ToBig(result), uint512ToBig(mulResult)) + } + } +} From 309029b30e11eb34f405f0376e7dd1d866544cc3 Mon Sep 17 00:00:00 2001 From: Dave Collins Date: Sun, 26 Jul 2026 01:30:44 -0500 Subject: [PATCH 17/17] secp256k1: Add arith README.md. --- dcrec/secp256k1/internal/arith/README.md | 23 +++++++++++++++++++++++ dcrec/secp256k1/internal/arith/arith.go | 2 ++ 2 files changed, 25 insertions(+) create mode 100644 dcrec/secp256k1/internal/arith/README.md diff --git a/dcrec/secp256k1/internal/arith/README.md b/dcrec/secp256k1/internal/arith/README.md new file mode 100644 index 0000000000..09660b20b7 --- /dev/null +++ b/dcrec/secp256k1/internal/arith/README.md @@ -0,0 +1,23 @@ +arith +===== + +[![Build Status](https://github.com/decred/dcrd/workflows/Build%20and%20Test/badge.svg)](https://github.com/decred/dcrd/actions) +[![ISC License](https://img.shields.io/badge/license-ISC-blue.svg)](http://copyfree.org) +[![Doc](https://img.shields.io/badge/doc-reference-blue.svg)](https://pkg.go.dev/github.com/decred/dcrd/dcrec/secp256k1/v4/internal/arith) + +## Overview + +This provides low-level constant-time primitives and modulus-agnostic arithmetic +shared by multiple implementations in this module. + +See [internal/proofs/README.md](../proofs/README.md) for formal verification of +the 512-bit multiplication and 512-bit squaring. + +## Installation and Updating + +This package is internal and therefore is neither directly installed nor needs +to be manually updated. + +## License + +Package arith is licensed under the [copyfree](http://copyfree.org) ISC License. diff --git a/dcrec/secp256k1/internal/arith/arith.go b/dcrec/secp256k1/internal/arith/arith.go index 13b35d6737..b6e75cd2f2 100644 --- a/dcrec/secp256k1/internal/arith/arith.go +++ b/dcrec/secp256k1/internal/arith/arith.go @@ -2,6 +2,8 @@ // Use of this source code is governed by an ISC // license that can be found in the LICENSE file. +// Package arith provides low-level constant-time primitives and +// modulus-agnostic arithmetic. package arith import "math/bits"