From 90b69c13eed20c8edf734756b8cb24e3ff721523 Mon Sep 17 00:00:00 2001 From: Isis Lovecruft Date: Mon, 16 Oct 2017 21:41:52 +0000 Subject: [PATCH] Revert "Merge remote-tracking branch 'floodyberry/optimized_scalar' into develop" MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit This reverts commit 804dab89247df02c90ed5ecfdc30be338a644a04, reversing changes made to 5d15ca77fffe2b9fb636fb67dee35aca05c82efc. This is due to a (previously undocumented) contract on the behaviours of (potentially unreduced mod \ell) "packed" scalars w.r.t. to the manner in which their bytes are interpreted. Upon documentation fixes and corresponding fixes being made on top of the floodyberry/optimized_scalar branch, this revert will again be reverted and then the additional changes merged (à la file:///usr/share/doc/git/html/howto/revert-a-faulty-merge.html). Signed-off-by: Isis Lovecruft --- src/constants_32bit.rs | 21 -- src/constants_64bit.rs | 13 - src/lib.rs | 7 - src/scalar.rs | 429 ++++++++++++++++++++++++++----- src/scalar_32bit.rs | 559 ----------------------------------------- src/scalar_64bit.rs | 475 ---------------------------------- 6 files changed, 361 insertions(+), 1143 deletions(-) delete mode 100644 src/scalar_32bit.rs delete mode 100644 src/scalar_64bit.rs diff --git a/src/constants_32bit.rs b/src/constants_32bit.rs index 004add7..bf46fe5 100644 --- a/src/constants_32bit.rs +++ b/src/constants_32bit.rs @@ -19,7 +19,6 @@ #![allow(non_snake_case)] use field_32bit::FieldElement32; -use scalar_32bit::Scalar32; use edwards::ExtendedPoint; use edwards::AffineNielsPoint; use edwards::EdwardsBasepointTable; @@ -97,26 +96,6 @@ pub const SQRT_MINUS_HALF: FieldElement32 = FieldElement32([ // sqrtMinusHalf -17256545, 3971863, 28865457, -1750208, 27359696, -16640980, 12573105, 1002827, -163343, 11073975, ]); -/// `L` is the order of base point, i.e. 2^252 + -/// 27742317777372353535851937790883648493 -pub const L: Scalar32 = Scalar32([ 0x1cf5d3ed, 0x009318d2, 0x1de73596, 0x1df3bd45, - 0x0000014d, 0x00000000, 0x00000000, 0x00000000, - 0x00100000 ]); - -/// `L` * `LFACTOR` = -1 (mod 2^29) -pub const LFACTOR: u32 = 0x12547e1b; - -/// `R` = R % L where R = 2^261 -pub const R: Scalar32 = Scalar32([ 0x114df9ed, 0x1a617303, 0x0f7c098c, 0x16793167, - 0x1ffd656e, 0x1fffffff, 0x1fffffff, 0x1fffffff, - 0x000fffff ]); - -/// `RR` = (R^2) % L where R = 2^261 -pub const RR: Scalar32 = Scalar32([ 0x0b5f9d12, 0x1e141b17, 0x158d7f3d, 0x143f3757, - 0x1972d781, 0x042feb7c, 0x1ceec73d, 0x1e184d1e, - 0x0005046d ]); - - /// Basepoint has y = 4/5. This is called `_POINT` to distinguish it from `_TABLE`, which should /// be used for scalar multiplication (it's much faster). pub const ED25519_BASEPOINT_POINT: ExtendedPoint = ExtendedPoint{ diff --git a/src/constants_64bit.rs b/src/constants_64bit.rs index e650757..ba7b2d9 100644 --- a/src/constants_64bit.rs +++ b/src/constants_64bit.rs @@ -19,7 +19,6 @@ #![allow(non_snake_case)] use field_64bit::FieldElement64; -use scalar_64bit::Scalar64; use edwards::ExtendedPoint; use edwards::AffineNielsPoint; use edwards::EdwardsBasepointTable; @@ -70,18 +69,6 @@ pub const SQRT_MINUS_APLUS2: FieldElement64 = FieldElement64([1693982333959686, /// `SQRT_MINUS_HALF` is sqrt(-1/2) pub const SQRT_MINUS_HALF: FieldElement64 = FieldElement64([266547196637087, 2134345371906993, 1135042577398223, 67298593331632, 743161882051057]); -/// `L` is the order of base point, i.e. 2^252 + 27742317777372353535851937790883648493 -pub const L: Scalar64 = Scalar64([ 0x0002631a5cf5d3ed, 0x000dea2f79cd6581, 0x000000000014def9, 0x0000000000000000, 0x0000100000000000 ]); - -/// `L` * `LFACTOR` = -1 (mod 2^51) -pub const LFACTOR: u64 = 0x51da312547e1b; - -/// `R` = R % L where R = 2^260 -pub const R: Scalar64 = Scalar64([ 0x000f48bd6721e6ed, 0x0003bab5ac67e45a, 0x000fffffeb35e51b, 0x000fffffffffffff, 0x00000fffffffffff ]); - -/// `RR` = (R^2) % L where R = 2^260 -pub const RR: Scalar64 = Scalar64([ 0x0009d265e952d13b, 0x000d63c715bea69f, 0x0005be65cb687604, 0x0003dceec73d217f, 0x000009411b7c309a ]); - /// Basepoint has y = 4/5. This is called `_POINT` to distinguish it from `_TABLE`, which should /// be used for scalar multiplication (it's much faster). pub const ED25519_BASEPOINT_POINT: ExtendedPoint = ExtendedPoint{ diff --git a/src/lib.rs b/src/lib.rs index be018f2..31c18ff 100644 --- a/src/lib.rs +++ b/src/lib.rs @@ -41,8 +41,6 @@ extern crate test; #[cfg(test)] extern crate sha2; -// this appears to only be used for serde support right now? -#[cfg(feature = "serde")] #[macro_use] extern crate arrayref; @@ -73,11 +71,6 @@ mod field_32bit; mod field_64bit; pub mod scalar; -#[cfg(not(feature="radix_51"))] -mod scalar_32bit; -#[cfg(feature="radix_51")] -mod scalar_64bit; - pub mod edwards; pub mod montgomery; diff --git a/src/scalar.rs b/src/scalar.rs index 399d00d..00124c1 100644 --- a/src/scalar.rs +++ b/src/scalar.rs @@ -22,11 +22,12 @@ //! //! The `Scalar` struct represents an element in ℤ/lℤ. //! -//! In contrast to `FieldElement`s, `Scalar`s are stored in +//! Arithmetic operations on `Scalar`s are done using 12 21-bit limbs. +//! However, in contrast to `FieldElement`s, `Scalar`s are stored in //! memory as bytes, allowing easy access to the bits of the `Scalar` //! when multiplying a point by a scalar. For efficient arithmetic -//! between two scalars, the `UnpackedScalar` struct (internally -//! either `Scalar32` or `Scalar64`) is stored as limbs. +//! between two scalars, the `UnpackedScalar` struct is stored as +//! limbs. use core::fmt::Debug; use core::ops::Neg; @@ -42,6 +43,9 @@ use rand::Rng; use digest::Digest; use generic_array::typenum::U64; +use constants; +use utils::{load3, load4}; + use subtle::slices_equal; use subtle::ConditionallyAssignable; use subtle::Equal; @@ -104,47 +108,50 @@ impl IndexMut for Scalar { impl<'b> MulAssign<&'b Scalar> for Scalar { fn mul_assign(&mut self, _rhs: &'b Scalar) { - *self = Scalar::mul(self, _rhs) + let result = (self as &Scalar) * _rhs; + self.0 = result.0; } } impl<'a, 'b> Mul<&'b Scalar> for &'a Scalar { type Output = Scalar; fn mul(self, _rhs: &'b Scalar) -> Scalar { - Scalar::mul(self, _rhs) + Scalar::multiply_add(self, _rhs, &Scalar::zero()) } } impl<'b> AddAssign<&'b Scalar> for Scalar { fn add_assign(&mut self, _rhs: &'b Scalar) { - *self = Scalar::add(self, _rhs); + *self = Scalar::multiply_add(&Scalar::one(), self, _rhs); } } impl<'a, 'b> Add<&'b Scalar> for &'a Scalar { type Output = Scalar; fn add(self, _rhs: &'b Scalar) -> Scalar { - Scalar::add(self, _rhs) + Scalar::multiply_add(&Scalar::one(), self, _rhs) } } impl<'b> SubAssign<&'b Scalar> for Scalar { fn sub_assign(&mut self, _rhs: &'b Scalar) { - *self = Scalar::sub(self, _rhs); + // (l-1)*_rhs + self = self - _rhs + *self = Scalar::multiply_add(&constants::l_minus_1, _rhs, self); } } impl<'a, 'b> Sub<&'b Scalar> for &'a Scalar { type Output = Scalar; fn sub(self, _rhs: &'b Scalar) -> Scalar { - Scalar::sub(self, _rhs) + // (l-1)*_rhs + self = self - _rhs + Scalar::multiply_add(&constants::l_minus_1, _rhs, self) } } impl<'a> Neg for &'a Scalar { type Output = Scalar; fn neg(self) -> Scalar { - Scalar::sub(&Scalar::zero(), self) + self * &constants::l_minus_1 } } @@ -227,18 +234,6 @@ impl<'de> Deserialize<'de> for Scalar { } } -/// An `UnpackedScalar` represents an element of the field GF(l), optimized for speed. -#[cfg(feature="radix_51")] -type UnpackedScalar = Scalar64; -#[cfg(feature="radix_51")] -use scalar_64bit::*; - -/// An `UnpackedScalar` represents an element of the field GF(l), optimized for speed. -#[cfg(not(feature="radix_51"))] -type UnpackedScalar = Scalar32; -#[cfg(not(feature="radix_51"))] -use scalar_32bit::*; - impl Scalar { /// Return a `Scalar` chosen uniformly at random using a user-provided RNG. /// @@ -389,6 +384,26 @@ impl Scalar { naf } + // Unpack a scalar into 12 21-bit limbs. + fn unpack(&self) -> UnpackedScalar { + let mask_21bits: i64 = (1 << 21) - 1; + let mut a = UnpackedScalar([0i64; 12]); + a[ 0] = mask_21bits & load3(&self.0[ 0..]) ; + a[ 1] = mask_21bits & (load4(&self.0[ 2..]) >> 5); + a[ 2] = mask_21bits & (load3(&self.0[ 5..]) >> 2); + a[ 3] = mask_21bits & (load4(&self.0[ 7..]) >> 7); + a[ 4] = mask_21bits & (load4(&self.0[10..]) >> 4); + a[ 5] = mask_21bits & (load3(&self.0[13..]) >> 1); + a[ 6] = mask_21bits & (load4(&self.0[15..]) >> 6); + a[ 7] = mask_21bits & (load3(&self.0[18..]) >> 3); + a[ 8] = mask_21bits & load3(&self.0[21..]) ; + a[ 9] = mask_21bits & (load4(&self.0[23..]) >> 5); + a[10] = mask_21bits & (load3(&self.0[26..]) >> 2); + a[11] = load4(&self.0[28..]) >> 7 ; + + a + } + /// Write this scalar in radix 16, with coefficients in `[-8,8)`, /// i.e., compute `a_i` such that /// @@ -427,41 +442,127 @@ impl Scalar { output } - /// Unpack this `Scalar` to an `UnpackedScalar` - pub fn unpack(&self) -> UnpackedScalar { - UnpackedScalar::from_bytes(&self.0) - } - - /// Compute `a + b` (mod l) - pub fn add(a: &Scalar, b: &Scalar) -> Scalar { - UnpackedScalar::add(&a.unpack(), &b.unpack()).pack() - } - - /// Compute `a - b` (mod l). - pub fn sub(a: &Scalar, b: &Scalar) -> Scalar { - UnpackedScalar::sub(&a.unpack(), &b.unpack()).pack() - } - - /// Compute `a * b` (mod l). - pub fn mul(a: &Scalar, b: &Scalar) -> Scalar { - UnpackedScalar::mul(&a.unpack(), &b.unpack()).pack() - } - - /// Compute `(a * b) + c` (mod l). + /// Compute `ab+c (mod l)`. + /// XXX should this exist, or should we just have Mul, Add etc impls + /// that unpack and then call UnpackedScalar::multiply_add ? pub fn multiply_add(a: &Scalar, b: &Scalar, c: &Scalar) -> Scalar { - UnpackedScalar::add(&UnpackedScalar::mul(&a.unpack(), &b.unpack()), &c.unpack()).pack() + // Unpack scalars into limbs + let al = a.unpack(); + let bl = b.unpack(); + let cl = c.unpack(); + + // Multiply and repack + UnpackedScalar::multiply_add(&al, &bl, &cl).pack() } /// Reduce a 512-bit little endian number mod l pub fn reduce(input: &[u8; 64]) -> Scalar { - UnpackedScalar::from_bytes_wide(input).pack() + let mut s = [0i64; 24]; + + // XXX express this as two unpack_limbs + // some issues re: masking with the top byte of the 32byte input + let mask_21bits: i64 = (1 << 21) -1; + s[0] = mask_21bits & load3(&input[ 0..]) ; + s[1] = mask_21bits & (load4(&input[ 2..]) >> 5); + s[2] = mask_21bits & (load3(&input[ 5..]) >> 2); + s[3] = mask_21bits & (load4(&input[ 7..]) >> 7); + s[4] = mask_21bits & (load4(&input[10..]) >> 4); + s[5] = mask_21bits & (load3(&input[13..]) >> 1); + s[6] = mask_21bits & (load4(&input[15..]) >> 6); + s[7] = mask_21bits & (load3(&input[18..]) >> 3); + s[8] = mask_21bits & load3(&input[21..]) ; + s[9] = mask_21bits & (load4(&input[23..]) >> 5); + s[10] = mask_21bits & (load3(&input[26..]) >> 2); + s[11] = mask_21bits & (load4(&input[28..]) >> 7); + s[12] = mask_21bits & (load4(&input[31..]) >> 4); + s[13] = mask_21bits & (load3(&input[34..]) >> 1); + s[14] = mask_21bits & (load4(&input[36..]) >> 6); + s[15] = mask_21bits & (load3(&input[39..]) >> 3); + s[16] = mask_21bits & load3(&input[42..]) ; + s[17] = mask_21bits & (load4(&input[44..]) >> 5); + s[18] = mask_21bits & (load3(&input[47..]) >> 2); + s[19] = mask_21bits & (load4(&input[49..]) >> 7); + s[20] = mask_21bits & (load4(&input[52..]) >> 4); + s[21] = mask_21bits & (load3(&input[55..]) >> 1); + s[22] = mask_21bits & (load4(&input[57..]) >> 6); + s[23] = load4(&input[60..]) >> 3 ; + + // XXX replacing the previous code in this function with the + // call to reduce_limbs adds two extra carry passes (the ones + // at the top of the reduce_limbs function). Otherwise they + // are identical. The test seems to work OK but it would be + // good to check that this really is OK to add. + UnpackedScalar::reduce_limbs(&mut s).pack() + } +} + +/// The `UnpackedScalar` struct represents an element in ℤ/lℤ as 12 +/// 21-bit limbs. +#[derive(Copy,Clone)] +pub struct UnpackedScalar(pub [i64; 12]); + +impl Index for UnpackedScalar { + type Output = i64; + + fn index(&self, _index: usize) -> &i64 { + &(self.0[_index]) + } +} + +impl IndexMut for UnpackedScalar { + fn index_mut(&mut self, _index: usize) -> &mut i64 { + &mut (self.0[_index]) } } impl UnpackedScalar { /// Pack the limbs of this `UnpackedScalar` into a `Scalar`. fn pack(&self) -> Scalar { - Scalar(self.to_bytes()) + let mut s = Scalar::zero(); + s[0] = (self.0[ 0] >> 0) as u8; + s[1] = (self.0[ 0] >> 8) as u8; + s[2] = ((self.0[ 0] >> 16) | (self.0[ 1] << 5)) as u8; + s[3] = (self.0[ 1] >> 3) as u8; + s[4] = (self.0[ 1] >> 11) as u8; + s[5] = ((self.0[ 1] >> 19) | (self.0[ 2] << 2)) as u8; + s[6] = (self.0[ 2] >> 6) as u8; + s[7] = ((self.0[ 2] >> 14) | (self.0[ 3] << 7)) as u8; + s[8] = (self.0[ 3] >> 1) as u8; + s[9] = (self.0[ 3] >> 9) as u8; + s[10] = ((self.0[ 3] >> 17) | (self.0[ 4] << 4)) as u8; + s[11] = (self.0[ 4] >> 4) as u8; + s[12] = (self.0[ 4] >> 12) as u8; + s[13] = ((self.0[ 4] >> 20) | (self.0[ 5] << 1)) as u8; + s[14] = (self.0[ 5] >> 7) as u8; + s[15] = ((self.0[ 5] >> 15) | (self.0[ 6] << 6)) as u8; + s[16] = (self.0[ 6] >> 2) as u8; + s[17] = (self.0[ 6] >> 10) as u8; + s[18] = ((self.0[ 6] >> 18) | (self.0[ 7] << 3)) as u8; + s[19] = (self.0[ 7] >> 5) as u8; + s[20] = (self.0[ 7] >> 13) as u8; + s[21] = (self.0[ 8] >> 0) as u8; + s[22] = (self.0[ 8] >> 8) as u8; + s[23] = ((self.0[ 8] >> 16) | (self.0[ 9] << 5)) as u8; + s[24] = (self.0[ 9] >> 3) as u8; + s[25] = (self.0[ 9] >> 11) as u8; + s[26] = ((self.0[ 9] >> 19) | (self.0[10] << 2)) as u8; + s[27] = (self.0[10] >> 6) as u8; + s[28] = ((self.0[10] >> 14) | (self.0[11] << 7)) as u8; + s[29] = (self.0[11] >> 1) as u8; + s[30] = (self.0[11] >> 9) as u8; + s[31] = (self.0[11] >> 17) as u8; + + s + } + + /// Return the zero scalar. + pub fn zero() -> UnpackedScalar { + UnpackedScalar([0,0,0,0,0,0,0,0,0,0,0,0]) + } + + /// Return the one scalar. + pub fn one() -> UnpackedScalar { + UnpackedScalar([1,0,0,0,0,0,0,0,0,0,0,0]) } /// Compute the multiplicative inverse of this scalar. @@ -470,25 +571,25 @@ impl UnpackedScalar { // https://briansmith.org/ecc-inversion-addition-chains-01#curve25519_scalar_inversion // as it was published on 2017-09-03. - let _1 = self.to_montgomery(); - let _10 = _1.montgomery_square(); - let _100 = _10.montgomery_square(); - let _11 = UnpackedScalar::montgomery_mul(&_10, &_1); - let _101 = UnpackedScalar::montgomery_mul(&_10, &_11); - let _111 = UnpackedScalar::montgomery_mul(&_10, &_101); - let _1001 = UnpackedScalar::montgomery_mul(&_10, &_111); - let _1011 = UnpackedScalar::montgomery_mul(&_10, &_1001); - let _1111 = UnpackedScalar::montgomery_mul(&_100, &_1011); + let _1 = *self; + let _10 = _1.square(); + let _100 = _10.square(); + let _11 = UnpackedScalar::multiply_add(&_10, &_1, &UnpackedScalar::zero()); + let _101 = UnpackedScalar::multiply_add(&_10, &_11, &UnpackedScalar::zero()); + let _111 = UnpackedScalar::multiply_add(&_10, &_101, &UnpackedScalar::zero()); + let _1001 = UnpackedScalar::multiply_add(&_10, &_111, &UnpackedScalar::zero()); + let _1011 = UnpackedScalar::multiply_add(&_10, &_1001, &UnpackedScalar::zero()); + let _1111 = UnpackedScalar::multiply_add(&_100, &_1011, &UnpackedScalar::zero()); // _10000 - let mut y = UnpackedScalar::montgomery_mul(&_1111, &_1); + let mut y = UnpackedScalar::multiply_add(&_1111, &_1, &UnpackedScalar::zero()); #[inline] fn square_multiply(y: &mut UnpackedScalar, squarings: usize, x: &UnpackedScalar) { for _ in 0..squarings { - *y = y.montgomery_square(); + *y = y.square(); } - *y = UnpackedScalar::montgomery_mul(y, x); + *y = UnpackedScalar::multiply_add(y, x, &UnpackedScalar::zero()); } square_multiply(&mut y, 123 + 3, &_101); @@ -519,14 +620,194 @@ impl UnpackedScalar { square_multiply(&mut y, 3, &_101); square_multiply(&mut y, 1 + 2, &_11); - y.from_montgomery() + y + } + + /// Compute `a^2 (mod l)`. + pub fn square(&self) -> UnpackedScalar { + let a = self.0; + let mut result = [0i64; 24]; + + result[0] = a[0]*a[0]; + result[1] = 2i64 * a[0]*a[1]; + result[2] = 2i64 * (a[0]*a[2]) + a[1]*a[1]; + result[3] = 2i64 * (a[0]*a[3] + a[1]*a[2]); + result[4] = 2i64 * (a[0]*a[4] + a[1]*a[3]) + a[2]*a[2]; + result[5] = 2i64 * (a[0]*a[5] + a[1]*a[4] + a[2]*a[3]); + result[6] = 2i64 * (a[0]*a[6] + a[1]*a[5] + a[2]*a[4]) + a[3]*a[3]; + result[7] = 2i64 * (a[0]*a[7] + a[1]*a[6] + a[2]*a[5] + a[3]*a[4]); + result[8] = 2i64 * (a[0]*a[8] + a[1]*a[7] + a[2]*a[6] + a[3]*a[5]) + a[4]*a[4]; + result[9] = 2i64 * (a[0]*a[9] + a[1]*a[8] + a[2]*a[7] + a[3]*a[6] + a[4]*a[5]); + result[10] = 2i64 * (a[0]*a[10] + a[1]*a[9] + a[2]*a[8] + a[3]*a[7] + a[4]*a[6]) + a[5]*a[5]; + result[11] = 2i64 * (a[0]*a[11] + a[1]*a[10] + a[2]*a[9] + a[3]*a[8] + a[4]*a[7] + a[5]*a[6]); + result[12] = 2i64 * (a[1]*a[11] + a[2]*a[10] + a[3]*a[9] + a[4]*a[8] + a[5]*a[7]) + a[6]*a[6]; + result[13] = 2i64 * (a[2]*a[11] + a[3]*a[10] + a[4]*a[9] + a[5]*a[8] + a[6]*a[7]); + result[14] = 2i64 * (a[3]*a[11] + a[4]*a[10] + a[5]*a[9] + a[6]*a[8]) + a[7]*a[7]; + result[15] = 2i64 * (a[4]*a[11] + a[5]*a[10] + a[6]*a[9] + a[7]*a[8]); + result[16] = 2i64 * (a[5]*a[11] + a[6]*a[10] + a[7]*a[9]) + a[8]*a[8]; + result[17] = 2i64 * (a[6]*a[11] + a[7]*a[10] + a[8]*a[9]); + result[18] = 2i64 * (a[7]*a[11] + a[8]*a[10]) + a[9]*a[9]; + result[19] = 2i64 * (a[8]*a[11] + a[9]*a[10]); + result[20] = 2i64 * (a[9]*a[11]) + a[10]*a[10]; + result[21] = 2i64 * (a[10]*a[11]); + result[22] = a[11]*a[11]; + result[23] = 0i64; + + // Reduce limbs + UnpackedScalar::reduce_limbs(&mut result) + } + + /// Compute `ab+c (mod l)`. + pub fn multiply_add(a: &UnpackedScalar, + b: &UnpackedScalar, + c: &UnpackedScalar) -> UnpackedScalar { + let mut result = [0i64; 24]; + + // Multiply a and b, and add c + result[0] = c[0] + a[0]*b[0]; + result[1] = c[1] + a[0]*b[1] + a[1]*b[0]; + result[2] = c[2] + a[0]*b[2] + a[1]*b[1] + a[2]*b[0]; + result[3] = c[3] + a[0]*b[3] + a[1]*b[2] + a[2]*b[1] + a[3]*b[0]; + result[4] = c[4] + a[0]*b[4] + a[1]*b[3] + a[2]*b[2] + a[3]*b[1] + a[4]*b[0]; + result[5] = c[5] + a[0]*b[5] + a[1]*b[4] + a[2]*b[3] + a[3]*b[2] + a[4]*b[1] + a[5]*b[0]; + result[6] = c[6] + a[0]*b[6] + a[1]*b[5] + a[2]*b[4] + a[3]*b[3] + a[4]*b[2] + a[5]*b[1] + a[6]*b[0]; + result[7] = c[7] + a[0]*b[7] + a[1]*b[6] + a[2]*b[5] + a[3]*b[4] + a[4]*b[3] + a[5]*b[2] + a[6]*b[1] + a[7]*b[0]; + result[8] = c[8] + a[0]*b[8] + a[1]*b[7] + a[2]*b[6] + a[3]*b[5] + a[4]*b[4] + a[5]*b[3] + a[6]*b[2] + a[7]*b[1] + a[8]*b[0]; + result[9] = c[9] + a[0]*b[9] + a[1]*b[8] + a[2]*b[7] + a[3]*b[6] + a[4]*b[5] + a[5]*b[4] + a[6]*b[3] + a[7]*b[2] + a[8]*b[1] + a[9]*b[0]; + result[10] = c[10] + a[0]*b[10] + a[1]*b[9] + a[2]*b[8] + a[3]*b[7] + a[4]*b[6] + a[5]*b[5] + a[6]*b[4] + a[7]*b[3] + a[8]*b[2] + a[9]*b[1] + a[10]*b[0]; + result[11] = c[11] + a[0]*b[11] + a[1]*b[10] + a[2]*b[9] + a[3]*b[8] + a[4]*b[7] + a[5]*b[6] + a[6]*b[5] + a[7]*b[4] + a[8]*b[3] + a[9]*b[2] + a[10]*b[1] + a[11]*b[0]; + result[12] = a[1]*b[11] + a[2]*b[10] + a[3]*b[9] + a[4]*b[8] + a[5]*b[7] + a[6]*b[6] + a[7]*b[5] + a[8]*b[4] + a[9]*b[3] + a[10]*b[2] + a[11]*b[1]; + result[13] = a[2]*b[11] + a[3]*b[10] + a[4]*b[9] + a[5]*b[8] + a[6]*b[7] + a[7]*b[6] + a[8]*b[5] + a[9]*b[4] + a[10]*b[3] + a[11]*b[2]; + result[14] = a[3]*b[11] + a[4]*b[10] + a[5]*b[9] + a[6]*b[8] + a[7]*b[7] + a[8]*b[6] + a[9]*b[5] + a[10]*b[4] + a[11]*b[3]; + result[15] = a[4]*b[11] + a[5]*b[10] + a[6]*b[9] + a[7]*b[8] + a[8]*b[7] + a[9]*b[6] + a[10]*b[5] + a[11]*b[4]; + result[16] = a[5]*b[11] + a[6]*b[10] + a[7]*b[9] + a[8]*b[8] + a[9]*b[7] + a[10]*b[6] + a[11]*b[5]; + result[17] = a[6]*b[11] + a[7]*b[10] + a[8]*b[9] + a[9]*b[8] + a[10]*b[7] + a[11]*b[6]; + result[18] = a[7]*b[11] + a[8]*b[10] + a[9]*b[9] + a[10]*b[8] + a[11]*b[7]; + result[19] = a[8]*b[11] + a[9]*b[10] + a[10]*b[9] + a[11]*b[8]; + result[20] = a[9]*b[11] + a[10]*b[10] + a[11]*b[9]; + result[21] = a[10]*b[11] + a[11]*b[10]; + result[22] = a[11]*b[11]; + result[23] = 0i64; + + // Reduce limbs + UnpackedScalar::reduce_limbs(&mut result) + } + + /// Reduce 24 limbs to 12, consuming the input. Reduction is mod + /// + /// l = 2^252 + 27742317777372353535851937790883648493, + /// + /// so + /// + /// 2^252 = -27742317777372353535851937790883648493 (mod l). + /// + /// We can write the right-hand side in 21-bit limbs as + /// + /// rhs = 666643 * 2^0 + /// + 470296 * 2^21 + /// + 654183 * 2^42 + /// - 997805 * 2^63 + /// + 136657 * 2^84 + /// - 683901 * 2^105 + /// + /// The (12+k)-th limb of `limbs` is the coefficient of + /// + /// 2^(252 + 21*k) + /// + /// since 12*21 = 252. By the above, we have that + /// + /// c * 2^(252 + 21*k) = c * 666643 * 2^(21*k) + /// + c * 470296 * 2^(42*k) + ... + /// + /// so we can eliminate it by adding those values to the lower + /// limbs. Reduction mod l amounts to eliminating all of the + /// high limbs while carrying as appropriate to prevent + /// overflows in the lower limbs. + fn reduce_limbs(mut limbs: &mut [i64; 24]) -> UnpackedScalar { + #[inline] + #[allow(dead_code)] + fn do_reduction(limbs: &mut [i64; 24], i: usize) { + limbs[i - 12] += limbs[i] * 666643; + limbs[i - 11] += limbs[i] * 470296; + limbs[i - 10] += limbs[i] * 654183; + limbs[i - 9] -= limbs[i] * 997805; + limbs[i - 8] += limbs[i] * 136657; + limbs[i - 7] -= limbs[i] * 683901; + limbs[i] = 0; + } + /// Carry excess from the `i`-th limb into the `(i+1)`-th limb. + /// Postcondition: `0 <= limbs[i] < 2^21`. + #[inline] + #[allow(dead_code)] + fn do_carry_uncentered(limbs: &mut [i64; 24], i: usize) { + let carry: i64 = limbs[i] >> 21; + limbs[i+1] += carry; + limbs[i ] -= carry << 21; + } + #[inline] + #[allow(dead_code)] + /// Carry excess from the `i`-th limb into the `(i+1)`-th limb. + /// Postcondition: `-2^20 <= limbs[i] < 2^20`. + fn do_carry_centered(limbs: &mut [i64; 24], i: usize) { + let carry: i64 = (limbs[i] + (1<<20)) >> 21; + limbs[i+1] += carry; + limbs[i ] -= carry << 21; + } + + for i in 0..23 { + do_carry_centered(&mut limbs, i); + } + for i in (0..23).filter(|x| x % 2 == 1) { + do_carry_centered(&mut limbs, i); + } + + do_reduction(&mut limbs, 23); + do_reduction(&mut limbs, 22); + do_reduction(&mut limbs, 21); + do_reduction(&mut limbs, 20); + do_reduction(&mut limbs, 19); + do_reduction(&mut limbs, 18); + + for i in (6..18).filter(|x| x % 2 == 0) { + do_carry_centered(&mut limbs, i); + } + for i in (6..16).filter(|x| x % 2 == 1) { + do_carry_centered(&mut limbs, i); + } + + do_reduction(&mut limbs, 17); + do_reduction(&mut limbs, 16); + do_reduction(&mut limbs, 15); + do_reduction(&mut limbs, 14); + do_reduction(&mut limbs, 13); + do_reduction(&mut limbs, 12); + + for i in (0..12).filter(|x| x % 2 == 0) { + do_carry_centered(&mut limbs, i); + } + for i in (0..12).filter(|x| x % 2 == 1) { + do_carry_centered(&mut limbs, i); + } + + do_reduction(&mut limbs, 12); + + for i in 0..12 { + do_carry_uncentered(&mut limbs, i); + } + + do_reduction(&mut limbs, 12); + + for i in 0..11 { + do_carry_uncentered(&mut limbs, i); + } + + UnpackedScalar(*array_ref!(limbs, 0, 12)) } } #[cfg(test)] mod test { use super::*; - use constants; /// x = 2238329342913194256032495932344128051776374960164957527413114840482143558222 pub static X: Scalar = Scalar( @@ -534,12 +815,6 @@ mod test { 0x59, 0x13, 0xb4, 0x64, 0x1b, 0xc2, 0x7d, 0x52, 0x52, 0xa5, 0x85, 0x10, 0x1b, 0xcc, 0x42, 0x44, 0xd4, 0x49, 0xf4, 0xa8, 0x79, 0xd9, 0xf2, 0x04]); - /// 1/x = 6859937278830797291664592131120606308688036382723378951768035303146619657244 - pub static XINV: Scalar = Scalar( - [0x1c, 0xdc, 0x17, 0xfc, 0xe0, 0xe9, 0xa5, 0xbb, - 0xd9, 0x24, 0x7e, 0x56, 0xbb, 0x01, 0x63, 0x47, - 0xbb, 0xba, 0x31, 0xed, 0xd5, 0xa9, 0xbb, 0x96, - 0xd5, 0x0b, 0xcd, 0x7a, 0x3f, 0x96, 0x2a, 0x0f]); /// y = 2592331292931086675770238855846338635550719849568364935475441891787804997264 pub static Y: Scalar = Scalar( [0x90, 0x76, 0x33, 0xfe, 0x1c, 0x4b, 0x66, 0xa4, @@ -677,7 +952,6 @@ mod test { #[test] fn invert() { let inv_X = X.invert(); - assert_eq!(inv_X, XINV); let should_be_one = &inv_X * &X; assert_eq!(should_be_one, Scalar::one()); } @@ -710,7 +984,7 @@ mod bench { use test::Bencher; use super::*; - use super::test::{X}; + use super::test::{X, Y, Z}; #[bench] fn scalar_random(b: &mut Bencher) { @@ -719,9 +993,28 @@ mod bench { b.iter(|| Scalar::random(&mut csprng)); } + #[bench] + fn scalar_multiply_add(b: &mut Bencher) { + b.iter(|| Scalar::multiply_add(&X, &Y, &Z)); + } + #[bench] fn invert(b: &mut Bencher) { let x = X.unpack(); b.iter(|| x.invert()); } + + #[bench] + fn square(b: &mut Bencher) { + let x = X.unpack(); + b.iter(|| x.square()); + } + + #[bench] + fn scalar_unpacked_multiply_add(b: &mut Bencher) { + let x = X.unpack(); + let y = Y.unpack(); + let z = Z.unpack(); + b.iter(|| UnpackedScalar::multiply_add(&x, &y, &z)); + } } diff --git a/src/scalar_32bit.rs b/src/scalar_32bit.rs deleted file mode 100644 index b0d8b57..0000000 --- a/src/scalar_32bit.rs +++ /dev/null @@ -1,559 +0,0 @@ -//! Arithmetic mod 2^252 + 27742317777372353535851937790883648493 -//! with 9 29-bit unsigned limbs -//! -//! To see that this is safe for intermediate results, note that -//! the largest limb in a 9 by 9 product of 29-bit limbs will be -//! (0x1fffffff^2) * 9 = 0x23fffffdc0000009 (62 bits). -//! -//! For a one level Karatsuba decomposition, the specific ranges -//! depend on how the limbs are combined, but will stay within -//! -0x1ffffffe00000008 (62 bits with sign bit) to -//! 0x43fffffbc0000011 (63 bits), which is still safe. -//! -//! (the 9th limb will never exceed 21 bits, so the actual -//! ranges are slightly smaller) - -use core::fmt::Debug; -use core::ops::{Index, IndexMut}; - -use constants; - -/// The `Scalar32` struct represents an element in ℤ/lℤ as 9 29-bit limbs -#[derive(Copy,Clone)] -pub struct Scalar32(pub [u32; 9]); - -impl Debug for Scalar32 { - fn fmt(&self, f: &mut ::core::fmt::Formatter) -> ::core::fmt::Result { - write!(f, "Scalar32: {:?}", &self.0[..]) - } -} - -impl Index for Scalar32 { - type Output = u32; - fn index(&self, _index: usize) -> &u32 { - &(self.0[_index]) - } -} - -impl IndexMut for Scalar32 { - fn index_mut(&mut self, _index: usize) -> &mut u32 { - &mut (self.0[_index]) - } -} - -/// u32 * u32 = u64 multiply helper -#[inline(always)] -fn m(x: u32, y: u32) -> u64 { - (x as u64) * (y as u64) -} - -impl Scalar32 { - /// Return the zero scalar. - pub fn zero() -> Scalar32 { - Scalar32([0,0,0,0,0,0,0,0,0]) - } - - /// Unpack a 32 byte / 512 bit scalar into 9 29-bit limbs, ignoring the upper 3 bits. - pub fn from_bytes(bytes: &[u8; 32]) -> Scalar32 { - let mut words = [0u32; 8]; - for i in 0..8 { - for j in 0..4 { - words[i] |= (bytes[(i * 4) + j] as u32) << (j * 8); - } - } - - let mask = (1u32 << 29) - 1; - let top_mask = (1u32 << 21) - 1; - let mut s = Scalar32::zero(); - - s[ 0] = words[0] & mask; - s[ 1] = ((words[0] >> 29) | (words[1] << 3)) & mask; - s[ 2] = ((words[1] >> 26) | (words[2] << 6)) & mask; - s[ 3] = ((words[2] >> 23) | (words[3] << 9)) & mask; - s[ 4] = ((words[3] >> 20) | (words[4] << 12)) & mask; - s[ 5] = ((words[4] >> 17) | (words[5] << 15)) & mask; - s[ 6] = ((words[5] >> 14) | (words[6] << 18)) & mask; - s[ 7] = ((words[6] >> 11) | (words[7] << 21)) & mask; - s[ 8] = (words[7] >> 8) & top_mask; - - s - } - - /// Reduce a 64 byte / 512 bit scalar mod l. - pub fn from_bytes_wide(bytes: &[u8; 64]) -> Scalar32 { - let mut words = [0u32; 16]; - for i in 0..16 { - for j in 0..4 { - words[i] |= (bytes[(i * 4) + j] as u32) << (j * 8); - } - } - - let mask = (1u32 << 29) - 1; - let mut lo = Scalar32::zero(); - let mut hi = Scalar32::zero(); - - lo[0] = words[ 0] & mask; - lo[1] = ((words[ 0] >> 29) | (words[ 1] << 3)) & mask; - lo[2] = ((words[ 1] >> 26) | (words[ 2] << 6)) & mask; - lo[3] = ((words[ 2] >> 23) | (words[ 3] << 9)) & mask; - lo[4] = ((words[ 3] >> 20) | (words[ 4] << 12)) & mask; - lo[5] = ((words[ 4] >> 17) | (words[ 5] << 15)) & mask; - lo[6] = ((words[ 5] >> 14) | (words[ 6] << 18)) & mask; - lo[7] = ((words[ 6] >> 11) | (words[ 7] << 21)) & mask; - lo[8] = ((words[ 7] >> 8) | (words[ 8] << 24)) & mask; - hi[0] = ((words[ 8] >> 5) | (words[ 9] << 27)) & mask; - hi[1] = (words[ 9] >> 2) & mask; - hi[2] = ((words[ 9] >> 31) | (words[10] << 1)) & mask; - hi[3] = ((words[10] >> 28) | (words[11] << 4)) & mask; - hi[4] = ((words[11] >> 25) | (words[12] << 7)) & mask; - hi[5] = ((words[12] >> 22) | (words[13] << 10)) & mask; - hi[6] = ((words[13] >> 19) | (words[14] << 13)) & mask; - hi[7] = ((words[14] >> 16) | (words[15] << 16)) & mask; - hi[8] = (words[15] >> 13) & mask; - - lo = Scalar32::montgomery_mul(&lo, &constants::R); // (lo * R) / R = lo - hi = Scalar32::montgomery_mul(&hi, &constants::RR); // (hi * R^2) / R = hi * R - - Scalar32::add(&hi, &lo) // (hi * R) + lo - } - - /// Pack the limbs of this `Scalar32` into 32 bytes. - pub fn to_bytes(&self) -> [u8; 32] { - let mut s = [0u8; 32]; - - s[0] = (self.0[ 0] >> 0) as u8; - s[1] = (self.0[ 0] >> 8) as u8; - s[2] = (self.0[ 0] >> 16) as u8; - s[3] = ((self.0[ 0] >> 24) | (self.0[ 1] << 5)) as u8; - s[4] = (self.0[ 1] >> 3) as u8; - s[5] = (self.0[ 1] >> 11) as u8; - s[6] = (self.0[ 1] >> 19) as u8; - s[7] = ((self.0[ 1] >> 27) | (self.0[ 2] << 2)) as u8; - s[8] = (self.0[ 2] >> 6) as u8; - s[9] = (self.0[ 2] >> 14) as u8; - s[10] = ((self.0[ 2] >> 22) | (self.0[ 3] << 7)) as u8; - s[11] = (self.0[ 3] >> 1) as u8; - s[12] = (self.0[ 3] >> 9) as u8; - s[13] = (self.0[ 3] >> 17) as u8; - s[14] = ((self.0[ 3] >> 25) | (self.0[ 4] << 4)) as u8; - s[15] = (self.0[ 4] >> 4) as u8; - s[16] = (self.0[ 4] >> 12) as u8; - s[17] = (self.0[ 4] >> 20) as u8; - s[18] = ((self.0[ 4] >> 28) | (self.0[ 5] << 1)) as u8; - s[19] = (self.0[ 5] >> 7) as u8; - s[20] = (self.0[ 5] >> 15) as u8; - s[21] = ((self.0[ 5] >> 23) | (self.0[ 6] << 6)) as u8; - s[22] = (self.0[ 6] >> 2) as u8; - s[23] = (self.0[ 6] >> 10) as u8; - s[24] = (self.0[ 6] >> 18) as u8; - s[25] = ((self.0[ 6] >> 26) | (self.0[ 7] << 3)) as u8; - s[26] = (self.0[ 7] >> 5) as u8; - s[27] = (self.0[ 7] >> 13) as u8; - s[28] = (self.0[ 7] >> 21) as u8; - s[29] = (self.0[ 8] >> 0) as u8; - s[30] = (self.0[ 8] >> 8) as u8; - s[31] = (self.0[ 8] >> 16) as u8; - - s - } - - /// Compute `a + b` (mod l). - pub fn add(a: &Scalar32, b: &Scalar32) -> Scalar32 { - let mut sum = Scalar32::zero(); - let mask = (1u32 << 29) - 1; - - // a + b - let mut carry: u32 = 0; - for i in 0..9 { - carry = a[i] + b[i] + (carry >> 29); - sum[i] = carry & mask; - } - - // subtract l if the sum is >= l - Scalar32::sub(&sum, &constants::L) - } - - /// Compute `a - b` (mod l). - pub fn sub(a: &Scalar32, b: &Scalar32) -> Scalar32 { - let mut difference = Scalar32::zero(); - let mask = (1u32 << 29) - 1; - - // a - b - let mut borrow: u32 = 0; - for i in 0..9 { - borrow = a[i].wrapping_sub(b[i] + (borrow >> 31)); - difference[i] = borrow & mask; - } - - // conditionally add l if the difference is negative - let underflow_mask = ((borrow >> 31) ^ 1).wrapping_sub(1); - let mut carry: u32 = 0; - for i in 0..9 { - carry = (carry >> 29) + difference[i] + (constants::L[i] & underflow_mask); - difference[i] = carry & mask; - } - - difference - } - - /// Compute `a * b`. - /// - /// This is implemented with a one-level refined Karatsuba decomposition - #[inline(always)] - fn mul_internal(a: &Scalar32, b: &Scalar32) -> [u64; 17] { - let mut z = [0u64; 17]; - - z[0] = m(a[0],b[0]); // c00 - z[1] = m(a[0],b[1]) + m(a[1],b[0]); // c01 - z[2] = m(a[0],b[2]) + m(a[1],b[1]) + m(a[2],b[0]); // c02 - z[3] = m(a[0],b[3]) + m(a[1],b[2]) + m(a[2],b[1]) + m(a[3],b[0]); // c03 - z[4] = m(a[0],b[4]) + m(a[1],b[3]) + m(a[2],b[2]) + m(a[3],b[1]) + m(a[4],b[0]); // c04 - z[5] = m(a[1],b[4]) + m(a[2],b[3]) + m(a[3],b[2]) + m(a[4],b[1]); // c05 - z[6] = m(a[2],b[4]) + m(a[3],b[3]) + m(a[4],b[2]); // c06 - z[7] = m(a[3],b[4]) + m(a[4],b[3]); // c07 - z[8] = (m(a[4],b[4])).wrapping_sub(z[3]); // c08 - c03 - - z[10] = z[5].wrapping_sub(m(a[5],b[5])); // c05mc10 - z[11] = z[6].wrapping_sub(m(a[5],b[6]) + m(a[6],b[5])); // c06mc11 - z[12] = z[7].wrapping_sub(m(a[5],b[7]) + m(a[6],b[6]) + m(a[7],b[5])); // c07mc12 - z[13] = m(a[5],b[8]) + m(a[6],b[7]) + m(a[7],b[6]) + m(a[8],b[5]); // c13 - z[14] = m(a[6],b[8]) + m(a[7],b[7]) + m(a[8],b[6]); // c14 - z[15] = m(a[7],b[8]) + m(a[8],b[7]); // c15 - z[16] = m(a[8],b[8]); // c16 - - z[ 5] = z[10].wrapping_sub(z[ 0]); // c05mc10 - c00 - z[ 6] = z[11].wrapping_sub(z[ 1]); // c06mc11 - c01 - z[ 7] = z[12].wrapping_sub(z[ 2]); // c07mc12 - c02 - z[ 8] = z[ 8].wrapping_sub(z[13]); // c08mc13 - c03 - z[ 9] = z[14].wrapping_add(z[ 4]); // c14 + c04 - z[10] = z[15].wrapping_add(z[10]); // c15 + c05mc10 - z[11] = z[16].wrapping_add(z[11]); // c16 + c06mc11 - - let aa = [ - a[0]+a[5], - a[1]+a[6], - a[2]+a[7], - a[3]+a[8] - ]; - - let bb = [ - b[0]+b[5], - b[1]+b[6], - b[2]+b[7], - b[3]+b[8] - ]; - - z[ 5] = (m(aa[0],bb[0])) .wrapping_add(z[ 5]); // c20 + c05mc10 - c00 - z[ 6] = (m(aa[0],bb[1]) + m(aa[1],bb[0])) .wrapping_add(z[ 6]); // c21 + c06mc11 - c01 - z[ 7] = (m(aa[0],bb[2]) + m(aa[1],bb[1]) + m(aa[2],bb[0])) .wrapping_add(z[ 7]); // c22 + c07mc12 - c02 - z[ 8] = (m(aa[0],bb[3]) + m(aa[1],bb[2]) + m(aa[2],bb[1]) + m(aa[3],bb[0])) .wrapping_add(z[ 8]); // c23 + c08mc13 - c03 - z[ 9] = (m(aa[0], b[4]) + m(aa[1],bb[3]) + m(aa[2],bb[2]) + m(aa[3],bb[1]) + m(a[4],bb[0])).wrapping_sub(z[ 9]); // c24 - c14 - c04 - z[10] = ( m(aa[1], b[4]) + m(aa[2],bb[3]) + m(aa[3],bb[2]) + m(a[4],bb[1])).wrapping_sub(z[10]); // c25 - c15 - c05mc10 - z[11] = ( m(aa[2], b[4]) + m(aa[3],bb[3]) + m(a[4],bb[2])).wrapping_sub(z[11]); // c26 - c16 - c06mc11 - z[12] = ( m(aa[3], b[4]) + m(a[4],bb[3])).wrapping_sub(z[12]); // c27 - c07mc12 - - z - } - - /// Compute `a^2`. - #[inline(always)] - fn square_internal(a: &Scalar32) -> [u64; 17] { - let aa = [ - a[0]*2, - a[1]*2, - a[2]*2, - a[3]*2, - a[4]*2, - a[5]*2, - a[6]*2, - a[7]*2 - ]; - - [ - m( a[0],a[0]), - m(aa[0],a[1]), - m(aa[0],a[2]) + m( a[1],a[1]), - m(aa[0],a[3]) + m(aa[1],a[2]), - m(aa[0],a[4]) + m(aa[1],a[3]) + m( a[2],a[2]), - m(aa[0],a[5]) + m(aa[1],a[4]) + m(aa[2],a[3]), - m(aa[0],a[6]) + m(aa[1],a[5]) + m(aa[2],a[4]) + m( a[3],a[3]), - m(aa[0],a[7]) + m(aa[1],a[6]) + m(aa[2],a[5]) + m(aa[3],a[4]), - m(aa[0],a[8]) + m(aa[1],a[7]) + m(aa[2],a[6]) + m(aa[3],a[5]) + m( a[4],a[4]), - m(aa[1],a[8]) + m(aa[2],a[7]) + m(aa[3],a[6]) + m(aa[4],a[5]), - m(aa[2],a[8]) + m(aa[3],a[7]) + m(aa[4],a[6]) + m( a[5],a[5]), - m(aa[3],a[8]) + m(aa[4],a[7]) + m(aa[5],a[6]), - m(aa[4],a[8]) + m(aa[5],a[7]) + m( a[6],a[6]), - m(aa[5],a[8]) + m(aa[6],a[7]), - m(aa[6],a[8]) + m( a[7],a[7]), - m(aa[7],a[8]), - m( a[8],a[8]), - ] - } - - /// Compute `limbs/R` (mod l), where R is the Montgomery modulus 2^261 - #[inline(always)] - fn montgomery_reduce(limbs: &[u64; 17]) -> Scalar32 { - - #[inline(always)] - fn part1(sum: u64) -> (u64, u32) { - let p = (sum as u32).wrapping_mul(constants::LFACTOR) & ((1u32 << 29) - 1); - ((sum + m(p,constants::L[0])) >> 29, p) - } - - #[inline(always)] - fn part2(sum: u64) -> (u64, u32) { - let w = (sum as u32) & ((1u32 << 29) - 1); - (sum >> 29, w) - } - - // note: l5,l6,l7 are zero, so their multiplies can be skipped - let l = &constants::L; - - // the first half computes the Montgomery adjustment factor n, and begins adding n*l to make limbs divisible by R - let (carry, n0) = part1( limbs[ 0]); - let (carry, n1) = part1(carry + limbs[ 1] + m(n0,l[1])); - let (carry, n2) = part1(carry + limbs[ 2] + m(n0,l[2]) + m(n1,l[1])); - let (carry, n3) = part1(carry + limbs[ 3] + m(n0,l[3]) + m(n1,l[2]) + m(n2,l[1])); - let (carry, n4) = part1(carry + limbs[ 4] + m(n0,l[4]) + m(n1,l[3]) + m(n2,l[2]) + m(n3,l[1])); - let (carry, n5) = part1(carry + limbs[ 5] + m(n1,l[4]) + m(n2,l[3]) + m(n3,l[2]) + m(n4,l[1])); - let (carry, n6) = part1(carry + limbs[ 6] + m(n2,l[4]) + m(n3,l[3]) + m(n4,l[2]) + m(n5,l[1])); - let (carry, n7) = part1(carry + limbs[ 7] + m(n3,l[4]) + m(n4,l[3]) + m(n5,l[2]) + m(n6,l[1])); - let (carry, n8) = part1(carry + limbs[ 8] + m(n0,l[8]) + m(n4,l[4]) + m(n5,l[3]) + m(n6,l[2]) + m(n7,l[1])); - - // limbs is divisible by R now, so we can divide by R by simply storing the upper half as the result - let (carry, r0) = part2(carry + limbs[ 9] + m(n1,l[8]) + m(n5,l[4]) + m(n6,l[3]) + m(n7,l[2]) + m(n8,l[1])); - let (carry, r1) = part2(carry + limbs[10] + m(n2,l[8]) + m(n6,l[4]) + m(n7,l[3]) + m(n8,l[2])); - let (carry, r2) = part2(carry + limbs[11] + m(n3,l[8]) + m(n7,l[4]) + m(n8,l[3])); - let (carry, r3) = part2(carry + limbs[12] + m(n4,l[8]) + m(n8,l[4])); - let (carry, r4) = part2(carry + limbs[13] + m(n5,l[8]) ); - let (carry, r5) = part2(carry + limbs[14] + m(n6,l[8]) ); - let (carry, r6) = part2(carry + limbs[15] + m(n7,l[8]) ); - let (carry, r7) = part2(carry + limbs[16] + m(n8,l[8])); - let r8 = carry as u32; - - // result may be >= l, so attempt to subtract l - Scalar32::sub(&Scalar32([r0,r1,r2,r3,r4,r5,r6,r7,r8]), l) - } - - /// Compute `a * b` (mod l). - #[inline(never)] - pub fn mul(a: &Scalar32, b: &Scalar32) -> Scalar32 { - let ab = Scalar32::montgomery_reduce(&Scalar32::mul_internal(a, b)); - Scalar32::montgomery_reduce(&Scalar32::mul_internal(&ab, &constants::RR)) - } - - /// Compute `a^2` (mod l). - #[inline(never)] - pub fn square(&self) -> Scalar32 { - let aa = Scalar32::montgomery_reduce(&Scalar32::square_internal(self)); - Scalar32::montgomery_reduce(&Scalar32::mul_internal(&aa, &constants::RR)) - } - - /// Compute `(a * b) / R` (mod l), where R is the Montgomery modulus 2^261 - #[inline(never)] - pub fn montgomery_mul(a: &Scalar32, b: &Scalar32) -> Scalar32 { - Scalar32::montgomery_reduce(&Scalar32::mul_internal(a, b)) - } - - /// Compute `(a^2) / R` (mod l) in Montgomery form, where R is the Montgomery modulus 2^261 - #[inline(never)] - pub fn montgomery_square(&self) -> Scalar32 { - Scalar32::montgomery_reduce(&Scalar32::square_internal(self)) - } - - /// Puts a Scalar32 in to Montgomery form, i.e. computes `a*R (mod l)` - #[inline(never)] - pub fn to_montgomery(&self) -> Scalar32 { - Scalar32::montgomery_mul(self, &constants::RR) - } - - /// Takes a Scalar32 out of Montgomery form, i.e. computes `a/R (mod l)` - pub fn from_montgomery(&self) -> Scalar32 { - let mut limbs = [0u64; 17]; - for i in 0..9 { - limbs[i] = self[i] as u64; - } - Scalar32::montgomery_reduce(&limbs) - } -} - - -#[cfg(test)] -mod test { - use super::*; - - /// Note: x is 2^253-1 which is slightly larger than the largest scalar produced by - /// this implementation (l-1), and should verify there are no overflows for valid scalars - /// - /// x = 2^253-1 = 14474011154664524427946373126085988481658748083205070504932198000989141204991 - /// x = 7237005577332262213973186563042994240801631723825162898930247062703686954002 mod l - /// x = 5147078182513738803124273553712992179887200054963030844803268920753008712037*R mod l in Montgomery form - pub static X: Scalar32 = Scalar32( - [0x1fffffff, 0x1fffffff, 0x1fffffff, 0x1fffffff, - 0x1fffffff, 0x1fffffff, 0x1fffffff, 0x1fffffff, - 0x001fffff]); - - /// x^2 = 3078544782642840487852506753550082162405942681916160040940637093560259278169 mod l - pub static XX: Scalar32 = Scalar32( - [0x00217559, 0x000b3401, 0x103ff43b, 0x1462a62c, - 0x1d6f9f38, 0x18e7a42f, 0x09a3dcee, 0x008dbe18, - 0x0006ce65]); - - /// x^2 = 2912514428060642753613814151688322857484807845836623976981729207238463947987*R mod l in Montgomery form - pub static XX_MONT: Scalar32 = Scalar32( - [0x152b4d2e, 0x0571d53b, 0x1da6d964, 0x188663b6, - 0x1d1b5f92, 0x19d50e3f, 0x12306c29, 0x0c6f26fe, - 0x00030edb]); - - /// y = 6145104759870991071742105800796537629880401874866217824609283457819451087098 - pub static Y: Scalar32 = Scalar32( - [0x1e1458fa, 0x165ba838, 0x1d787b36, 0x0e577f3a, - 0x1d2baf06, 0x1d689a19, 0x1fff3047, 0x117704ab, - 0x000d9601]); - - /// x*y = 36752150652102274958925982391442301741 - pub static XY: Scalar32 = Scalar32( - [0x0ba7632d, 0x017736bb, 0x15c76138, 0x0c69daa1, - 0x000001ba, 0x00000000, 0x00000000, 0x00000000, - 0x00000000]); - - /// x*y = 3783114862749659543382438697751927473898937741870308063443170013240655651591*R mod l in Montgomery form - pub static XY_MONT: Scalar32 = Scalar32( - [0x077b51e1, 0x1c64e119, 0x02a19ef5, 0x18d2129e, - 0x00de0430, 0x045a7bc8, 0x04cfc7c9, 0x1c002681, - 0x000bdc1c]); - - /// a = 2351415481556538453565687241199399922945659411799870114962672658845158063753 - pub static A: Scalar32 = Scalar32( - [0x07b3be89, 0x02291b60, 0x14a99f03, 0x07dc3787, - 0x0a782aae, 0x16262525, 0x0cfdb93f, 0x13f5718d, - 0x000532da]); - - /// b = 4885590095775723760407499321843594317911456947580037491039278279440296187236 - pub static B: Scalar32 = Scalar32( - [0x15421564, 0x1e69fd72, 0x093d9692, 0x161785be, - 0x1587d69f, 0x09d9dada, 0x130246c0, 0x0c0a8e72, - 0x000acd25]); - - /// a+b = 0 - /// a-b = 4702830963113076907131374482398799845891318823599740229925345317690316127506 - pub static AB: Scalar32 = Scalar32( - [0x0f677d12, 0x045236c0, 0x09533e06, 0x0fb86f0f, - 0x14f0555c, 0x0c4c4a4a, 0x19fb727f, 0x07eae31a, - 0x000a65b5]); - - // c = (2^512 - 1) % l = 1627715501170711445284395025044413883736156588369414752970002579683115011840 - pub static C: Scalar32 = Scalar32( - [0x049c0f00, 0x00308f1a, 0x0164d1e9, 0x1c374ed1, - 0x1be65d00, 0x19e90bfa, 0x08f73bb1, 0x036f8613, - 0x00039941]); - - #[test] - fn mul_max() { - let res = Scalar32::mul(&X, &X); - for i in 0..9 { - assert!(res[i] == XX[i]); - } - } - - #[test] - fn square_max() { - let res = X.square(); - for i in 0..9 { - assert!(res[i] == XX[i]); - } - } - - #[test] - fn montgomery_mul_max() { - let res = Scalar32::montgomery_mul(&X, &X); - for i in 0..9 { - assert!(res[i] == XX_MONT[i]); - } - } - - #[test] - fn montgomery_square_max() { - let res = X.montgomery_square(); - for i in 0..9 { - assert!(res[i] == XX_MONT[i]); - } - } - - #[test] - fn mul() { - let res = Scalar32::mul(&X, &Y); - for i in 0..9 { - assert!(res[i] == XY[i]); - } - } - - #[test] - fn montgomery_mul() { - let res = Scalar32::montgomery_mul(&X, &Y); - for i in 0..9 { - assert!(res[i] == XY_MONT[i]); - } - } - - #[test] - fn add() { - let res = Scalar32::add(&A, &B); - let zero = Scalar32::zero(); - for i in 0..9 { - assert!(res[i] == zero[i]); - } - } - - #[test] - fn sub() { - let res = Scalar32::sub(&A, &B); - for i in 0..9 { - assert!(res[i] == AB[i]); - } - } - - #[test] - fn from_bytes_wide() { - let bignum = [255u8; 64]; // 2^512 - 1 - let reduced = Scalar32::from_bytes_wide(&bignum); - for i in 0..9 { - assert!(reduced[i] == C[i]); - } - } -} - - -#[cfg(all(test, feature = "bench"))] -mod bench { - use test::Bencher; - - use super::*; - use super::test::{X, Y}; - - #[bench] - fn square(b: &mut Bencher) { - b.iter(|| X.square()); - } - - #[bench] - fn mul(b: &mut Bencher) { - b.iter(|| Scalar32::mul(&X, &Y)); - } - - #[bench] - fn montgomery_square(b: &mut Bencher) { - b.iter(|| X.montgomery_square()); - } - - #[bench] - fn montgomery_mul(b: &mut Bencher) { - b.iter(|| Scalar32::montgomery_mul(&X, &Y)); - } - - #[bench] - fn from_bytes_wide(b: &mut Bencher) { - let bignum = [255u8; 64]; // 2^512 - 1 - b.iter(|| Scalar32::from_bytes_wide(&bignum)); - } -} diff --git a/src/scalar_64bit.rs b/src/scalar_64bit.rs deleted file mode 100644 index 7076b56..0000000 --- a/src/scalar_64bit.rs +++ /dev/null @@ -1,475 +0,0 @@ -//! Arithmetic mod 2^252 + 27742317777372353535851937790883648493 -//! with 5 52-bit unsigned limbs. 51-bit limbs would cover the -//! desired bit range (253 bits), but isn't large enough to reduce -//! a 512 bit number with Montgomery multiplication, so 52 bits is -//! used instead -//! -//! To see that this is safe for intermediate results, note that -//! the largest limb in a 5 by 5 product of 52-bit limbs will be -//! (0xfffffffffffff^2) * 5 = 0x4ffffffffffff60000000000005 (107 bits). -//! -//! (the 5th limb will never exceed 45 bits, so the actual -//! ranges are slightly smaller) - - -use core::fmt::Debug; -use core::ops::{Index, IndexMut}; - -use constants; - -/// The `Scalar64` struct represents an element in ℤ/lℤ as 5 52-bit limbs -#[derive(Copy,Clone)] -pub struct Scalar64(pub [u64; 5]); - -impl Debug for Scalar64 { - fn fmt(&self, f: &mut ::core::fmt::Formatter) -> ::core::fmt::Result { - write!(f, "Scalar64: {:?}", &self.0[..]) - } -} - -impl Index for Scalar64 { - type Output = u64; - fn index(&self, _index: usize) -> &u64 { - &(self.0[_index]) - } -} - -impl IndexMut for Scalar64 { - fn index_mut(&mut self, _index: usize) -> &mut u64 { - &mut (self.0[_index]) - } -} - -/// u64 * u64 = u128 multiply helper -#[inline(always)] -fn m(x: u64, y: u64) -> u128 { - (x as u128) * (y as u128) -} - -impl Scalar64 { - /// Return the zero scalar - pub fn zero() -> Scalar64 { - Scalar64([0,0,0,0,0]) - } - - /// Unpack a 32 byte / 256 bit scalar into 5 52-bit limbs, ignoring the upper 3 bits - pub fn from_bytes(bytes: &[u8; 32]) -> Scalar64 { - let mut words = [0u64; 8]; - for i in 0..4 { - for j in 0..8 { - words[i] |= (bytes[(i * 8) + j] as u64) << (j * 8); - } - } - - let mask = (1u64 << 52) - 1; - let top_mask = (1u64 << 45) - 1; - let mut s = Scalar64::zero(); - - s[ 0] = words[0] & mask; - s[ 1] = ((words[0] >> 52) | (words[1] << 12)) & mask; - s[ 2] = ((words[1] >> 40) | (words[2] << 24)) & mask; - s[ 3] = ((words[2] >> 28) | (words[3] << 36)) & mask; - s[ 4] = (words[3] >> 16) & top_mask; - - s - } - - /// Reduce a 64 byte / 512 bit scalar mod l - pub fn from_bytes_wide(bytes: &[u8; 64]) -> Scalar64 { - let mut words = [064; 16]; - for i in 0..8 { - for j in 0..8 { - words[i] |= (bytes[(i * 8) + j] as u64) << (j * 8); - } - } - - let mask = (1u64 << 52) - 1; - let mut lo = Scalar64::zero(); - let mut hi = Scalar64::zero(); - - lo[0] = words[ 0] & mask; - lo[1] = ((words[ 0] >> 52) | (words[ 1] << 12)) & mask; - lo[2] = ((words[ 1] >> 40) | (words[ 2] << 24)) & mask; - lo[3] = ((words[ 2] >> 28) | (words[ 3] << 36)) & mask; - lo[4] = ((words[ 3] >> 16) | (words[ 4] << 48)) & mask; - hi[0] = (words[ 4] >> 4) & mask; - hi[1] = ((words[ 4] >> 56) | (words[ 5] << 8)) & mask; - hi[2] = ((words[ 5] >> 44) | (words[ 6] << 20)) & mask; - hi[3] = ((words[ 6] >> 32) | (words[ 7] << 32)) & mask; - hi[4] = words[ 7] >> 20 ; - - lo = Scalar64::montgomery_mul(&lo, &constants::R); // (lo * R) / R = lo - hi = Scalar64::montgomery_mul(&hi, &constants::RR); // (hi * R^2) / R = hi * R - - Scalar64::add(&hi, &lo) - } - - /// Pack the limbs of this `Scalar64` into 32 bytes - pub fn to_bytes(&self) -> [u8; 32] { - let mut s = [0u8; 32]; - - s[0] = (self.0[ 0] >> 0) as u8; - s[1] = (self.0[ 0] >> 8) as u8; - s[2] = (self.0[ 0] >> 16) as u8; - s[3] = (self.0[ 0] >> 24) as u8; - s[4] = (self.0[ 0] >> 32) as u8; - s[5] = (self.0[ 0] >> 40) as u8; - s[6] = ((self.0[ 0] >> 48) | (self.0[ 1] << 4)) as u8; - s[7] = (self.0[ 1] >> 4) as u8; - s[8] = (self.0[ 1] >> 12) as u8; - s[9] = (self.0[ 1] >> 20) as u8; - s[10] = (self.0[ 1] >> 28) as u8; - s[11] = (self.0[ 1] >> 36) as u8; - s[12] = (self.0[ 1] >> 44) as u8; - s[13] = (self.0[ 2] >> 0) as u8; - s[14] = (self.0[ 2] >> 8) as u8; - s[15] = (self.0[ 2] >> 16) as u8; - s[16] = (self.0[ 2] >> 24) as u8; - s[17] = (self.0[ 2] >> 32) as u8; - s[18] = (self.0[ 2] >> 40) as u8; - s[19] = ((self.0[ 2] >> 48) | (self.0[ 3] << 4)) as u8; - s[20] = (self.0[ 3] >> 4) as u8; - s[21] = (self.0[ 3] >> 12) as u8; - s[22] = (self.0[ 3] >> 20) as u8; - s[23] = (self.0[ 3] >> 28) as u8; - s[24] = (self.0[ 3] >> 36) as u8; - s[25] = (self.0[ 3] >> 44) as u8; - s[26] = (self.0[ 4] >> 0) as u8; - s[27] = (self.0[ 4] >> 8) as u8; - s[28] = (self.0[ 4] >> 16) as u8; - s[29] = (self.0[ 4] >> 24) as u8; - s[30] = (self.0[ 4] >> 32) as u8; - s[31] = (self.0[ 4] >> 40) as u8; - - s - } - - /// Compute `a + b` (mod l) - pub fn add(a: &Scalar64, b: &Scalar64) -> Scalar64 { - let mut sum = Scalar64::zero(); - let mask = (1u64 << 52) - 1; - - // a + b - let mut carry: u64 = 0; - for i in 0..5 { - carry = a[i] + b[i] + (carry >> 52); - sum[i] = carry & mask; - } - - // subtract l if the sum is >= l - Scalar64::sub(&sum, &constants::L) - } - - /// Compute `a - b` (mod l) - pub fn sub(a: &Scalar64, b: &Scalar64) -> Scalar64 { - let mut difference = Scalar64::zero(); - let mask = (1u64 << 52) - 1; - - // a - b - let mut borrow: u64 = 0; - for i in 0..5 { - borrow = a[i].wrapping_sub(b[i] + (borrow >> 63)); - difference[i] = borrow & mask; - } - - // conditionally add l if the difference is negative - let underflow_mask = ((borrow >> 63) ^ 1).wrapping_sub(1); - let mut carry: u64 = 0; - for i in 0..5 { - carry = (carry >> 52) + difference[i] + (constants::L[i] & underflow_mask); - difference[i] = carry & mask; - } - - difference - } - - /// Compute `a * b` - #[inline(always)] - fn mul_internal(a: &Scalar64, b: &Scalar64) -> [u128; 9] { - [ - m(a[0],b[0]), - m(a[0],b[1]) + m(a[1],b[0]), - m(a[0],b[2]) + m(a[1],b[1]) + m(a[2],b[0]), - m(a[0],b[3]) + m(a[1],b[2]) + m(a[2],b[1]) + m(a[3],b[0]), - m(a[0],b[4]) + m(a[1],b[3]) + m(a[2],b[2]) + m(a[3],b[1]) + m(a[4],b[0]), - m(a[1],b[4]) + m(a[2],b[3]) + m(a[3],b[2]) + m(a[4],b[1]), - m(a[2],b[4]) + m(a[3],b[3]) + m(a[4],b[2]), - m(a[3],b[4]) + m(a[4],b[3]), - m(a[4],b[4]) - ] - } - - /// Compute `a^2` - #[inline(always)] - fn square_internal(a: &Scalar64) -> [u128; 9] { - let aa = [ - a[0]*2, - a[1]*2, - a[2]*2, - a[3]*2, - ]; - - [ - m( a[0],a[0]), - m(aa[0],a[1]), - m(aa[0],a[2]) + m( a[1],a[1]), - m(aa[0],a[3]) + m(aa[1],a[2]), - m(aa[0],a[4]) + m(aa[1],a[3]) + m( a[2],a[2]), - m(aa[1],a[4]) + m(aa[2],a[3]), - m(aa[2],a[4]) + m( a[3],a[3]), - m(aa[3],a[4]), - m(a[4],a[4]) - ] - } - - /// Compute `limbs/R` (mod l), where R is the Montgomery modulus 2^260 - #[inline(always)] - fn montgomery_reduce(limbs: &[u128; 9]) -> Scalar64 { - - #[inline(always)] - fn part1(sum: u128) -> (u128, u64) { - let p = (sum as u64).wrapping_mul(constants::LFACTOR) & ((1u64 << 52) - 1); - ((sum + m(p,constants::L[0])) >> 52, p) - } - - #[inline(always)] - fn part2(sum: u128) -> (u128, u64) { - let w = (sum as u64) & ((1u64 << 52) - 1); - (sum >> 52, w) - } - - // note: l3 is zero, so its multiplies can be skipped - let l = &constants::L; - - // the first half computes the Montgomery adjustment factor n, and begins adding n*l to make limbs divisible by R - let (carry, n0) = part1( limbs[0]); - let (carry, n1) = part1(carry + limbs[1] + m(n0,l[1])); - let (carry, n2) = part1(carry + limbs[2] + m(n0,l[2]) + m(n1,l[1])); - let (carry, n3) = part1(carry + limbs[3] + m(n1,l[2]) + m(n2,l[1])); - let (carry, n4) = part1(carry + limbs[4] + m(n0,l[4]) + m(n2,l[2]) + m(n3,l[1])); - - // limbs is divisible by R now, so we can divide by R by simply storing the upper half as the result - let (carry, r0) = part2(carry + limbs[5] + m(n1,l[4]) + m(n3,l[2]) + m(n4,l[1])); - let (carry, r1) = part2(carry + limbs[6] + m(n2,l[4]) + m(n4,l[2])); - let (carry, r2) = part2(carry + limbs[7] + m(n3,l[4]) ); - let (carry, r3) = part2(carry + limbs[8] + m(n4,l[4])); - let r4 = carry as u64; - - // result may be >= l, so attempt to subtract l - Scalar64::sub(&Scalar64([r0,r1,r2,r3,r4]), l) - } - - /// Compute `a * b` (mod l) - #[inline(never)] - pub fn mul(a: &Scalar64, b: &Scalar64) -> Scalar64 { - let ab = Scalar64::montgomery_reduce(&Scalar64::mul_internal(a, b)); - Scalar64::montgomery_reduce(&Scalar64::mul_internal(&ab, &constants::RR)) - } - - /// Compute `a^2` (mod l) - #[inline(never)] - pub fn square(&self) -> Scalar64 { - let aa = Scalar64::montgomery_reduce(&Scalar64::square_internal(self)); - Scalar64::montgomery_reduce(&Scalar64::mul_internal(&aa, &constants::RR)) - } - - /// Compute `(a * b) / R` (mod l), where R is the Montgomery modulus 2^260 - #[inline(never)] - pub fn montgomery_mul(a: &Scalar64, b: &Scalar64) -> Scalar64 { - Scalar64::montgomery_reduce(&Scalar64::mul_internal(a, b)) - } - - /// Compute `(a^2) / R` (mod l) in Montgomery form, where R is the Montgomery modulus 2^260 - #[inline(never)] - pub fn montgomery_square(&self) -> Scalar64 { - Scalar64::montgomery_reduce(&Scalar64::square_internal(self)) - } - - /// Puts a Scalar64 in to Montgomery form, i.e. computes `a*R (mod l)` - #[inline(never)] - pub fn to_montgomery(&self) -> Scalar64 { - Scalar64::montgomery_mul(self, &constants::RR) - } - - /// Takes a Scalar64 out of Montgomery form, i.e. computes `a/R (mod l)` - #[inline(never)] - pub fn from_montgomery(&self) -> Scalar64 { - let mut limbs = [0u128; 9]; - for i in 0..5 { - limbs[i] = self[i] as u128; - } - Scalar64::montgomery_reduce(&limbs) - } -} - - -#[cfg(test)] -mod test { - use super::*; - - /// Note: x is 2^253-1 which is slightly larger than the largest scalar produced by - /// this implementation (l-1), and should show there are no overflows for valid scalars - /// - /// x = 14474011154664524427946373126085988481658748083205070504932198000989141204991 - /// x = 7237005577332262213973186563042994240801631723825162898930247062703686954002 mod l - /// x = 3057150787695215392275360544382990118917283750546154083604586903220563173085*R mod l in Montgomery form - pub static X: Scalar64 = Scalar64( - [0x000fffffffffffff, 0x000fffffffffffff, 0x000fffffffffffff, 0x000fffffffffffff, - 0x00001fffffffffff]); - - /// x^2 = 3078544782642840487852506753550082162405942681916160040940637093560259278169 mod l - pub static XX: Scalar64 = Scalar64( - [0x0001668020217559, 0x000531640ffd0ec0, 0x00085fd6f9f38a31, 0x000c268f73bb1cf4, - 0x000006ce65046df0]); - - /// x^2 = 4413052134910308800482070043710297189082115023966588301924965890668401540959*R mod l in Montgomery form - pub static XX_MONT: Scalar64 = Scalar64( - [0x000c754eea569a5c, 0x00063b6ed36cb215, 0x0008ffa36bf25886, 0x000e9183614e7543, - 0x0000061db6c6f26f]); - - /// y = 6145104759870991071742105800796537629880401874866217824609283457819451087098 - pub static Y: Scalar64 = Scalar64( - [0x000b75071e1458fa, 0x000bf9d75e1ecdac, 0x000433d2baf0672b, 0x0005fffcc11fad13, - 0x00000d96018bb825]); - - /// x*y = 36752150652102274958925982391442301741 mod l - pub static XY: Scalar64 = Scalar64( - [0x000ee6d76ba7632d, 0x000ed50d71d84e02, 0x00000000001ba634, 0x0000000000000000, - 0x0000000000000000]); - - /// x*y = 658448296334113745583381664921721413881518248721417041768778176391714104386*R mod l in Montgomery form - pub static XY_MONT: Scalar64 = Scalar64( - [0x0006d52bf200cfd5, 0x00033fb1d7021570, 0x000f201bc07139d8, 0x0001267e3e49169e, - 0x000007b839c00268]); - - /// a = 2351415481556538453565687241199399922945659411799870114962672658845158063753 - pub static A: Scalar64 = Scalar64( - [0x0005236c07b3be89, 0x0001bc3d2a67c0c4, 0x000a4aa782aae3ee, 0x0006b3f6e4fec4c4, - 0x00000532da9fab8c]); - - /// b = 4885590095775723760407499321843594317911456947580037491039278279440296187236 - pub static B: Scalar64 = Scalar64( - [0x000d3fae55421564, 0x000c2df24f65a4bc, 0x0005b5587d69fb0b, 0x00094c091b013b3b, - 0x00000acd25605473]); - - /// a+b = 0 - /// a-b = 4702830963113076907131374482398799845891318823599740229925345317690316127506 - pub static AB: Scalar64 = Scalar64( - [0x000a46d80f677d12, 0x0003787a54cf8188, 0x0004954f0555c7dc, 0x000d67edc9fd8989, - 0x00000a65b53f5718]); - - // c = (2^512 - 1) % l = 1627715501170711445284395025044413883736156588369414752970002579683115011840 - pub static C: Scalar64 = Scalar64( - [0x000611e3449c0f00, 0x000a768859347a40, 0x0007f5be65d00e1b, 0x0009a3dceec73d21, - 0x00000399411b7c30]); - - #[test] - fn mul_max() { - let res = Scalar64::mul(&X, &X); - for i in 0..5 { - assert!(res[i] == XX[i]); - } - } - - #[test] - fn square_max() { - let res = X.square(); - for i in 0..5 { - assert!(res[i] == XX[i]); - } - } - - #[test] - fn montgomery_mul_max() { - let res = Scalar64::montgomery_mul(&X, &X); - for i in 0..5 { - assert!(res[i] == XX_MONT[i]); - } - } - - #[test] - fn montgomery_square_max() { - let res = X.montgomery_square(); - for i in 0..5 { - assert!(res[i] == XX_MONT[i]); - } - } - - #[test] - fn mul() { - let res = Scalar64::mul(&X, &Y); - for i in 0..5 { - assert!(res[i] == XY[i]); - } - } - - #[test] - fn montgomery_mul() { - let res = Scalar64::montgomery_mul(&X, &Y); - for i in 0..5 { - assert!(res[i] == XY_MONT[i]); - } - } - - #[test] - fn add() { - let res = Scalar64::add(&A, &B); - let zero = Scalar64::zero(); - for i in 0..5 { - assert!(res[i] == zero[i]); - } - } - - #[test] - fn sub() { - let res = Scalar64::sub(&A, &B); - for i in 0..5 { - assert!(res[i] == AB[i]); - } - } - - #[test] - fn from_bytes_wide() { - let bignum = [255u8; 64]; // 2^512 - 1 - let reduced = Scalar64::from_bytes_wide(&bignum); - println!("{:?}", reduced); - for i in 0..5 { - assert!(reduced[i] == C[i]); - } - } -} - - -#[cfg(all(test, feature = "bench"))] -mod bench { - use test::Bencher; - - use super::*; - use super::test::{X, Y}; - - #[bench] - fn square(b: &mut Bencher) { - b.iter(|| X.square()); - } - - #[bench] - fn mul(b: &mut Bencher) { - b.iter(|| Scalar64::mul(&X, &Y)); - } - - #[bench] - fn montgomery_square(b: &mut Bencher) { - b.iter(|| X.montgomery_square()); - } - - #[bench] - fn montgomery_mul(b: &mut Bencher) { - b.iter(|| Scalar64::montgomery_mul(&X, &Y)); - } - - #[bench] - fn from_bytes_wide(b: &mut Bencher) { - let bignum = [255u8; 64]; // 2^512 - 1 - b.iter(|| Scalar64::from_bytes_wide(&bignum)); - } -}