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

Merge remote-tracking branch 'hdevalence/feature/64bit-multiplication_r3' into develop

Browse source
This commit is contained in:
Isis Lovecruft 2017-03-14 04:23:11 +00:00
commit 6991b4264b
Failed to extract signature
6 changed files with 2106 additions and 165 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

4
Cargo.toml
View file

@ -36,13 +36,15 @@ version = "^0.6"
version = "0.4"
[features]
nightly = ["basepoint_table_creation"]
nightly = ["basepoint_table_creation", "radix_51"]
default = ["std"]
std = ["rand"]
yolocrypto = []
# Needs nightly for placement new
basepoint_table_creation = []
bench = []
# Radix-51 arithmetic using u128
radix_51 = []
# The development profile, used for `cargo build`.
[profile.dev]

1576
src/constants.rs
View file

File diff suppressed because it is too large Load diff

35
src/curve.rs
View file

@ -1192,10 +1192,9 @@ mod test {
/// XXX what does Signal do here?
#[test]
fn u_minus_one_monty() {
let mut m1 = FieldElement::zero();
m1[0] = -1;
let m1_bytes = m1.to_bytes();
let div_by_zero_u = CompressedMontgomeryU(m1_bytes);
let minus_one = FieldElement::minus_one();
let minus_one_bytes = minus_one.to_bytes();
let div_by_zero_u = CompressedMontgomeryU(minus_one_bytes);
assert!(div_by_zero_u.decompress().is_none());
}
@ -1400,6 +1399,26 @@ mod test {
assert!( ExtendedPoint::identity().is_identity() == true);
assert!(constants::ED25519_BASEPOINT.is_identity() == false);
}
/// Rust's debug builds have overflow and underflow trapping,
/// and enable `debug_assert!()`. This performs many scalar
/// multiplications to attempt to trigger possible overflows etc.
///
/// For instance, the `radix_51` `Mul` implementation for
/// `FieldElements` requires the input `Limb`s to be bounded by
/// 2^54, but we cannot enforce this dynamically at runtime, or
/// statically at compile time (until Rust gets type-level
/// integers, at which point we can encode "bits of headroom" into
/// the type system and prove correctness).
#[test]
fn monte_carlo_overflow_underflow_debug_assert_test() {
let mut P = ExtendedPoint::basepoint();
// N.B. each scalar_mult does 1407 field mults, 1024 field squarings,
// so this does ~ 1M of each operation.
for _ in 0..1_000 {
P = P.scalar_mult(&A_SCALAR);
}
}
}
// ------------------------------------------------------------------------
@ -1436,7 +1455,7 @@ mod bench {
}
#[bench]
fn add_extended_and_cached_output_completed(b: &mut Bencher) {
fn add_extended_and_projective_niels_output_completed(b: &mut Bencher) {
let p1 = constants::ED25519_BASEPOINT;
let p2 = constants::ED25519_BASEPOINT.to_projective_niels();
@ -1444,7 +1463,7 @@ mod bench {
}
#[bench]
fn add_extended_and_cached_output_extended(b: &mut Bencher) {
fn add_extended_and_projective_niels_output_extended(b: &mut Bencher) {
let p1 = constants::ED25519_BASEPOINT;
let p2 = constants::ED25519_BASEPOINT.to_projective_niels();
@ -1452,7 +1471,7 @@ mod bench {
}
#[bench]
fn add_extended_and_precomputed_output_completed(b: &mut Bencher) {
fn add_extended_and_affine_niels_output_completed(b: &mut Bencher) {
let p1 = constants::ED25519_BASEPOINT;
let p2 = constants::ED25519_BASEPOINT.to_affine_niels();
@ -1460,7 +1479,7 @@ mod bench {
}
#[bench]
fn add_extended_and_precomputed_output_extended(b: &mut Bencher) {
fn add_extended_and_affine_niels_output_extended(b: &mut Bencher) {
let p1 = constants::ED25519_BASEPOINT;
let p2 = constants::ED25519_BASEPOINT.to_affine_niels();

641
src/field.rs
View file

@ -29,22 +29,37 @@ use subtle::byte_is_nonzero;
use subtle::CTAssignable;
use subtle::CTEq;
use utils::{load3, load4};
use utils::{load3, load4, load8};
use constants;
/// With the `radix51` feature enabled, `FieldElements` are represented
/// in radix 2^51 as five `u64`s.
#[cfg(feature="radix_51")]
pub type Limb = u64;
/// FieldElement represents an element of the field GF(2^255 - 19). An element
/// t, entries t[0]...t[9], represents the integer t[0]+2^26 t[1]+2^51 t[2]+2^77
/// t[3]+2^102 t[4]+...+2^230 t[9]. Bounds on each t[i] vary depending on
/// context.
#[cfg(feature="radix_51")]
#[derive(Copy, Clone)]
pub struct FieldElement(pub [u64; 5]);
/// FieldElements are represented as an array of ten "Limbs", which are radix
/// 25.5, that is, each Limb of a FieldElement alternates between being
/// represented as a factor of 2^25 or 2^26 more than the last corresponding
/// integer.
#[cfg(not(feature="radix_51"))]
pub type Limb = i32;
/// FieldElement represents an element of the field GF(2^255 - 19). An element
/// t, entries t[0]...t[9], represents the integer t[0]+2^26 t[1]+2^51 t[2]+2^77
/// t[3]+2^102 t[4]+...+2^230 t[9]. Bounds on each t[i] vary depending on
/// context.
#[cfg(not(feature="radix_51"))]
#[derive(Copy, Clone)]
pub struct FieldElement(pub [Limb; 10]);
pub struct FieldElement(pub [i32; 10]);
impl Eq for FieldElement {}
impl PartialEq for FieldElement {
@ -105,7 +120,7 @@ impl IndexMut<usize> for FieldElement {
impl<'b> AddAssign<&'b FieldElement> for FieldElement {
fn add_assign(&mut self, _rhs: &'b FieldElement) { // fsum()
for i in 0..10 {
for i in 0..self.0.len() {
self[i] += _rhs[i];
}
}
@ -121,32 +136,187 @@ impl<'a, 'b> Add<&'b FieldElement> for &'a FieldElement {
}
impl<'b> SubAssign<&'b FieldElement> for FieldElement {
#[cfg(not(feature="radix_51"))]
fn sub_assign(&mut self, _rhs: &'b FieldElement) { // fdifference()
for i in 0..10 {
self[i] -= _rhs[i];
}
}
#[cfg(feature="radix_51")]
fn sub_assign(&mut self, _rhs: &'b FieldElement) {
let result = (self as &FieldElement) - _rhs;
for i in 0..5 {
self.0[i] = result.0[i];
}
}
}
impl<'a, 'b> Sub<&'b FieldElement> for &'a FieldElement {
type Output = FieldElement;
#[cfg(not(feature="radix_51"))]
fn sub(self, _rhs: &'b FieldElement) -> FieldElement {
let mut output = self.clone();
output -= _rhs;
output
}
#[cfg(feature="radix_51")]
fn sub(self, _rhs: &'b FieldElement) -> FieldElement {
// To avoid underflow, first add a multiple of p.
// Choose 16*p = p << 4 to be larger than 54-bit _rhs.
//
// If we could statically track the bitlengths of the limbs
// of every FieldElement, we could choose a multiple of p
// just bigger than _rhs and avoid having to do a reduction.
//
// Since we don't yet have type-level integers to do this, we
// have to add an explicit reduction call here, which is a
// significant cost.
FieldElement::reduce([
(self.0[0] + 36028797018963664u64) - _rhs.0[0],
(self.0[1] + 36028797018963952u64) - _rhs.0[1],
(self.0[2] + 36028797018963952u64) - _rhs.0[2],
(self.0[3] + 36028797018963952u64) - _rhs.0[3],
(self.0[4] + 36028797018963952u64) - _rhs.0[4],
])
}
}
impl<'b> MulAssign<&'b FieldElement> for FieldElement {
fn mul_assign(&mut self, _rhs: &'b FieldElement) {
self.0 = self.multiply(_rhs).0;
let result = (self as &FieldElement) * _rhs;
self.0 = result.0;
}
}
impl<'a, 'b> Mul<&'b FieldElement> for &'a FieldElement {
type Output = FieldElement;
#[cfg(feature="radix_51")]
fn mul(self, _rhs: &'b FieldElement) -> FieldElement {
self.multiply(_rhs)
/// Multiply two 64-bit integers with 128 bits of output.
#[inline(always)]
fn m(x: u64, y: u64) -> u128 { (x as u128) * (y as u128) }
// Alias self, _rhs for more readable formulas
let a: &[u64; 5] = &self.0;
let b: &[u64; 5] = &_rhs.0;
// Multiply to get 128-bit coefficients of output
let c0: u128 = m(a[0],b[0]) + ( m(a[4],b[1]) + m(a[3],b[2]) + m(a[2],b[3]) + m(a[1],b[4]) )*19;
let mut c1: u128 = m(a[1],b[0]) + m(a[0],b[1]) + ( m(a[4],b[2]) + m(a[3],b[3]) + m(a[2],b[4]) )*19;
let mut c2: u128 = m(a[2],b[0]) + m(a[1],b[1]) + m(a[0],b[2]) + ( m(a[4],b[3]) + m(a[3],b[4]) )*19;
let mut c3: u128 = m(a[3],b[0]) + m(a[2],b[1]) + m(a[1],b[2]) + m(a[0],b[3]) + ( m(a[4],b[4]) )*19;
let mut c4: u128 = m(a[4],b[0]) + m(a[3],b[1]) + m(a[2],b[2]) + m(a[1],b[3]) + m(a[0],b[4]);
// Now c[i] < 2^2b * (1+i + (4-i)*19) < 2^(2b + lg(1+4*19)) < 2^(2b + 6.27)
// where b is the bitlength of the input limbs.
// The carry (c[i] >> 51) fits into a u64 iff 2b+6.27 < 64+51 iff b <= 54.
// After the first carry pass, all c[i] fit into u64.
debug_assert!(a[0] < (1 << 54)); debug_assert!(b[0] < (1 << 54));
debug_assert!(a[1] < (1 << 54)); debug_assert!(b[1] < (1 << 54));
debug_assert!(a[2] < (1 << 54)); debug_assert!(b[2] < (1 << 54));
debug_assert!(a[3] < (1 << 54)); debug_assert!(b[3] < (1 << 54));
debug_assert!(a[4] < (1 << 54)); debug_assert!(b[4] < (1 << 54));
// The 128-bit output limbs are stored in two 64-bit registers (low/high part).
// By rebinding the names after carrying, we free the upper registers for reuse.
let low_51_bit_mask = (1u64 << 51) - 1;
c1 += (c0 >> 51) as u128;
let mut c0: u64 = (c0 as u64) & low_51_bit_mask;
c2 += (c1 >> 51) as u128;
let c1: u64 = (c1 as u64) & low_51_bit_mask;
c3 += (c2 >> 51) as u128;
let c2: u64 = (c2 as u64) & low_51_bit_mask;
c4 += (c3 >> 51) as u128;
let c3: u64 = (c3 as u64) & low_51_bit_mask;
c0 += ((c4 >> 51) as u64) * 19;
let c4: u64 = (c4 as u64) & low_51_bit_mask;
FieldElement::reduce([c0,c1,c2,c3,c4])
}
#[cfg(not(feature="radix_51"))]
fn mul(self, _rhs: &'b FieldElement) -> FieldElement {
// Notes preserved from ed25519.go (presumably originally from ref10):
//
// Calculates h = f * g. Can overlap h with f or g.
//
// # Preconditions
//
// * |f[i]| bounded by 1.1*2^26, 1.1*2^25, 1.1*2^26, 1.1*2^25, etc.
// * |g[i]| bounded by 1.1*2^26, 1.1*2^25, 1.1*2^26, 1.1*2^25, etc.
//
// # Postconditions
//
// * |h| bounded by 1.1*2^25, 1.1*2^24, 1.1*2^25, 1.1*2^24, etc.
//
// ## Notes on implementation strategy
//
// * Using schoolbook multiplication.
// * Karatsuba would save a little in some cost models.
//
// * Most multiplications by 2 and 19 are 32-bit precomputations;
// cheaper than 64-bit postcomputations.
//
// * There is one remaining multiplication by 19 in the carry chain;
// one *19 precomputation can be merged into this,
// but the resulting data flow is considerably less clean.
//
// * There are 12 carries below.
// 10 of them are 2-way parallelizable and vectorizable.
// Can get away with 11 carries, but then data flow is much deeper.
//
// * With tighter constraints on inputs can squeeze carries into int32.
let f0 = self[0] as i64;
let f1 = self[1] as i64;
let f2 = self[2] as i64;
let f3 = self[3] as i64;
let f4 = self[4] as i64;
let f5 = self[5] as i64;
let f6 = self[6] as i64;
let f7 = self[7] as i64;
let f8 = self[8] as i64;
let f9 = self[9] as i64;
let f1_2 = (2 * self[1]) as i64;
let f3_2 = (2 * self[3]) as i64;
let f5_2 = (2 * self[5]) as i64;
let f7_2 = (2 * self[7]) as i64;
let f9_2 = (2 * self[9]) as i64;
let g0 = _rhs[0] as i64;
let g1 = _rhs[1] as i64;
let g2 = _rhs[2] as i64;
let g3 = _rhs[3] as i64;
let g4 = _rhs[4] as i64;
let g5 = _rhs[5] as i64;
let g6 = _rhs[6] as i64;
let g7 = _rhs[7] as i64;
let g8 = _rhs[8] as i64;
let g9 = _rhs[9] as i64;
let g1_19 = (19 * _rhs[1]) as i64; /* 1.4*2^29 */
let g2_19 = (19 * _rhs[2]) as i64; /* 1.4*2^30; still ok */
let g3_19 = (19 * _rhs[3]) as i64;
let g4_19 = (19 * _rhs[4]) as i64;
let g5_19 = (19 * _rhs[5]) as i64;
let g6_19 = (19 * _rhs[6]) as i64;
let g7_19 = (19 * _rhs[7]) as i64;
let g8_19 = (19 * _rhs[8]) as i64;
let g9_19 = (19 * _rhs[9]) as i64;
let h0 = f0*g0 + f1_2*g9_19 + f2*g8_19 + f3_2*g7_19 + f4*g6_19 + f5_2*g5_19 + f6*g4_19 + f7_2*g3_19 + f8*g2_19 + f9_2*g1_19;
let h1 = f0*g1 + f1*g0 + f2*g9_19 + f3*g8_19 + f4*g7_19 + f5*g6_19 + f6*g5_19 + f7*g4_19 + f8*g3_19 + f9*g2_19;
let h2 = f0*g2 + f1_2*g1 + f2*g0 + f3_2*g9_19 + f4*g8_19 + f5_2*g7_19 + f6*g6_19 + f7_2*g5_19 + f8*g4_19 + f9_2*g3_19;
let h3 = f0*g3 + f1*g2 + f2*g1 + f3*g0 + f4*g9_19 + f5*g8_19 + f6*g7_19 + f7*g6_19 + f8*g5_19 + f9*g4_19;
let h4 = f0*g4 + f1_2*g3 + f2*g2 + f3_2*g1 + f4*g0 + f5_2*g9_19 + f6*g8_19 + f7_2*g7_19 + f8*g6_19 + f9_2*g5_19;
let h5 = f0*g5 + f1*g4 + f2*g3 + f3*g2 + f4*g1 + f5*g0 + f6*g9_19 + f7*g8_19 + f8*g7_19 + f9*g6_19;
let h6 = f0*g6 + f1_2*g5 + f2*g4 + f3_2*g3 + f4*g2 + f5_2*g1 + f6*g0 + f7_2*g9_19 + f8*g8_19 + f9_2*g7_19;
let h7 = f0*g7 + f1*g6 + f2*g5 + f3*g4 + f4*g3 + f5*g2 + f6*g1 + f7*g0 + f8*g9_19 + f9*g8_19;
let h8 = f0*g8 + f1_2*g7 + f2*g6 + f3_2*g5 + f4*g4 + f5_2*g3 + f6*g2 + f7_2*g1 + f8*g0 + f9_2*g9_19;
let h9 = f0*g9 + f1*g8 + f2*g7 + f3*g6 + f4*g5 + f5*g4 + f6*g3 + f7*g2 + f8*g1 + f9*g0;
FieldElement::reduce(&[h0, h1, h2, h3, h4, h5, h6, h7, h8, h9])
}
}
@ -162,6 +332,9 @@ impl<'a> Neg for &'a FieldElement {
impl CTAssignable for FieldElement {
/// Conditionally assign another FieldElement to this one.
///
/// XXX fixup tests to avoid limb specs
/// XXX_radix_51
///
/// If `choice == 0`, replace `self` with `self`:
///
/// ```
@ -189,38 +362,97 @@ impl CTAssignable for FieldElement {
/// # Preconditions
///
/// * `choice` in {0,1}
#[cfg(not(feature="radix_51"))]
fn conditional_assign(&mut self, f: &FieldElement, choice: u8) {
let mask = -(choice as Limb);
for i in 0..10 {
self[i] ^= mask & (self[i] ^ f[i]);
}
}
#[cfg(feature="radix_51")]
fn conditional_assign(&mut self, f: &FieldElement, choice: u8) {
let mask = (-(choice as i64)) as u64;
for i in 0..5 {
self.0[i] ^= mask & (self.0[i] ^ f.0[i]);
}
}
}
impl FieldElement {
/// Invert the sign of this field element
#[cfg(not(feature="radix_51"))]
pub fn negate(&mut self) {
for i in 0..10 {
self[i] = -self[i];
}
}
/// Invert the sign of this field element
#[cfg(feature="radix_51")]
pub fn negate(&mut self) {
// See commentary in the Sub impl
let neg = FieldElement::reduce([
36028797018963664u64 - self.0[0],
36028797018963952u64 - self.0[1],
36028797018963952u64 - self.0[2],
36028797018963952u64 - self.0[3],
36028797018963952u64 - self.0[4],
]);
self.0 = neg.0;
}
/// Construct the additive identity
/// Construct zero.
#[cfg(not(feature="radix_51"))]
pub fn zero() -> FieldElement {
FieldElement([ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 ])
}
/// Construct zero.
#[cfg(feature="radix_51")]
pub fn zero() -> FieldElement {
FieldElement([ 0, 0, 0, 0, 0 ])
}
/// Construct the multiplicative identity
/// Construct one.
#[cfg(not(feature="radix_51"))]
pub fn one() -> FieldElement {
FieldElement([ 1, 0, 0, 0, 0, 0, 0, 0, 0, 0 ])
}
/// Construct one.
#[cfg(feature="radix_51")]
pub fn one() -> FieldElement {
FieldElement([ 1, 0, 0, 0, 0 ])
}
/// Construct -1.
#[cfg(not(feature="radix_51"))]
pub fn minus_one() -> FieldElement {
FieldElement([-1, 0, 0, 0, 0, 0, 0, 0, 0, 0 ])
}
/// Construct -1.
#[cfg(feature="radix_51")]
pub fn minus_one() -> FieldElement {
FieldElement([2251799813685228, 2251799813685247, 2251799813685247, 2251799813685247, 2251799813685247])
}
fn combine_coeffs(input: &[i64;10]) -> FieldElement { //FeCombine
/// Given 64-bit limbs, reduce to enforce the bound c_i < 2^51.
#[cfg(feature="radix_51")]
#[inline(always)]
fn reduce(mut limbs: [u64; 5]) -> FieldElement {
let low_51_bit_mask = (1u64 << 51) - 1;
limbs[1] += limbs[0] >> 51;
limbs[0] = limbs[0] & low_51_bit_mask;
limbs[2] += limbs[1] >> 51;
limbs[1] = limbs[1] & low_51_bit_mask;
limbs[3] += limbs[2] >> 51;
limbs[2] = limbs[2] & low_51_bit_mask;
limbs[4] += limbs[3] >> 51;
limbs[3] = limbs[3] & low_51_bit_mask;
limbs[0] += (limbs[4] >> 51) * 19;
limbs[4] = limbs[4] & low_51_bit_mask;
FieldElement(limbs)
}
#[cfg(not(feature="radix_51"))]
fn reduce(input: &[i64;10]) -> FieldElement { //FeCombine
let mut c = [0i64;10];
let mut h = input.clone();
@ -315,6 +547,7 @@ impl FieldElement {
/// Create a FieldElement by demarshalling an array of 32 bytes.
///
/// # Example
/// XXX eliminate limbs
///
/// ```
/// # use curve25519_dalek::field::FieldElement;
@ -331,6 +564,7 @@ impl FieldElement {
/// # Return
///
/// Returns a new FieldElement.
#[cfg(not(feature="radix_51"))]
pub fn from_bytes(data: &[u8; 32]) -> FieldElement { //FeFromBytes
let mut h = [0i64;10];
h[0] = load4(&data[ 0..]);
@ -344,41 +578,29 @@ impl FieldElement {
h[8] = load3(&data[26..]) << 4;
h[9] = (load3(&data[29..]) & 8388607) << 2;
FieldElement::combine_coeffs(&h)
FieldElement::reduce(&h)
}
/// Parse a `FieldElement` from 32 bytes.
#[cfg(feature="radix_51")]
pub fn from_bytes(bytes: &[u8;32]) -> FieldElement {
let low_51_bit_mask = (1u64 << 51) - 1;
FieldElement(
// load bits [ 0, 64), no shift
[ load8(&bytes[ 0..]) & low_51_bit_mask
// load bits [ 48,112), shift to [ 51,112)
, (load8(&bytes[ 6..]) >> 3) & low_51_bit_mask
// load bits [ 96,160), shift to [102,160)
, (load8(&bytes[12..]) >> 6) & low_51_bit_mask
// load bits [152,216), shift to [153,216)
, (load8(&bytes[19..]) >> 1) & low_51_bit_mask
// load bits [192,256), shift to [204,112)
, (load8(&bytes[24..]) >> 12) & low_51_bit_mask
])
}
/// Marshal this FieldElement into a 32-byte array.
///
/// # Preconditions
///
/// * `|h[i]|` bounded by 1.1*2^25, 1.1*2^24, 1.1*2^25, 1.1*2^24, etc.
///
/// # Lemma
///
/// Write p = 2^255 - 19 and q = floor(h/p).
///
/// Basic claim: q = floor(2^(-255)(h + 19 * 2^-25 h9 + 2^-1)).
///
/// # Proof
///
/// Have |h|<=p so |q|<=1 so |19^2 * 2^-255 * q| < 1/4.
///
/// Also have |h-2^230 * h9| < 2^230 so |19 * 2^-255 * (h-2^230 * h9)| < 1/4.
///
/// Write y=2^(-1)-19^2 2^(-255)q-19 2^(-255)(h-2^230 h9), then 0<y<1.
///
/// Write r = h - pq.
///
/// Have 0 <= r< = p-1 = 2^255 - 20.
///
/// Thus 0 <= r + 19 * 2^-255 * r < r + 19 * 2^-255 * 2^255 <= 2^255 - 1.
///
/// Write x = r + 19 * 2^-255 * r + y.
///
/// Then 0 < x < 2^255 so floor(2^(-255)x) = 0 so floor(q+2^(-255)x) = q.
///
/// Have q+2^(-255)x = 2^-255 * (h + 19 * 2^-25 * h9 + 2^-1),
/// so floor(2^-255 * (h + 19 * 2^-25 * h9 + 2^-1)) = q.
/// XXX eliminate limbs
///
/// # Example
///
@ -395,7 +617,41 @@ impl FieldElement {
/// let bytes: [u8; 32] = fe.to_bytes();
/// assert!(data == bytes);
/// ```
#[cfg(not(feature="radix_51"))]
pub fn to_bytes(&self) -> [u8;32] { //FeToBytes
// Comment preserved from ed25519.go (presumably originally from ref10):
//
// # Preconditions
//
// * `|h[i]|` bounded by 1.1*2^25, 1.1*2^24, 1.1*2^25, 1.1*2^24, etc.
//
// # Lemma
//
// Write p = 2^255 - 19 and q = floor(h/p).
//
// Basic claim: q = floor(2^(-255)(h + 19 * 2^-25 h9 + 2^-1)).
//
// # Proof
//
// Have |h|<=p so |q|<=1 so |19^2 * 2^-255 * q| < 1/4.
//
// Also have |h-2^230 * h9| < 2^230 so |19 * 2^-255 * (h-2^230 * h9)| < 1/4.
//
// Write y=2^(-1)-19^2 2^(-255)q-19 2^(-255)(h-2^230 h9), then 0<y<1.
//
// Write r = h - pq.
//
// Have 0 <= r< = p-1 = 2^255 - 20.
//
// Thus 0 <= r + 19 * 2^-255 * r < r + 19 * 2^-255 * 2^255 <= 2^255 - 1.
//
// Write x = r + 19 * 2^-255 * r + y.
//
// Then 0 < x < 2^255 so floor(2^(-255)x) = 0 so floor(q+2^(-255)x) = q.
//
// Have q+2^(-255)x = 2^-255 * (h + 19 * 2^-25 * h9 + 2^-1),
// so floor(2^-255 * (h + 19 * 2^-25 * h9 + 2^-1)) = q.
//
let mut carry = [0i32; 10];
let mut h = self.clone();
@ -486,10 +742,90 @@ impl FieldElement {
s[31] = (h[9] >> 18) as u8;
//Clear high bit
debug_assert!((s[31] & 0b1000_0000u8) == 0u8);
s[31] &= 127u8;
s
}
/// Serialize this `FieldElement` to bytes.
#[cfg(feature="radix_51")]
pub fn to_bytes(&self) -> [u8;32] {
// This reduces to the range [0,2^255), but we need [0,2^255-19)
let mut limbs = FieldElement::reduce(self.0).0;
// Let h = limbs[0] + limbs[1]*2^51 + ... + limbs[4]*2^204.
//
// Write h = pq + r with 0 <= r < p. We want to compute r = h mod p.
//
// Since h < 2^255, q = 0 or 1, with q = 0 when h < p and q = 1 when h >= p.
//
// Notice that h >= p <==> h + 19 >= p + 19 <==> h + 19 >= 2^255.
// Therefore q can be computed as the carry bit of h + 19.
let mut q = (limbs[0] + 19) >> 51;
q = (limbs[1] + q) >> 51;
q = (limbs[2] + q) >> 51;
q = (limbs[3] + q) >> 51;
q = (limbs[4] + q) >> 51;
// Now we can compute r as r = h - pq = r - (2^255-19)q = r + 19q - 2^255q
limbs[0] += 19*q;
// Now carry the result to compute r + 19q ...
let low_51_bit_mask = (1u64 << 51) - 1;
limbs[1] += limbs[0] >> 51;
limbs[0] = limbs[0] & low_51_bit_mask;
limbs[2] += limbs[1] >> 51;
limbs[1] = limbs[1] & low_51_bit_mask;
limbs[3] += limbs[2] >> 51;
limbs[2] = limbs[2] & low_51_bit_mask;
limbs[4] += limbs[3] >> 51;
limbs[3] = limbs[3] & low_51_bit_mask;
// ... but instead of carrying (limbs[4] >> 51) = 2^255q
// into another limb, discard it, subtracting the value
limbs[4] = limbs[4] & low_51_bit_mask;
// Now arrange the bits of the limbs.
let mut s = [0u8;32];
s[ 0] = limbs[0] as u8;
s[ 1] = (limbs[0] >> 8) as u8;
s[ 2] = (limbs[0] >> 16) as u8;
s[ 3] = (limbs[0] >> 24) as u8;
s[ 4] = (limbs[0] >> 32) as u8;
s[ 5] = (limbs[0] >> 40) as u8;
s[ 6] = ((limbs[0] >> 48) | (limbs[1] << 3)) as u8;
s[ 7] = (limbs[1] >> 5) as u8;
s[ 8] = (limbs[1] >> 13) as u8;
s[ 9] = (limbs[1] >> 21) as u8;
s[10] = (limbs[1] >> 29) as u8;
s[11] = (limbs[1] >> 37) as u8;
s[12] = ((limbs[1] >> 45) | (limbs[2] << 6)) as u8;
s[13] = (limbs[2] >> 2) as u8;
s[14] = (limbs[2] >> 10) as u8;
s[15] = (limbs[2] >> 18) as u8;
s[16] = (limbs[2] >> 26) as u8;
s[17] = (limbs[2] >> 34) as u8;
s[18] = (limbs[2] >> 42) as u8;
s[19] = ((limbs[2] >> 50) | (limbs[3] << 1)) as u8;
s[20] = (limbs[3] >> 7) as u8;
s[21] = (limbs[3] >> 15) as u8;
s[22] = (limbs[3] >> 23) as u8;
s[23] = (limbs[3] >> 31) as u8;
s[24] = (limbs[3] >> 39) as u8;
s[25] = ((limbs[3] >> 47) | (limbs[4] << 4)) as u8;
s[26] = (limbs[4] >> 4) as u8;
s[27] = (limbs[4] >> 12) as u8;
s[28] = (limbs[4] >> 20) as u8;
s[29] = (limbs[4] >> 28) as u8;
s[30] = (limbs[4] >> 36) as u8;
s[31] = (limbs[4] >> 44) as u8;
//Clear high bit
debug_assert!((s[31] & 0b1000_0000u8) == 0u8);
s[31] &= 127u8;
return s
}
/// XXX clarify documentation
/// Determine if this field element, represented as a byte array,
@ -585,87 +921,7 @@ impl FieldElement {
return byte_is_nonzero(x);
}
/// Calculates h = f * g. Can overlap h with f or g.
///
/// # Preconditions
///
/// * |f[i]| bounded by 1.1*2^26, 1.1*2^25, 1.1*2^26, 1.1*2^25, etc.
/// * |g[i]| bounded by 1.1*2^26, 1.1*2^25, 1.1*2^26, 1.1*2^25, etc.
///
/// # Postconditions
///
/// * |h| bounded by 1.1*2^25, 1.1*2^24, 1.1*2^25, 1.1*2^24, etc.
///
/// ## Notes on implementation strategy
///
/// * Using schoolbook multiplication.
/// * Karatsuba would save a little in some cost models.
///
/// * Most multiplications by 2 and 19 are 32-bit precomputations;
/// cheaper than 64-bit postcomputations.
///
/// * There is one remaining multiplication by 19 in the carry chain;
/// one *19 precomputation can be merged into this,
/// but the resulting data flow is considerably less clean.
///
/// * There are 12 carries below.
/// 10 of them are 2-way parallelizable and vectorizable.
/// Can get away with 11 carries, but then data flow is much deeper.
///
/// * With tighter constraints on inputs can squeeze carries into int32.
pub fn multiply(&self, _rhs: &FieldElement) -> FieldElement {
let f0 = self[0] as i64;
let f1 = self[1] as i64;
let f2 = self[2] as i64;
let f3 = self[3] as i64;
let f4 = self[4] as i64;
let f5 = self[5] as i64;
let f6 = self[6] as i64;
let f7 = self[7] as i64;
let f8 = self[8] as i64;
let f9 = self[9] as i64;
let f1_2 = (2 * self[1]) as i64;
let f3_2 = (2 * self[3]) as i64;
let f5_2 = (2 * self[5]) as i64;
let f7_2 = (2 * self[7]) as i64;
let f9_2 = (2 * self[9]) as i64;
let g0 = _rhs[0] as i64;
let g1 = _rhs[1] as i64;
let g2 = _rhs[2] as i64;
let g3 = _rhs[3] as i64;
let g4 = _rhs[4] as i64;
let g5 = _rhs[5] as i64;
let g6 = _rhs[6] as i64;
let g7 = _rhs[7] as i64;
let g8 = _rhs[8] as i64;
let g9 = _rhs[9] as i64;
let g1_19 = (19 * _rhs[1]) as i64; /* 1.4*2^29 */
let g2_19 = (19 * _rhs[2]) as i64; /* 1.4*2^30; still ok */
let g3_19 = (19 * _rhs[3]) as i64;
let g4_19 = (19 * _rhs[4]) as i64;
let g5_19 = (19 * _rhs[5]) as i64;
let g6_19 = (19 * _rhs[6]) as i64;
let g7_19 = (19 * _rhs[7]) as i64;
let g8_19 = (19 * _rhs[8]) as i64;
let g9_19 = (19 * _rhs[9]) as i64;
let h0 = f0*g0 + f1_2*g9_19 + f2*g8_19 + f3_2*g7_19 + f4*g6_19 + f5_2*g5_19 + f6*g4_19 + f7_2*g3_19 + f8*g2_19 + f9_2*g1_19;
let h1 = f0*g1 + f1*g0 + f2*g9_19 + f3*g8_19 + f4*g7_19 + f5*g6_19 + f6*g5_19 + f7*g4_19 + f8*g3_19 + f9*g2_19;
let h2 = f0*g2 + f1_2*g1 + f2*g0 + f3_2*g9_19 + f4*g8_19 + f5_2*g7_19 + f6*g6_19 + f7_2*g5_19 + f8*g4_19 + f9_2*g3_19;
let h3 = f0*g3 + f1*g2 + f2*g1 + f3*g0 + f4*g9_19 + f5*g8_19 + f6*g7_19 + f7*g6_19 + f8*g5_19 + f9*g4_19;
let h4 = f0*g4 + f1_2*g3 + f2*g2 + f3_2*g1 + f4*g0 + f5_2*g9_19 + f6*g8_19 + f7_2*g7_19 + f8*g6_19 + f9_2*g5_19;
let h5 = f0*g5 + f1*g4 + f2*g3 + f3*g2 + f4*g1 + f5*g0 + f6*g9_19 + f7*g8_19 + f8*g7_19 + f9*g6_19;
let h6 = f0*g6 + f1_2*g5 + f2*g4 + f3_2*g3 + f4*g2 + f5_2*g1 + f6*g0 + f7_2*g9_19 + f8*g8_19 + f9_2*g7_19;
let h7 = f0*g7 + f1*g6 + f2*g5 + f3*g4 + f4*g3 + f5*g2 + f6*g1 + f7*g0 + f8*g9_19 + f9*g8_19;
let h8 = f0*g8 + f1_2*g7 + f2*g6 + f3_2*g5 + f4*g4 + f5_2*g3 + f6*g2 + f7_2*g1 + f8*g0 + f9_2*g9_19;
let h9 = f0*g9 + f1*g8 + f2*g7 + f3*g6 + f4*g5 + f5*g4 + f6*g3 + f7*g2 + f8*g1 + f9*g0;
FieldElement::combine_coeffs(&[h0, h1, h2, h3, h4, h5, h6, h7, h8, h9])
}
#[cfg(not(feature="radix_51"))]
fn square_inner(&self) -> [i64;10] {
let f0 = self[0] as i64;
let f1 = self[1] as i64;
@ -705,9 +961,61 @@ impl FieldElement {
h
}
#[cfg(feature="radix_51")]
#[inline(always)]
fn square_inner(&self) -> [u64; 5] {
/// Multiply two 64-bit integers with 128 bits of output.
#[inline(always)]
fn m(x: u64, y: u64) -> u128 { (x as u128) * (y as u128) }
// Alias self, _rhs for more readable formulas
let a: &[u64; 5] = &self.0;
// Precomputation: 64-bit multiply by 19
let a3_19 = 19 * a[3];
let a4_19 = 19 * a[4];
// Multiply to get 128-bit coefficients of output
let c0: u128 = m(a[0], a[0]) + 2*( m(a[1], a4_19) + m(a[2], a3_19) );
let mut c1: u128 = m(a[3], a3_19) + 2*( m(a[0], a[1]) + m(a[2], a4_19) );
let mut c2: u128 = m(a[1], a[1]) + 2*( m(a[0], a[2]) + m(a[4], a3_19) );
let mut c3: u128 = m(a[4], a4_19) + 2*( m(a[0], a[3]) + m(a[1], a[2]) );
let mut c4: u128 = m(a[2], a[2]) + 2*( m(a[0], a[4]) + m(a[1], a[3]) );
// Same bound as in multiply:
// c[i] < 2^2b * (1+i + (4-i)*19) < 2^(2b + lg(1+4*19)) < 2^(2b + 6.27)
// where b is the bitlength of the input limbs.
//
// The carry (c[i] >> 51) fits into a u64 iff 2b+6.27 < 64+51 iff b <= 54.
// After the first carry pass, all c[i] fit into u64.
debug_assert!(a[0] < (1 << 54));
debug_assert!(a[1] < (1 << 54));
debug_assert!(a[2] < (1 << 54));
debug_assert!(a[3] < (1 << 54));
debug_assert!(a[4] < (1 << 54));
// The 128-bit output limbs are stored in two 64-bit registers (low/high part).
// By rebinding the names after carrying, we free the upper registers for reuse.
let low_51_bit_mask = (1u64 << 51) - 1;
c1 += (c0 >> 51) as u128;
let mut c0: u64 = (c0 as u64) & low_51_bit_mask;
c2 += (c1 >> 51) as u128;
let c1: u64 = (c1 as u64) & low_51_bit_mask;
c3 += (c2 >> 51) as u128;
let c2: u64 = (c2 as u64) & low_51_bit_mask;
c4 += (c3 >> 51) as u128;
let c3: u64 = (c3 as u64) & low_51_bit_mask;
c0 += ((c4 >> 51) as u64) * 19;
let c4: u64 = (c4 as u64) & low_51_bit_mask;
// Now c_i all fit into u64, but are not yet bounded by 2^51.
[c0,c1,c2,c3,c4]
}
/// Calculates h = f*f. Can overlap h with f.
///
/// XXX limbs: better to talk about headroom?
///
/// # Preconditions
///
/// * |f[i]| bounded by 1.1*2^26, 1.1*2^25, 1.1*2^26, 1.1*2^25, etc.
@ -715,12 +1023,20 @@ impl FieldElement {
/// # Postconditions
///
/// * |h[i]| bounded by 1.1*2^25, 1.1*2^24, 1.1*2^25, 1.1*2^24, etc.
#[cfg(not(feature="radix_51"))]
pub fn square(&self) -> FieldElement {
FieldElement::combine_coeffs(&self.square_inner())
FieldElement::reduce(&self.square_inner())
}
/// Compute `self^2`.
#[cfg(feature="radix_51")]
pub fn square(&self) -> FieldElement {
FieldElement::reduce( self.square_inner())
}
/// Square this field element and multiply the result by 2.
///
/// XXX explain why square2 exists vs square (overflow)
///
/// # Preconditions
///
/// * |f[i]| bounded by 1.65*2^26, 1.65*2^25, 1.65*2^26, 1.65*2^25, etc.
@ -733,12 +1049,26 @@ impl FieldElement {
///
/// See fe_mul.c in ref10 implementation for discussion of implementation
/// strategy.
#[cfg(not(feature="radix_51"))]
pub fn square2(&self) -> FieldElement {
let mut coeffs = self.square_inner();
for i in 0..10 {
for i in 0..self.0.len() {
coeffs[i] += coeffs[i];
}
FieldElement::combine_coeffs(&coeffs)
FieldElement::reduce(&coeffs)
}
/// Compute `2 * self^2`.
#[cfg(feature="radix_51")]
pub fn square2(&self) -> FieldElement {
let mut limbs = self.square_inner();
// For this to work, need to have 1 extra bit of headroom after carry
// --> max 53 bit inputs, not 54
limbs[0] *= 2;
limbs[1] *= 2;
limbs[2] *= 2;
limbs[3] *= 2;
limbs[4] *= 2;
FieldElement::reduce(limbs)
}
#[inline]
@ -954,43 +1284,55 @@ mod test {
0x15, 0x21, 0xf9, 0xe3, 0xe1, 0x61, 0x21, 0x55];
#[test]
fn fieldelement_a_mul_a() {
fn a_mul_a_vs_a_squared_constant() {
let a = FieldElement::from_bytes(&A_BYTES);
let asq = FieldElement::from_bytes(&ASQ_BYTES);
assert_eq!(asq, &a*&a);
}
#[test]
fn a_square_vs_a_squared_constant() {
let a = FieldElement::from_bytes(&A_BYTES);
let asq = FieldElement::from_bytes(&ASQ_BYTES);
assert_eq!(asq, a.square());
}
#[test]
fn fieldelement_a_square2() {
fn a_square2_vs_a_squared_constant() {
let a = FieldElement::from_bytes(&A_BYTES);
let asq = FieldElement::from_bytes(&ASQ_BYTES);
assert_eq!(a.square2(), &asq+&asq);
}
#[test]
fn fieldelement_a_inv() {
fn a_invert_vs_inverse_of_a_constant() {
let a = FieldElement::from_bytes(&A_BYTES);
let ainv = FieldElement::from_bytes(&AINV_BYTES);
assert_eq!(ainv, a.invert());
let should_be_inverse = a.invert();
assert_eq!(ainv, should_be_inverse);
assert_eq!(FieldElement::one(), &a * &should_be_inverse);
}
#[test]
fn fieldelement_a_p58() {
fn a_p58_vs_ap58_constant() {
let a = FieldElement::from_bytes(&A_BYTES);
let ap58 = FieldElement::from_bytes(&AP58_BYTES);
assert_eq!(ap58, a.pow_p58());
}
#[test]
fn fieldelement_a_chi() {
fn chi_on_square_and_nonsquare() {
let a = FieldElement::from_bytes(&A_BYTES);
// a is square
assert_eq!(a.chi(), FieldElement::one());
let mut two_bytes = [0u8; 32]; two_bytes[0] = 2;
let two = FieldElement::from_bytes(&two_bytes);
// 2 is nonsquare
assert_eq!(two.chi(), FieldElement::minus_one());
}
#[test]
fn fieldelement_eq() {
fn equality() {
let a = FieldElement::from_bytes(&A_BYTES);
let ainv = FieldElement::from_bytes(&AINV_BYTES);
assert!(a == a);
@ -1000,44 +1342,49 @@ mod test {
/// Notice that the last element has the high bit set, which
/// should be ignored
static B_BYTES: [u8;32] =
[113, 191, 169, 143, 91, 234, 121, 15, 241, 131, 217, 36, 230, 101, 92, 234, 8, 208, 170, 251, 97, 127, 70, 210, 58, 23, 166, 87, 240, 169, 184, 178];
static B_LIMBS: FieldElement = FieldElement(
[-5652623, 8034020, 8266223, -13556020, -5672552, -5582839, -12603138, 15161929, -16418207, 13296296]);
[113, 191, 169, 143, 91, 234, 121, 15,
241, 131, 217, 36, 230, 101, 92, 234,
8, 208, 170, 251, 97, 127, 70, 210,
58, 23, 166, 87, 240, 169, 184, 178];
#[test]
fn fieldelement_frombytes_highbit_is_ignored() {
fn from_bytes_highbit_is_ignored() {
let mut cleared_bytes = B_BYTES.clone();
cleared_bytes[31] &= 127u8;
let orig_elt = FieldElement::from_bytes(&B_BYTES);
let cleared_elt = FieldElement::from_bytes(&cleared_bytes);
for i in 0..10 {
assert!(orig_elt[i] == cleared_elt[i]);
}
let with_highbit_set = FieldElement::from_bytes(&B_BYTES);
let without_highbit_set = FieldElement::from_bytes(&cleared_bytes);
assert_eq!(without_highbit_set, with_highbit_set);
}
#[cfg(not(feature="radix_51"))]
static B_LIMBS_RADIX_25_5: FieldElement = FieldElement(
[-5652623, 8034020, 8266223, -13556020, -5672552,
-5582839, -12603138, 15161929, -16418207, 13296296]);
#[cfg(not(feature="radix_51"))]
#[test]
fn fieldelement_to_bytes() {
fn from_bytes_vs_radix_25_5_limb_constants() {
let test_elt = FieldElement::from_bytes(&B_BYTES);
for i in 0..10 {
assert!(test_elt[i] == B_LIMBS[i]);
assert!(test_elt[i] == B_LIMBS_RADIX_25_5[i]);
}
}
#[cfg(not(feature="radix_51"))]
#[test]
fn fieldelement_from_bytes() {
let test_bytes = B_LIMBS.to_bytes();
fn radix_25_5_limb_constants_to_bytes_vs_byte_constants() {
let test_bytes = B_LIMBS_RADIX_25_5.to_bytes();
for i in 0..31 {
assert!(test_bytes[i] == B_BYTES[i]);
}
// high bit is set to zero in to_bytes
// Check that high bit is set to zero in to_bytes
assert!(test_bytes[31] == (B_BYTES[31] & 127u8));
}
#[test]
fn conditional_negate() {
let one = FieldElement([ 1,0,0,0,0,0,0,0,0,0]);
let minus_one = FieldElement([-1,0,0,0,0,0,0,0,0,0]);
let one = FieldElement::one();
let minus_one = FieldElement::minus_one();
let mut x = one;
x.conditional_negate(1u8);
assert_eq!(x, minus_one);

1
src/lib.rs
View file

@ -12,6 +12,7 @@
#![cfg_attr(not(feature = "std"), no_std)]
#![cfg_attr(not(feature = "std"), feature(collections))]
#![cfg_attr(feature = "nightly", feature(box_syntax))]
#![cfg_attr(feature = "nightly", feature(i128_type))]
#![allow(unused_features)]
#![cfg_attr(feature = "bench", feature(test))]
#![deny(missing_docs)] // refuse to compile if documentation is missing

14
src/utils.rs
View file

@ -29,3 +29,17 @@ pub fn load4(input: &[u8]) -> i64 {
| ((input[2] as i64) << 16)
| ((input[3] as i64) << 24)
}
/// Convert an array of (at least) eight bytes into a u64.
#[inline]
//#[allow(dead_code)]
pub fn load8(input: &[u8]) -> u64 {
(input[0] as u64)
| ((input[1] as u64) << 8)
| ((input[2] as u64) << 16)
| ((input[3] as u64) << 24)
| ((input[4] as u64) << 32)
| ((input[5] as u64) << 40)
| ((input[6] as u64) << 48)
| ((input[7] as u64) << 56)
}