saymrwulf/curve25519-dalek-source
1
0
Fork 0
mirror of https://github.com/saymrwulf/curve25519-dalek-source.git synced 2026-09-04 20:24:10 +00:00
Code Issues Projects Releases Packages Wiki Activity Actions

Revert "Revert "Merge remote-tracking branch 'floodyberry/optimized_scalar' into develop""

Browse source 
This reverts commit 90b69c13ee.

Signed-off-by: Isis Lovecruft <isis@torproject.org>
This commit is contained in:
Isis Lovecruft 2017-10-30 19:27:44 +00:00
parent 3f52aa4df2
commit 6079b0269f
Failed to extract signature
6 changed files with 1143 additions and 361 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

21
src/constants_32bit.rs
View file

@ -19,6 +19,7 @@
#![allow(non_snake_case)]
use field_32bit::FieldElement32;
use scalar_32bit::Scalar32;
use edwards::ExtendedPoint;
use edwards::AffineNielsPoint;
use edwards::EdwardsBasepointTable;
@ -96,6 +97,26 @@ 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{

13
src/constants_64bit.rs
View file

@ -19,6 +19,7 @@
#![allow(non_snake_case)]
use field_64bit::FieldElement64;
use scalar_64bit::Scalar64;
use edwards::ExtendedPoint;
use edwards::AffineNielsPoint;
use edwards::EdwardsBasepointTable;
@ -69,6 +70,18 @@ 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{

7
src/lib.rs
View file

@ -41,6 +41,8 @@ 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;
@ -71,6 +73,11 @@ 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;

429
src/scalar.rs
View file

@ -22,12 +22,11 @@
//!
//! The `Scalar` struct represents an element in ℤ/lℤ.
//!
//! Arithmetic operations on `Scalar`s are done using 12 21-bit limbs.
//! However, in contrast to `FieldElement`s, `Scalar`s are stored in
//! 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 is stored as
//! limbs.
//! between two scalars, the `UnpackedScalar` struct (internally
//! either `Scalar32` or `Scalar64`) is stored as limbs.
use core::fmt::Debug;
use core::ops::Neg;
@ -43,9 +42,6 @@ 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;
@ -108,50 +104,47 @@ impl IndexMut<usize> for Scalar {
impl<'b> MulAssign<&'b Scalar> for Scalar {
fn mul_assign(&mut self, _rhs: &'b Scalar) {
let result = (self as &Scalar) * _rhs;
self.0 = result.0;
*self = Scalar::mul(self, _rhs)
}
}
impl<'a, 'b> Mul<&'b Scalar> for &'a Scalar {
type Output = Scalar;
fn mul(self, _rhs: &'b Scalar) -> Scalar {
Scalar::multiply_add(self, _rhs, &Scalar::zero())
Scalar::mul(self, _rhs)
}
}
impl<'b> AddAssign<&'b Scalar> for Scalar {
fn add_assign(&mut self, _rhs: &'b Scalar) {
*self = Scalar::multiply_add(&Scalar::one(), self, _rhs);
*self = Scalar::add(self, _rhs);
}
}
impl<'a, 'b> Add<&'b Scalar> for &'a Scalar {
type Output = Scalar;
fn add(self, _rhs: &'b Scalar) -> Scalar {
Scalar::multiply_add(&Scalar::one(), self, _rhs)
Scalar::add(self, _rhs)
}
}
impl<'b> SubAssign<&'b Scalar> for Scalar {
fn sub_assign(&mut self, _rhs: &'b Scalar) {
// (l-1)*_rhs + self = self - _rhs
*self = Scalar::multiply_add(&constants::l_minus_1, _rhs, self);
*self = Scalar::sub(self, _rhs);
}
}
impl<'a, 'b> Sub<&'b Scalar> for &'a Scalar {
type Output = Scalar;
fn sub(self, _rhs: &'b Scalar) -> Scalar {
// (l-1)*_rhs + self = self - _rhs
Scalar::multiply_add(&constants::l_minus_1, _rhs, self)
Scalar::sub(self, _rhs)
}
}
impl<'a> Neg for &'a Scalar {
type Output = Scalar;
fn neg(self) -> Scalar {
self * &constants::l_minus_1
Scalar::sub(&Scalar::zero(), self)
}
}
@ -234,6 +227,18 @@ 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.
///
@ -384,26 +389,6 @@ 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
///
@ -442,127 +427,41 @@ impl Scalar {
output
}
/// 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 {
// Unpack scalars into limbs
let al = a.unpack();
let bl = b.unpack();
let cl = c.unpack();
/// Unpack this `Scalar` to an `UnpackedScalar`
pub fn unpack(&self) -> UnpackedScalar {
UnpackedScalar::from_bytes(&self.0)
}
// Multiply and repack
UnpackedScalar::multiply_add(&al, &bl, &cl).pack()
/// 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).
pub fn multiply_add(a: &Scalar, b: &Scalar, c: &Scalar) -> Scalar {
UnpackedScalar::add(&UnpackedScalar::mul(&a.unpack(), &b.unpack()), &c.unpack()).pack()
}
/// Reduce a 512-bit little endian number mod l
pub fn reduce(input: &[u8; 64]) -> Scalar {
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<usize> for UnpackedScalar {
type Output = i64;
fn index(&self, _index: usize) -> &i64 {
&(self.0[_index])
}
}
impl IndexMut<usize> for UnpackedScalar {
fn index_mut(&mut self, _index: usize) -> &mut i64 {
&mut (self.0[_index])
UnpackedScalar::from_bytes_wide(input).pack()
}
}
impl UnpackedScalar {
/// Pack the limbs of this `UnpackedScalar` into a `Scalar`.
fn pack(&self) -> Scalar {
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])
Scalar(self.to_bytes())
}
/// Compute the multiplicative inverse of this scalar.
@ -571,25 +470,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;
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());
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);
// _10000
let mut y = UnpackedScalar::multiply_add(&_1111, &_1, &UnpackedScalar::zero());
let mut y = UnpackedScalar::montgomery_mul(&_1111, &_1);
#[inline]
fn square_multiply(y: &mut UnpackedScalar, squarings: usize, x: &UnpackedScalar) {
for _ in 0..squarings {
*y = y.square();
*y = y.montgomery_square();
}
*y = UnpackedScalar::multiply_add(y, x, &UnpackedScalar::zero());
*y = UnpackedScalar::montgomery_mul(y, x);
}
square_multiply(&mut y, 123 + 3, &_101);
@ -620,194 +519,14 @@ impl UnpackedScalar {
square_multiply(&mut y, 3, &_101);
square_multiply(&mut y, 1 + 2, &_11);
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))
y.from_montgomery()
}
}
#[cfg(test)]
mod test {
use super::*;
use constants;
/// x = 2238329342913194256032495932344128051776374960164957527413114840482143558222
pub static X: Scalar = Scalar(
@ -815,6 +534,12 @@ 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,
@ -952,6 +677,7 @@ 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());
}
@ -984,7 +710,7 @@ mod bench {
use test::Bencher;
use super::*;
use super::test::{X, Y, Z};
use super::test::{X};
#[bench]
fn scalar_random(b: &mut Bencher) {
@ -993,28 +719,9 @@ 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));
}
}

559
src/scalar_32bit.rs Normal file
View file

@ -0,0 +1,559 @@
//! 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<usize> for Scalar32 {
type Output = u32;
fn index(&self, _index: usize) -> &u32 {
&(self.0[_index])
}
}
impl IndexMut<usize> 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));
}
}

475
src/scalar_64bit.rs Normal file
View file

@ -0,0 +1,475 @@
//! 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<usize> for Scalar64 {
type Output = u64;
fn index(&self, _index: usize) -> &u64 {
&(self.0[_index])
}
}
impl IndexMut<usize> 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));
}
}