diff --git a/lib/crypto/src/lib.rs b/lib/crypto/src/lib.rs index d7ec89880e..bf75d0b084 100644 --- a/lib/crypto/src/lib.rs +++ b/lib/crypto/src/lib.rs @@ -36,6 +36,7 @@ pub mod kdf; mod keccak; pub mod mlkem; mod montgomery; +mod naf; mod nist; mod poly1305; pub mod rsa; diff --git a/lib/crypto/src/montgomery.rs b/lib/crypto/src/montgomery.rs index 71934f8a4d..d8450a434e 100644 --- a/lib/crypto/src/montgomery.rs +++ b/lib/crypto/src/montgomery.rs @@ -14,6 +14,8 @@ use core::marker::PhantomData; /// constants for R = 2^(64N), evaluated at compile time. pub(crate) trait Modulus: 'static { const P: [u64; N]; + /// How a reduction step forms q p. + const FORM: Form = Form::Generic; /// -p^-1 mod 2^64. const N0: u64 = n0(Self::P[0]); /// R mod p, the Montgomery form of one. @@ -22,6 +24,19 @@ pub(crate) trait Modulus: 'static { const R2: [u64; N] = shift_mod(&one(), &Self::P, 128 * N); } +/// The shape of a modulus that its reduction step exploits. +#[derive(Clone, Copy, PartialEq, Eq)] +pub(crate) enum Form { + Generic, + /// p = 2^256 - 2^224 + 2^192 + 2^96 - 1. Its low word is 2^64 - 1, so + /// -p^-1 = 1 mod 2^64 and a step's quotient digit is the low word q + /// itself; t0 + q (2^64 - 1) + q (2^32 - 1) 2^64 leaves q 2^32 at word + /// one, and only q p3 needs a multiplication (Gueron and Krasnov, "Fast + /// prime field elliptic-curve cryptography with 256-bit primes", J. + /// Cryptogr. Eng. 2015). + P256, +} + /// x * R mod p, the Montgomery form of an integer x < p, for constants. pub(crate) const fn encode(x: &[u64; N], p: &[u64; N]) -> [u64; N] { shift_mod(x, p, 64 * N) @@ -99,24 +114,7 @@ impl, const N: usize> Residue { /// reduction step (CIOS; Koc, Acar and Kaliski, "Analyzing and comparing /// Montgomery multiplication algorithms", IEEE Micro 1996). pub(crate) fn mul(&self, b: &Self) -> Self { - let mut t = [0u64; N]; - let mut high = 0u64; - for i in 0..N { - let mut c = 0; - for j in 0..N { - (t[j], c) = mac(t[j], self.0[i], b.0[j], c); - } - let (top, overflow) = high.overflowing_add(c); - let q = t[0].wrapping_mul(M::N0); - let (_, mut c) = mac(t[0], q, M::P[0], 0); - for j in 1..N { - (t[j - 1], c) = mac(t[j], q, M::P[j], c); - } - let (top, carry) = top.overflowing_add(c); - t[N - 1] = top; - // t < 2p after each step, so at most one of the two is set. - high = u64::from(overflow) | u64::from(carry); - } + let (t, high) = product::(&self.0, &b.0); Self(reduce(&t, high, &M::P), PhantomData) } pub(crate) fn square(&self) -> Self { @@ -193,6 +191,48 @@ pub(crate) fn to_be(words: &[u64; N], bytes: &mut [u8]) { } } +/// a b / R as N words and a top bit, below 2p: each product row followed by +/// its reduction step (CIOS). +#[inline(always)] +const fn product, const N: usize>(a: &[u64; N], b: &[u64; N]) -> ([u64; N], u64) { + let mut t = [0u64; N]; + let mut high = 0u64; + let mut i = 0; + while i < N { + let mut c = 0; + let mut j = 0; + while j < N { + (t[j], c) = mac(t[j], a[i], b[j], c); + j += 1; + } + let (top, overflow) = high.overflowing_add(c); + let q = t[0].wrapping_mul(M::N0); + let mut c; + match M::FORM { + Form::Generic => { + (_, c) = mac(t[0], q, M::P[0], 0); + let mut j = 1; + while j < N { + (t[j - 1], c) = mac(t[j], q, M::P[j], c); + j += 1; + } + } + Form::P256 => { + let carry; + (t[0], carry) = t[1].overflowing_add(q << 32); + (t[1], c) = mac(t[2], q >> 32, 1, carry as u64); + (t[2], c) = mac(t[3], q, M::P[3], c); + } + } + let (top, carry) = top.overflowing_add(c); + t[N - 1] = top; + // t < 2p after each step, so at most one of the two is set. + high = (overflow | carry) as u64; + i += 1; + } + (t, high) +} + /// v mod p for v + high * 2^(64N) < 2p: p subtracted under a mask. fn reduce(v: &[u64; N], high: u64, p: &[u64; N]) -> [u64; N] { let mut d = [0; N]; @@ -210,11 +250,11 @@ fn reduce(v: &[u64; N], high: u64, p: &[u64; N]) -> [u64; N] { /// (lo, hi) of t + a * b + c, which cannot exceed 2^128 - 1. #[inline(always)] -pub(crate) fn mac(t: u64, a: u64, b: u64, c: u64) -> (u64, u64) { - let v = u128::from(a) - .wrapping_mul(u128::from(b)) - .wrapping_add(u128::from(t)) - .wrapping_add(u128::from(c)); +pub(crate) const fn mac(t: u64, a: u64, b: u64, c: u64) -> (u64, u64) { + let v = (a as u128) + .wrapping_mul(b as u128) + .wrapping_add(t as u128) + .wrapping_add(c as u128); (v as u64, (v >> 64) as u64) } diff --git a/lib/crypto/src/naf.rs b/lib/crypto/src/naf.rs new file mode 100644 index 0000000000..e2c83c7c7f --- /dev/null +++ b/lib/crypto/src/naf.rs @@ -0,0 +1,39 @@ +//! Width-w non-adjacent forms of public scalars, for variable-time +//! verification (Hankerson, Menezes and Vanstone, Guide to Elliptic Curve +//! Cryptography, algorithm 3.35). + +/// The width-w NAF of k into `digits`, w in 2..=8: odd digits below 2^(w-1) +/// in magnitude, any two nonzero ones at least w positions apart. k must be +/// below 2^(64N) - 2^7, so that adding a digit never carries out, and +/// `digits` must be at least one longer than k's bit length. Variable time: +/// k is public. +pub(crate) fn naf(k: &[u64; N], w: u32, digits: &mut [i8]) { + debug_assert!((2..=8).contains(&w)); + let mut k = *k; + digits.fill(0); + let mut i = 0; + while k != [0; N] { + if k[0] & 1 == 1 { + let window = (k[0] & ((1 << w) - 1)) as i32; + let d = if window >= 1 << (w - 1) { window - (1 << w) } else { window }; + digits[i] = d as i8; + if d >= 0 { + k[0] -= d as u64; + } else { + let (s, mut carry) = k[0].overflowing_add((-d) as u64); + k[0] = s; + for word in &mut k[1..] { + if !carry { + break; + } + (*word, carry) = word.overflowing_add(1); + } + } + } + for j in 0..N - 1 { + k[j] = k[j] >> 1 | k[j + 1] << 63; + } + k[N - 1] >>= 1; + i += 1; + } +} diff --git a/lib/crypto/src/nist.rs b/lib/crypto/src/nist.rs index f068049d49..3924148948 100644 --- a/lib/crypto/src/nist.rs +++ b/lib/crypto/src/nist.rs @@ -1,7 +1,10 @@ -//! SEC 1 NIST P-256/P-384 arithmetic. ECDH uses fixed-width scalars and a -//! double-and-add schedule with complete masked handling of exceptional points. -//! SEC 1 decoding, DER decoding, ECDSA verification and key-validity refusals -//! depend on public data. Only uncompressed SEC 1 points are accepted. +//! SEC 1 NIST P-256/P-384 arithmetic. Secret scalar multiplication (ECDH +//! and key generation) runs Booth's signed 5-bit windows on a fixed schedule +//! with masked table reads and complete masked handling of exceptional +//! additions. SEC 1 decoding, DER decoding, ECDSA verification and +//! key-validity refusals depend on public data; verification multiplies in +//! variable time, through `PublicScalar` only. Only uncompressed SEC 1 +//! points are accepted. #![allow( clippy::needless_range_loop, reason = "masked selects index parallel fixed-size arrays, as `montgomery` does" @@ -9,7 +12,7 @@ use crate::{ Error, - montgomery::{self, Modulus, Residue, below, from_be, to_be}, + montgomery::{self, Form, Modulus, Residue, below, from_be, to_be}, secret::{Workspace, wipe64}, }; @@ -73,11 +76,45 @@ impl, const N: usize> Scalar { } Ok(scalar) } - fn bit(&self, i: usize) -> u64 { - (self.0[i / 64] >> (i % 64)) & 1 + /// Booth digit i in [-16, 16]: with u the six bits 5i - 1 ..= 5i + 4 + /// (bit -1 and bits past the top are zero), d = (u + 1) / 2 - 32 u_5, so + /// that the scalar is the sum of d_i 2^(5i). The index is public; the + /// arithmetic does not branch. + fn booth(&self, i: usize) -> i8 { + let bits = |from: usize| -> u64 { + // The six bits from `from`, a public position. + let (word, shift) = (from / 64, from % 64); + let low = if word < N { self.0[word] >> shift } else { 0 }; + let high = if shift > 58 && word + 1 < N { self.0[word + 1] << (64 - shift) } else { 0 }; + (low | high) & 63 + }; + let u = if i == 0 { bits(0) << 1 & 63 } else { bits(5 * i - 1) }; + (((u + 1) >> 1) as i8).wrapping_sub(((u >> 5) << 5) as i8) + } +} + +/// A public scalar below n, for variable-time verification only; a secret +/// [`Scalar`] does not convert to it. +struct PublicScalar([u64; N]); + +impl PublicScalar { + /// The width-w NAF; k < n < 2^(64N) - 2^7, and 385 digits hold N = 6. + fn naf(&self, w: u32) -> [i8; 385] { + let mut digits = [0i8; 385]; + crate::naf::naf(&self.0, w, &mut digits); + digits } } +/// The windows of [`Point::mul`]: 5-bit Booth windows over 64N bits, plus +/// the one that absorbs the top window's carry. +const fn windows(n: usize) -> usize { + (64 * n) / 5 + 1 +} + +/// The odd multiples P, 3P, ..., 15P of a public point, for verification. +type OddMultiples = [Point; 8]; + /// Jacobian coordinates (X : Y : Z) for x = X/Z^2, y = Y/Z^3; Z = 0 is the /// point at infinity. pub(crate) struct Point, const N: usize> { @@ -151,13 +188,100 @@ impl, const N: usize> Point { let out = Self::select(&out, self, q.z.is_zero()); Self::select(&out, q, self.z.is_zero()) } - /// scalar * self on a fixed schedule. + /// scalar * self on a fixed schedule: Booth's signed 5-bit windows over + /// the multiples P, 2P, ..., 16P, each window five doublings and one + /// complete addition of a masked table read. self is public; neither + /// branches nor addresses depend on the scalar. pub(crate) fn mul(&self, scalar: &Scalar) -> Self { + let mut table = [*self; 16]; + table[1] = self.double(); + for j in 2..16 { + table[j] = table[j - 1].add(self); + } + let top = windows(N) - 1; + let mut r = Self::lookup(&table, scalar.booth(top)); + for i in (0..top).rev() { + for _ in 0..5 { + r = r.double(); + } + r = r.add(&Self::lookup(&table, scalar.booth(i))); + } + r + } + /// d P from table[j] = (j + 1) P, for d in [-16, 16], reading every + /// entry: the identity for zero, and Y negated by mask for d < 0. + fn lookup(table: &[Self; 16], d: i8) -> Self { + let negative = u64::from((d as u8) >> 7); + let magnitude = ((d ^ 0i8.wrapping_sub(negative as i8)).wrapping_add(negative as i8)) as u64; + let mut chosen = Self::identity(); + for (j, entry) in table.iter().enumerate() { + chosen = Self::select(&chosen, entry, 1 ^ crate::ct::nonzero(magnitude ^ (j as u64 + 1))); + } + let negated = Residue::zero().sub(&chosen.y); + chosen.y = Residue::select(&chosen.y, &negated, negative); + chosen + } + /// self + q for public points, branching on their values: no masked + /// doubling or selects. + fn add_public(&self, q: &Self) -> Self { + if self.z.is_zero() == 1 { + return *q; + } + if q.z.is_zero() == 1 { + return *self; + } + let z1 = self.z.square(); + let z2 = q.z.square(); + let u1 = self.x.mul(&z2); + let u2 = q.x.mul(&z1); + let s1 = self.y.mul(&q.z.mul(&z2)); + let s2 = q.y.mul(&self.z.mul(&z1)); + let h = u2.sub(&u1); + let r = s2.sub(&s1); + if h.is_zero() == 1 { + return if r.is_zero() == 1 { self.double() } else { Self::identity() }; + } + let h2 = h.square(); + let h3 = h.mul(&h2); + let u = u1.mul(&h2); + let x = r.square().sub(&h3).sub(&u.double()); + let y = r.mul(&u.sub(&x)).sub(&s1.mul(&h3)); + let z = h.mul(&self.z.mul(&q.z)); + Self { x, y, z } + } + fn neg(&self) -> Self { + Self { + x: self.x, + y: Residue::zero().sub(&self.y), + z: self.z, + } + } + fn odd_multiples(&self) -> OddMultiples { + let twice = self.double(); + let mut odd = [*self; 8]; + for j in 1..8 { + odd[j] = odd[j - 1].add_public(&twice); + } + odd + } + /// u1 G + u2 self in one doubling chain over both scalars' width-5 NAFs + /// (Straus's method), in variable time: every input is public. + fn joint(&self, u1: &PublicScalar, u2: &PublicScalar) -> Self { + let (g, q) = (Self::base().odd_multiples(), self.odd_multiples()); + let (d1, d2) = (u1.naf(5), u2.naf(5)); + let Some(top) = (0..64 * N + 1).rev().find(|&i| d1[i] != 0 || d2[i] != 0) else { + return Self::identity(); + }; let mut r = Self::identity(); - for i in (0..N * 64).rev() { - let doubled = r.double(); - let added = doubled.add(self); - r = Self::select(&doubled, &added, scalar.bit(i)) + for i in (0..=top).rev() { + r = r.double(); + for (d, table) in [(d1[i], &g), (d2[i], &q)] { + if d > 0 { + r = r.add_public(&table[d as usize / 2]); + } else if d < 0 { + r = r.add_public(&table[(-d) as usize / 2].neg()); + } + } } r } @@ -237,16 +361,30 @@ impl, const N: usize> Point { let s = Scalar::::from_bytes(&s[..N * 8]).map_err(|_| Error::InvalidSignature)?; let z = ModN::::from_wide(&from_be(digest)); let inverse = ModN::::from_wide(&s.0).invert(); - let u = Scalar(z.mul(&inverse).integer(), core::marker::PhantomData); - let v = Scalar(ModN::::from_wide(&r.0).mul(&inverse).integer(), core::marker::PhantomData); - let result = Self::base().mul(&u).add(&self.mul(&v)); - let (mut x, mut y) = ([0; N], [0; N]); - result.affine(&mut x, &mut y).map_err(|_| Error::InvalidSignature)?; - let x = ModN::::from_wide(&x).integer(); - if x == r.0 { - Ok(()) - } else { - Err(Error::InvalidSignature) + let u1 = PublicScalar(z.mul(&inverse).integer()); + let u2 = PublicScalar(ModN::::from_wide(&r.0).mul(&inverse).integer()); + let result = self.joint(&u1, &u2); + if result.z.is_zero() == 1 { + return Err(Error::InvalidSignature); + } + // x(R) = X / Z^2 is r mod n exactly when X = r' Z^2 for r' = r, or + // r' = r + n where that is still below p: no inversion. + let z2 = result.z.square(); + let r_field = Coord::::from_integer(&r.0).expect("r < n < p"); + if result.x.eq(&r_field.mul(&z2)) == 1 { + return Ok(()); + } + let mut r_plus_n = [0; N]; + let mut carry = false; + for i in 0..N { + let (sum, c1) = r.0[i].overflowing_add(C::Order::P[i]); + let (sum, c2) = sum.overflowing_add(u64::from(carry)); + r_plus_n[i] = sum; + carry = c1 | c2; + } + match (carry, Coord::::from_integer(&r_plus_n)) { + (false, Some(r_field)) if result.x.eq(&r_field.mul(&z2)) == 1 => Ok(()), + _ => Err(Error::InvalidSignature), } } } @@ -296,7 +434,12 @@ macro_rules! modulus { } }; } -modulus!(P256Field, 4, [0xffffffffffffffff, 0x00000000ffffffff, 0, 0xffffffff00000001]); +/// P-256's prime, with its special reduction step. +pub(crate) struct P256Field; +impl Modulus<4> for P256Field { + const P: [u64; 4] = [0xffffffffffffffff, 0x00000000ffffffff, 0, 0xffffffff00000001]; + const FORM: Form = Form::P256; +} modulus!(P256Order, 4, [0xf3b9cac2fc632551, 0xbce6faada7179e84, 0xffffffffffffffff, 0xffffffff00000000]); modulus!( P384Field, @@ -398,4 +541,119 @@ mod tests { let invalid = [n, [0; 32], [0xff; 32], n].concat(); assert!(matches!(Scalar::::first_valid(&invalid), Err(Error::Entropy))); } + + fn scalars(order: &[u64; N]) -> std::vec::Vec<[u64; N]> { + let mut out = std::vec![[0; N], [u64::MAX; N]]; + let mut n_minus_1 = *order; + n_minus_1[0] -= 1; + out.push(n_minus_1); + let mut one = [0; N]; + one[0] = 1; + out.push(one); + // Every window's six bits 0b011111 or 0b100000: digits +16 and -16. + out.push([0x7bdef7bdef7bdef7; N]); + out.push([0x8421084210842108; N]); + let mut seed = 0x0123456789abcdefu64; + for _ in 0..40 { + out.push(core::array::from_fn(|_| { + seed ^= seed << 13; + seed ^= seed >> 7; + seed ^= seed << 17; + seed + })); + } + out + } + + fn integer(words: &[u64; N]) -> num_bigint::BigInt { + num_bigint::BigInt::from_bytes_le( + num_bigint::Sign::Plus, + &words.iter().flat_map(|w| w.to_le_bytes()).collect::>(), + ) + } + + /// Booth digits lie in [-16, 16] and recompose the scalar, for both word + /// counts and at the patterns that give digits of +16 and -16. + fn check_booth, const N: usize>() { + use num_bigint::BigInt; + let mut extremes = [false; 2]; + for k in scalars(&C::Order::P) { + let scalar = Scalar::(k, core::marker::PhantomData); + let mut sum = BigInt::from(0); + for i in 0..windows(N) { + let d = scalar.booth(i); + assert!((-16..=16).contains(&d)); + extremes[0] |= d == 16; + extremes[1] |= d == -16; + sum += BigInt::from(d) << (5 * i); + } + assert_eq!(sum, integer(&k)); + } + assert_eq!(extremes, [true, true]); + } + + #[test] + fn booth_digits_recompose_the_scalar() { + check_booth::(); + check_booth::(); + } + + /// The width-w NAF of 384-bit scalars recomposes them. + #[test] + fn naf_digits_recompose_384_bit_scalars() { + use num_bigint::BigInt; + for k in scalars(&P384Order::P) { + let k = Coord::::from_wide(&k).integer(); + for w in [2u32, 5, 8] { + let digits = PublicScalar(k).naf(w); + let sum = digits.iter().enumerate().fold(BigInt::from(0), |s, (i, &d)| s + (BigInt::from(d) << i)); + assert_eq!(sum, integer(&k)); + } + } + } + + /// The fixed-schedule Booth multiplication equals the public joint + /// multiplication, a different recoding and addition, for scalars below + /// n including those with extreme digits. + fn check_mul, const N: usize>() { + let zero = PublicScalar([0; N]); + for k in scalars(&C::Order::P) { + let k = ModN::::from_wide(&k).integer(); + let ours = Point::::base().mul(&Scalar(k, core::marker::PhantomData)); + let theirs = Point::::identity().joint(&PublicScalar(k), &zero); + let (mut a, mut b) = ([0; N], [0; N]); + let (mut c, mut d) = ([0; N], [0; N]); + match (ours.affine(&mut a, &mut b), theirs.affine(&mut c, &mut d)) { + (Ok(()), Ok(())) => assert_eq!((a, b), (c, d)), + (Err(_), Err(_)) => assert_eq!(k, [0; N]), + _ => panic!("one side is the identity"), + } + } + } + + /// The public path's exceptional sums: equal summands double, and + /// opposite ones cancel to the identity. + fn check_exceptional_sums, const N: usize>() { + let mut one = [0; N]; + one[0] = 1; + let (mut a, mut b) = ([0; N], [0; N]); + let (mut c, mut d) = ([0; N], [0; N]); + let g = Point::::base(); + g.joint(&PublicScalar(one), &PublicScalar(one)).affine(&mut a, &mut b).unwrap(); + g.double().affine(&mut c, &mut d).unwrap(); + assert_eq!((a, b), (c, d)); + assert_eq!(g.add_public(&g.neg()).z.is_zero(), 1); + } + + #[test] + fn public_sums_of_equal_and_opposite_points() { + check_exceptional_sums::(); + check_exceptional_sums::(); + } + + #[test] + fn booth_multiplication_matches_the_public_path() { + check_mul::(); + check_mul::(); + } } diff --git a/lib/crypto/src/scalar25519.rs b/lib/crypto/src/scalar25519.rs index a1e08cd9b8..9f826c14a0 100644 --- a/lib/crypto/src/scalar25519.rs +++ b/lib/crypto/src/scalar25519.rs @@ -184,37 +184,10 @@ impl PublicScalar { bytes(&self.0) } - /// The width-w non-adjacent form, w in 2..=8: odd digits below 2^(w-1) in - /// magnitude, any two nonzero ones at least w positions apart (Hankerson, - /// Menezes and Vanstone, Guide to Elliptic Curve Cryptography, algorithm - /// 3.35). Variable time: the scalar is public. + /// The width-w NAF (`naf::naf`); k < L < 2^253, so 256 digits hold it. pub(crate) fn naf(&self, w: u32) -> [i8; 256] { - debug_assert!((2..=8).contains(&w)); - let mut k = self.0; let mut digits = [0i8; 256]; - let mut i = 0; - while k != [0; 4] { - if k[0] & 1 == 1 { - let window = (k[0] & ((1 << w) - 1)) as i32; - let d = if window >= 1 << (w - 1) { window - (1 << w) } else { window }; - digits[i] = d as i8; - // k - d, where |d| < 2^7 and k < 2^253 + 2^7: no word overflows. - if d >= 0 { - k[0] -= d as u64; - } else { - let (s, mut carry) = k[0].overflowing_add((-d) as u64); - k[0] = s; - for word in &mut k[1..] { - if !carry { - break; - } - (*word, carry) = word.overflowing_add(1); - } - } - } - k = [k[0] >> 1 | k[1] << 63, k[1] >> 1 | k[2] << 63, k[2] >> 1 | k[3] << 63, k[3] >> 1]; - i += 1; - } + crate::naf::naf(&self.0, w, &mut digits); digits } }