saymrwulf/betrusted-curve25519-dalek-source
1
0
Fork 0
mirror of https://github.com/saymrwulf/betrusted-curve25519-dalek-source.git synced 2026-09-12 21:40:32 +00:00
Code Issues Projects Releases Packages Wiki Activity Actions

Merge remote-tracking branch 'chain/feature/refactor-32bit_r2' into develop

Browse source
This commit is contained in:
Isis Lovecruft 2017-11-18 00:40:07 +00:00
commit 0d59f6e709
Failed to extract signature
13 changed files with 506 additions and 3279 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

15
Cargo.toml
View file

@ -15,6 +15,7 @@ exclude = [
".gitignore", ".gitignore",
".travis.yml", ".travis.yml",
] ]
build = "build.rs"
[package.metadata.docs.rs] [package.metadata.docs.rs]
rustdoc-args = ["--html-in-header", ".cargo/registry/src/github.com-1ecc6299db9ec823/curve25519-dalek-0.13.2/rustdoc-include-katex-header.html"] rustdoc-args = ["--html-in-header", ".cargo/registry/src/github.com-1ecc6299db9ec823/curve25519-dalek-0.13.2/rustdoc-include-katex-header.html"]
@ -50,6 +51,17 @@ version = "0.6"
[dev-dependencies.serde_cbor] [dev-dependencies.serde_cbor]
version = "0.6" version = "0.6"
[build-dependencies]
subtle = "^0.3"
rand = "0.3"
generic-array = "^0.8"
digest = "0.6"
arrayref = "0.3.4"
[build-dependencies.serde]
version = "1.0"
optional = true
[features] [features]
nightly = ["radix_51", "subtle/nightly"] nightly = ["radix_51", "subtle/nightly"]
default = ["std"] default = ["std"]
@ -60,3 +72,6 @@ yolocrypto = []
bench = [] bench = []
# Radix-51 arithmetic using u128 # Radix-51 arithmetic using u128
radix_51 = [] radix_51 = []
# Include precomputed basepoint tables. This is off by default so that build.rs can generate the tables, and then re-enabled by build.rs in the main-stage compilation.
precomputed_tables = []

103
build.rs Normal file
View file

@ -0,0 +1,103 @@
#![cfg_attr(feature = "nightly", feature(i128_type))]
#![allow(unused_variables)]
#![allow(non_snake_case)]
#![allow(dead_code)]
extern crate core;
extern crate subtle;
extern crate rand;
extern crate digest;
extern crate generic_array;
#[macro_use]
extern crate arrayref;
use std::env;
use std::fs::File;
use std::io::Write;
use std::path::Path;
// Replicate lib.rs in the build.rs, since we're effectively building the whole crate twice.
//
// This should be fixed up by refactoring our code to seperate the "minimal" parts from the rest.
//
// For instance, this shouldn't exist here at all, but it does.
#[cfg(feature = "serde")]
extern crate serde;
#[path="src/field.rs"]
mod field;
#[cfg(not(feature="radix_51"))]
#[path="src/field_32bit.rs"]
mod field_32bit;
#[cfg(feature="radix_51")]
#[path="src/field_64bit.rs"]
mod field_64bit;
#[path="src/scalar.rs"]
mod scalar;
#[cfg(not(feature="radix_51"))]
#[path="src/scalar_32bit.rs"]
mod scalar_32bit;
#[cfg(feature="radix_51")]
#[path="src/scalar_64bit.rs"]
mod scalar_64bit;
#[path="src/montgomery.rs"]
mod montgomery;
#[path="src/edwards.rs"]
mod edwards;
#[path="src/ristretto.rs"]
mod ristretto;
#[path="src/constants.rs"]
mod constants;
#[cfg(not(feature="radix_51"))]
#[path="src/constants_32bit.rs"]
mod constants_32bit;
#[cfg(feature="radix_51")]
#[path="src/constants_64bit.rs"]
mod constants_64bit;
use edwards::EdwardsBasepointTable;
fn main() {
// Enable the "precomputed_tables" feature in the main build stage
println!("cargo:rustc-cfg=feature=\"precomputed_tables\"\n");
let out_dir = env::var("OUT_DIR").unwrap();
let dest_path = Path::new(&out_dir).join("basepoint_table.rs");
let mut f = File::create(&dest_path).unwrap();
// Generate a table of precomputed multiples of the basepoint
let table = EdwardsBasepointTable::create(&constants::ED25519_BASEPOINT_POINT);
f.write_all(format!("\n
#[cfg(feature=\"radix_51\")]
use field_64bit::FieldElement64;
#[cfg(not(feature=\"radix_51\"))]
use field_32bit::FieldElement32;
use edwards::AffineNielsPoint;
use edwards::EdwardsBasepointTable;
/// Table containing precomputed multiples of the basepoint `B = (x,4/5)`.
///
/// The table is defined so `constants::base[i][j-1] = j*(16^2i)*B`,
/// for `0 ≤ i < 32`, `1 ≤ j < 9`.
pub const ED25519_BASEPOINT_TABLE: EdwardsBasepointTable = {:?};
\n\n", &table).as_bytes()).unwrap();
// Now generate AFFINE_ODD_MULTIPLES_OF_BASEPOINT
let B = &constants::ED25519_BASEPOINT_POINT;
let B2 = B.double();
let mut odd_multiples = [B.to_affine_niels(); 8];
for i in 0..7 {
odd_multiples[i+1] = (&B2 + &odd_multiples[i]).to_extended().to_affine_niels();
}
f.write_all(format!("\n
/// Odd multiples of the basepoint `[B, 3B, 5B, 7B, 9B, 11B, 13B, 15B]`.
pub(crate) const AFFINE_ODD_MULTIPLES_OF_BASEPOINT: [AffineNielsPoint; 8] = {:?};
\n\n", &odd_multiples).as_bytes()).unwrap();
}

17
src/constants.rs
View file

@ -30,7 +30,8 @@
#![allow(non_snake_case)] #![allow(non_snake_case)]
use edwards::CompressedEdwardsY; use edwards::CompressedEdwardsY;
use ristretto::{RistrettoPoint, RistrettoBasepointTable}; use ristretto::RistrettoPoint;
use montgomery::CompressedMontgomeryU; use montgomery::CompressedMontgomeryU;
use scalar::Scalar; use scalar::Scalar;
@ -86,7 +87,15 @@ pub const BASEPOINT_ORDER_MINUS_2: Scalar = Scalar([
0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x10, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x10,
]); ]);
// Precomputed basepoint table is generated into a file by build.rs
#[cfg(feature="precomputed_tables")]
include!(concat!(env!("OUT_DIR"), "/basepoint_table.rs"));
#[cfg(feature="precomputed_tables")]
use ristretto::RistrettoBasepointTable;
/// The Ed25519 basepoint, as a RistrettoPoint /// The Ed25519 basepoint, as a RistrettoPoint
#[cfg(feature="precomputed_tables")]
pub const RISTRETTO_BASEPOINT_TABLE: RistrettoBasepointTable pub const RISTRETTO_BASEPOINT_TABLE: RistrettoBasepointTable
= RistrettoBasepointTable(ED25519_BASEPOINT_TABLE); = RistrettoBasepointTable(ED25519_BASEPOINT_TABLE);
@ -142,7 +151,7 @@ mod test {
#[cfg(not(feature="radix_51"))] #[cfg(not(feature="radix_51"))]
fn sqrt_minus_aplus2() { fn sqrt_minus_aplus2() {
use field_32bit::FieldElement32; use field_32bit::FieldElement32;
let minus_aplus2 = FieldElement32([-486664,0,0,0,0,0,0,0,0,0]); let minus_aplus2 = -&FieldElement32([486664,0,0,0,0,0,0,0,0,0]);
let sqrt = constants::SQRT_MINUS_APLUS2; let sqrt = constants::SQRT_MINUS_APLUS2;
let sq = &sqrt * &sqrt; let sq = &sqrt * &sqrt;
assert_eq!(sq, minus_aplus2); assert_eq!(sq, minus_aplus2);
@ -172,8 +181,8 @@ mod test {
#[test] #[test]
fn test_d_vs_ratio() { fn test_d_vs_ratio() {
use field_32bit::FieldElement32; use field_32bit::FieldElement32;
let a = FieldElement32([-121665,0,0,0,0,0,0,0,0,0]); let a = -&FieldElement32([121665,0,0,0,0,0,0,0,0,0]);
let b = FieldElement32([ 121666,0,0,0,0,0,0,0,0,0]); let b = FieldElement32([121666,0,0,0,0,0,0,0,0,0]);
let d = &a * &b.invert(); let d = &a * &b.invert();
let d2 = &d + &d; let d2 = &d + &d;
assert_eq!(d, constants::EDWARDS_D); assert_eq!(d, constants::EDWARDS_D);

1475
src/constants_32bit.rs
View file

File diff suppressed because it is too large Load diff

1437
src/constants_64bit.rs
View file

File diff suppressed because it is too large Load diff

37
src/edwards.rs
View file

@ -846,7 +846,7 @@ impl<'a, 'b> Mul<&'b ExtendedPoint> for &'a Scalar {
#[cfg(any(feature = "alloc", feature = "std"))] #[cfg(any(feature = "alloc", feature = "std"))]
pub fn multiscalar_mult<'a, 'b, I, J>(scalars: I, points: J) -> ExtendedPoint pub fn multiscalar_mult<'a, 'b, I, J>(scalars: I, points: J) -> ExtendedPoint
where I: IntoIterator<Item = &'a Scalar>, where I: IntoIterator<Item = &'a Scalar>,
J: IntoIterator<Item = &'b ExtendedPoint> J: IntoIterator<Item = &'b ExtendedPoint>
{ {
//assert_eq!(scalars.len(), points.len()); //assert_eq!(scalars.len(), points.len());
@ -1111,39 +1111,50 @@ impl ExtendedPoint {
impl Debug for ExtendedPoint { impl Debug for ExtendedPoint {
fn fmt(&self, f: &mut ::core::fmt::Formatter) -> ::core::fmt::Result { fn fmt(&self, f: &mut ::core::fmt::Formatter) -> ::core::fmt::Result {
write!(f, "ExtendedPoint(\n\tX: {:?},\n\tY: {:?},\n\tZ: {:?},\n\tT: {:?}\n)", write!(f, "ExtendedPoint{{\n\tX: {:?},\n\tY: {:?},\n\tZ: {:?},\n\tT: {:?}\n}}",
&self.X, &self.Y, &self.Z, &self.T) &self.X, &self.Y, &self.Z, &self.T)
} }
} }
impl Debug for ProjectivePoint { impl Debug for ProjectivePoint {
fn fmt(&self, f: &mut ::core::fmt::Formatter) -> ::core::fmt::Result { fn fmt(&self, f: &mut ::core::fmt::Formatter) -> ::core::fmt::Result {
write!(f, "ProjectivePoint(\n\tX: {:?},\n\tY: {:?},\n\tZ: {:?}\n)", write!(f, "ProjectivePoint{{\n\tX: {:?},\n\tY: {:?},\n\tZ: {:?}\n}}",
&self.X, &self.Y, &self.Z) &self.X, &self.Y, &self.Z)
} }
} }
impl Debug for CompletedPoint { impl Debug for CompletedPoint {
fn fmt(&self, f: &mut ::core::fmt::Formatter) -> ::core::fmt::Result { fn fmt(&self, f: &mut ::core::fmt::Formatter) -> ::core::fmt::Result {
write!(f, "CompletedPoint(\n\tX: {:?},\n\tY: {:?},\n\tZ: {:?},\n\tT: {:?}\n)", write!(f, "CompletedPoint{{\n\tX: {:?},\n\tY: {:?},\n\tZ: {:?},\n\tT: {:?}\n}}",
&self.X, &self.Y, &self.Z, &self.T) &self.X, &self.Y, &self.Z, &self.T)
} }
} }
impl Debug for AffineNielsPoint { impl Debug for AffineNielsPoint {
fn fmt(&self, f: &mut ::core::fmt::Formatter) -> ::core::fmt::Result { fn fmt(&self, f: &mut ::core::fmt::Formatter) -> ::core::fmt::Result {
write!(f, "AffineNielsPoint(\n\ty_plus_x: {:?},\n\ty_minus_x: {:?},\n\txy2d: {:?}\n)", write!(f, "AffineNielsPoint{{\n\ty_plus_x: {:?},\n\ty_minus_x: {:?},\n\txy2d: {:?}\n}}",
&self.y_plus_x, &self.y_minus_x, &self.xy2d) &self.y_plus_x, &self.y_minus_x, &self.xy2d)
} }
} }
impl Debug for ProjectiveNielsPoint { impl Debug for ProjectiveNielsPoint {
fn fmt(&self, f: &mut ::core::fmt::Formatter) -> ::core::fmt::Result { fn fmt(&self, f: &mut ::core::fmt::Formatter) -> ::core::fmt::Result {
write!(f, "ProjectiveNielsPoint(\n\tY_plus_X: {:?},\n\tY_minus_X: {:?},\n\tZ: {:?},\n\tT2d: {:?}\n)", write!(f, "ProjectiveNielsPoint{{\n\tY_plus_X: {:?},\n\tY_minus_X: {:?},\n\tZ: {:?},\n\tT2d: {:?}\n}}",
&self.Y_plus_X, &self.Y_minus_X, &self.Z, &self.T2d) &self.Y_plus_X, &self.Y_minus_X, &self.Z, &self.T2d)
} }
} }
impl Debug for EdwardsBasepointTable {
fn fmt(&self, f: &mut ::core::fmt::Formatter) -> ::core::fmt::Result {
write!(f, "EdwardsBasepointTable([\n")?;
for i in 0..32 {
write!(f, "\t{:?},\n", &self.0[i])?;
}
write!(f, "])")
}
}
// ------------------------------------------------------------------------ // ------------------------------------------------------------------------
// Variable-time functions // Variable-time functions
// ------------------------------------------------------------------------ // ------------------------------------------------------------------------
@ -1217,6 +1228,7 @@ pub mod vartime {
/// Given a point `A` and scalars `a` and `b`, compute the point /// Given a point `A` and scalars `a` and `b`, compute the point
/// `aA+bB`, where `B` is the Ed25519 basepoint (i.e., `B = (x,4/5)` /// `aA+bB`, where `B` is the Ed25519 basepoint (i.e., `B = (x,4/5)`
/// with x positive). /// with x positive).
#[cfg(feature="precomputed_tables")]
pub fn double_scalar_mult_basepoint(a: &Scalar, pub fn double_scalar_mult_basepoint(a: &Scalar,
A: &ExtendedPoint, A: &ExtendedPoint,
b: &Scalar) -> ExtendedPoint { b: &Scalar) -> ExtendedPoint {
@ -1354,6 +1366,7 @@ mod test {
/// Test that computing 1*basepoint gives the correct basepoint. /// Test that computing 1*basepoint gives the correct basepoint.
#[test] #[test]
#[cfg(feature="precomputed_tables")]
fn basepoint_mult_one_vs_basepoint() { fn basepoint_mult_one_vs_basepoint() {
let bp = &constants::ED25519_BASEPOINT_TABLE * &Scalar::one(); let bp = &constants::ED25519_BASEPOINT_TABLE * &Scalar::one();
let compressed = bp.compress(); let compressed = bp.compress();
@ -1362,6 +1375,7 @@ mod test {
/// Test that `EdwardsBasepointTable::basepoint()` gives the correct basepoint. /// Test that `EdwardsBasepointTable::basepoint()` gives the correct basepoint.
#[test] #[test]
#[cfg(feature="precomputed_tables")]
fn basepoint_table_basepoint_function_correct() { fn basepoint_table_basepoint_function_correct() {
let bp = constants::ED25519_BASEPOINT_TABLE.basepoint(); let bp = constants::ED25519_BASEPOINT_TABLE.basepoint();
assert_eq!(bp.compress(), constants::BASE_CMPRSSD); assert_eq!(bp.compress(), constants::BASE_CMPRSSD);
@ -1412,6 +1426,7 @@ mod test {
/// Sanity check for conversion to precomputed points /// Sanity check for conversion to precomputed points
#[test] #[test]
#[cfg(feature="precomputed_tables")]
fn to_affine_niels_clears_denominators() { fn to_affine_niels_clears_denominators() {
// construct a point as aB so it has denominators (ie. Z != 1) // construct a point as aB so it has denominators (ie. Z != 1)
let aB = &constants::ED25519_BASEPOINT_TABLE * &A_SCALAR; let aB = &constants::ED25519_BASEPOINT_TABLE * &A_SCALAR;
@ -1423,6 +1438,7 @@ mod test {
/// Test basepoint_mult versus a known scalar multiple from ed25519.py /// Test basepoint_mult versus a known scalar multiple from ed25519.py
#[test] #[test]
#[cfg(feature="precomputed_tables")]
fn basepoint_mult_vs_ed25519py() { fn basepoint_mult_vs_ed25519py() {
let aB = &constants::ED25519_BASEPOINT_TABLE * &A_SCALAR; let aB = &constants::ED25519_BASEPOINT_TABLE * &A_SCALAR;
assert_eq!(aB.compress(), A_TIMES_BASEPOINT); assert_eq!(aB.compress(), A_TIMES_BASEPOINT);
@ -1430,6 +1446,7 @@ mod test {
/// Test that multiplication by the basepoint order kills the basepoint /// Test that multiplication by the basepoint order kills the basepoint
#[test] #[test]
#[cfg(feature="precomputed_tables")]
fn basepoint_mult_by_basepoint_order() { fn basepoint_mult_by_basepoint_order() {
let B = &constants::ED25519_BASEPOINT_TABLE; let B = &constants::ED25519_BASEPOINT_TABLE;
let should_be_id = B * &constants::BASEPOINT_ORDER; let should_be_id = B * &constants::BASEPOINT_ORDER;
@ -1462,6 +1479,7 @@ mod test {
/// Test that computing 2*basepoint is the same as basepoint.double() /// Test that computing 2*basepoint is the same as basepoint.double()
#[test] #[test]
#[cfg(feature="precomputed_tables")]
fn basepoint_mult_two_vs_basepoint2() { fn basepoint_mult_two_vs_basepoint2() {
let mut two_bytes = [0u8; 32]; two_bytes[0] = 2; let mut two_bytes = [0u8; 32]; two_bytes[0] = 2;
let bp2 = &constants::ED25519_BASEPOINT_TABLE * &Scalar(two_bytes); let bp2 = &constants::ED25519_BASEPOINT_TABLE * &Scalar(two_bytes);
@ -1555,6 +1573,7 @@ mod test {
/// Test double_scalar_mult_vartime vs ed25519.py /// Test double_scalar_mult_vartime vs ed25519.py
#[test] #[test]
#[cfg(feature="precomputed_tables")]
fn double_scalar_mult_basepoint_vs_ed25519py() { fn double_scalar_mult_basepoint_vs_ed25519py() {
let A = A_TIMES_BASEPOINT.decompress().unwrap(); let A = A_TIMES_BASEPOINT.decompress().unwrap();
let result = vartime::double_scalar_mult_basepoint(&A_SCALAR, &A, &B_SCALAR); let result = vartime::double_scalar_mult_basepoint(&A_SCALAR, &A, &B_SCALAR);
@ -1635,6 +1654,7 @@ mod bench {
} }
#[bench] #[bench]
#[cfg(feature="precomputed_tables")]
fn basepoint_mult(b: &mut Bencher) { fn basepoint_mult(b: &mut Bencher) {
let B = &constants::ED25519_BASEPOINT_TABLE; let B = &constants::ED25519_BASEPOINT_TABLE;
b.iter(|| B * &A_SCALAR); b.iter(|| B * &A_SCALAR);
@ -1647,6 +1667,7 @@ mod bench {
} }
#[bench] #[bench]
#[cfg(feature="precomputed_tables")]
fn bench_select_precomputed_point(b: &mut Bencher) { fn bench_select_precomputed_point(b: &mut Bencher) {
b.iter(|| select_precomputed_point(0, &constants::ED25519_BASEPOINT_TABLE.0[0])); b.iter(|| select_precomputed_point(0, &constants::ED25519_BASEPOINT_TABLE.0[0]));
} }
@ -1704,14 +1725,15 @@ mod bench {
b.iter(|| p1.mult_by_cofactor()); b.iter(|| p1.mult_by_cofactor());
} }
#[cfg(feature="basepoint_table_creation")]
#[bench] #[bench]
#[cfg(feature="basepoint_table_creation")]
fn create_basepoint_table(b: &mut Bencher) { fn create_basepoint_table(b: &mut Bencher) {
let aB = &constants::ED25519_BASEPOINT_TABLE * &A_SCALAR; let aB = &constants::ED25519_BASEPOINT_TABLE * &A_SCALAR;
b.iter(|| EdwardsBasepointTable::create(&aB)); b.iter(|| EdwardsBasepointTable::create(&aB));
} }
#[bench] #[bench]
#[cfg(feature="basepoint_table_creation")]
fn ten_fold_scalar_mult(b: &mut Bencher) { fn ten_fold_scalar_mult(b: &mut Bencher) {
let mut csprng: OsRng = OsRng::new().unwrap(); let mut csprng: OsRng = OsRng::new().unwrap();
// Create 10 random scalars // Create 10 random scalars
@ -1735,6 +1757,7 @@ mod bench {
} }
#[bench] #[bench]
#[cfg(feature="basepoint_table_creation")]
fn ten_fold_scalar_mult(b: &mut Bencher) { fn ten_fold_scalar_mult(b: &mut Bencher) {
let mut csprng: OsRng = OsRng::new().unwrap(); let mut csprng: OsRng = OsRng::new().unwrap();
// Create 10 random scalars // Create 10 random scalars

23
src/field.rs
View file

@ -404,29 +404,6 @@ mod test {
assert_eq!(without_highbit_set, with_highbit_set); assert_eq!(without_highbit_set, with_highbit_set);
} }
#[cfg(not(feature="radix_51"))]
static B_LIMBS_RADIX_25_5: FieldElement32 = FieldElement32(
[-5652623, 8034020, 8266223, -13556020, -5672552,
-5582839, -12603138, 15161929, -16418207, 13296296]);
#[cfg(not(feature="radix_51"))]
#[test]
fn from_bytes_vs_radix_25_5_limb_constants() {
let test_elt = FieldElement::from_bytes(&B_BYTES);
assert_eq!(test_elt.0, B_LIMBS_RADIX_25_5.0);
}
#[cfg(not(feature="radix_51"))]
#[test]
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]);
}
// Check that high bit is set to zero in to_bytes
assert!(test_bytes[31] == (B_BYTES[31] & 127u8));
}
#[test] #[test]
fn conditional_negate() { fn conditional_negate() {
let one = FieldElement::one(); let one = FieldElement::one();

598
src/field_32bit.rs
View file

@ -11,8 +11,8 @@
//! Field arithmetic for ℤ/(2²⁵⁵-19), using 32-bit arithmetic with //! Field arithmetic for ℤ/(2²⁵⁵-19), using 32-bit arithmetic with
//! 64-bit products. //! 64-bit products.
//! //!
//! Based on Adam Langley's curve25519-donna and (Golang) ed25519 //! This code was originally derived from Adam Langley's
//! implementations. //! curve25519-donna and (Golang) ed25519 implementations.
//! //!
//! This implementation is intended for platforms that can multiply //! This implementation is intended for platforms that can multiply
//! 32-bit inputs to produce 64-bit outputs. //! 32-bit inputs to produce 64-bit outputs.
@ -30,14 +30,11 @@ use core::ops::Neg;
use subtle::ConditionallyAssignable; use subtle::ConditionallyAssignable;
use utils::{load3, load4};
/// A `FieldElement32` represents an element of the field GF(2^255 - 19). /// A `FieldElement32` represents an element of the field GF(2^255 - 19).
/// ///
/// In the 32-bit implementation, a `FieldElement32` is represented in /// In the 32-bit implementation, a `FieldElement32` is represented in
/// radix 2^25.5 as ten `i32`s, so that an element t, entries /// radix 2^25.5 as ten `u32`s, so that an element t, entries
/// t[0],...,t[9], represents the integer t[0]+2^26 t[1]+2^51 /// t[0],...,t[9], represents `sum(t[i]*2^ceil(i*51/2))`.
/// t[2]+2^77 t[3]+2^102 t[4]+...+2^230 t[9].
/// ///
/// The coefficients t[i] are allowed to grow between multiplications. /// The coefficients t[i] are allowed to grow between multiplications.
/// ///
@ -56,11 +53,11 @@ use utils::{load3, load4};
/// faster. However, the `FieldElement64` implementation requires Rust's /// faster. However, the `FieldElement64` implementation requires Rust's
/// `u128`, which is not yet stable. /// `u128`, which is not yet stable.
#[derive(Copy, Clone)] #[derive(Copy, Clone)]
pub struct FieldElement32(pub (crate) [i32; 10]); pub struct FieldElement32(pub (crate) [u32; 10]);
impl Debug for FieldElement32 { impl Debug for FieldElement32 {
fn fmt(&self, f: &mut ::core::fmt::Formatter) -> ::core::fmt::Result { fn fmt(&self, f: &mut ::core::fmt::Formatter) -> ::core::fmt::Result {
write!(f, "FieldElement32: {:?}", &self.0[..]) write!(f, "FieldElement32({:?})", &self.0[..])
} }
} }
@ -83,9 +80,22 @@ impl<'a, 'b> Add<&'b FieldElement32> for &'a FieldElement32 {
impl<'b> SubAssign<&'b FieldElement32> for FieldElement32 { impl<'b> SubAssign<&'b FieldElement32> for FieldElement32 {
fn sub_assign(&mut self, _rhs: &'b FieldElement32) { fn sub_assign(&mut self, _rhs: &'b FieldElement32) {
for i in 0..10 { // See comment in FieldElement64::Sub
self.0[i] -= _rhs.0[i]; //
} // Compute a - b as ((a + 2^4 * p) - b) to avoid underflow.
let b = &_rhs.0;
self.0 = FieldElement32::reduce([
((self.0[0] + (0x3ffffed << 4)) - b[0]) as u64,
((self.0[1] + (0x1ffffff << 4)) - b[1]) as u64,
((self.0[2] + (0x3ffffff << 4)) - b[2]) as u64,
((self.0[3] + (0x1ffffff << 4)) - b[3]) as u64,
((self.0[4] + (0x3ffffff << 4)) - b[4]) as u64,
((self.0[5] + (0x1ffffff << 4)) - b[5]) as u64,
((self.0[6] + (0x3ffffff << 4)) - b[6]) as u64,
((self.0[7] + (0x1ffffff << 4)) - b[7]) as u64,
((self.0[8] + (0x3ffffff << 4)) - b[8]) as u64,
((self.0[9] + (0x1ffffff << 4)) - b[9]) as u64,
]).0;
} }
} }
@ -108,86 +118,101 @@ impl<'b> MulAssign<&'b FieldElement32> for FieldElement32 {
impl<'a, 'b> Mul<&'b FieldElement32> for &'a FieldElement32 { impl<'a, 'b> Mul<&'b FieldElement32> for &'a FieldElement32 {
type Output = FieldElement32; type Output = FieldElement32;
fn mul(self, _rhs: &'b FieldElement32) -> FieldElement32 { fn mul(self, _rhs: &'b FieldElement32) -> FieldElement32 {
// Notes preserved from ed25519.go (presumably originally from ref10): /// Helper function to multiply two 32-bit integers with 64 bits
// /// of output.
// Calculates h = f * g. Can overlap h with f or g. #[inline(always)]
// fn m(x: u32, y: u32) -> u64 { (x as u64) * (y as u64) }
// # 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[0] as i64;
let f1 = self.0[1] as i64;
let f2 = self.0[2] as i64;
let f3 = self.0[3] as i64;
let f4 = self.0[4] as i64;
let f5 = self.0[5] as i64;
let f6 = self.0[6] as i64;
let f7 = self.0[7] as i64;
let f8 = self.0[8] as i64;
let f9 = self.0[9] as i64;
let f1_2 = (2 * self.0[1]) as i64; // Alias self, _rhs for more readable formulas
let f3_2 = (2 * self.0[3]) as i64; let x: &[u32;10] = &self.0; let y: &[u32;10] = &_rhs.0;
let f5_2 = (2 * self.0[5]) as i64;
let f7_2 = (2 * self.0[7]) as i64;
let f9_2 = (2 * self.0[9]) as i64;
let g0 = _rhs.0[0] as i64; // We assume that the input limbs x[i], y[i] are bounded by:
let g1 = _rhs.0[1] as i64; //
let g2 = _rhs.0[2] as i64; // x[i], y[i] < 2^(26 + b) if i even
let g3 = _rhs.0[3] as i64; // x[i], y[i] < 2^(25 + b) if i odd
let g4 = _rhs.0[4] as i64; //
let g5 = _rhs.0[5] as i64; // where b is a (real) parameter representing the excess bits of
let g6 = _rhs.0[6] as i64; // the limbs. We track the bitsizes of all variables through
let g7 = _rhs.0[7] as i64; // the computation and solve at the end for the allowable
let g8 = _rhs.0[8] as i64; // headroom bitsize b (which determines how many additions we
let g9 = _rhs.0[9] as i64; // can perform between reductions or multiplications).
let g1_19 = (19 * _rhs.0[1]) as i64; /* 1.4*2^29 */ let y1_19 = 19 * y[1]; // This fits in a u32
let g2_19 = (19 * _rhs.0[2]) as i64; /* 1.4*2^30; still ok */ let y2_19 = 19 * y[2]; // iff 26 + b + lg(19) < 32
let g3_19 = (19 * _rhs.0[3]) as i64; let y3_19 = 19 * y[3]; // if b < 32 - 26 - 4.248 = 1.752
let g4_19 = (19 * _rhs.0[4]) as i64; let y4_19 = 19 * y[4];
let g5_19 = (19 * _rhs.0[5]) as i64; let y5_19 = 19 * y[5]; // below, b<2.5: this is a bottleneck,
let g6_19 = (19 * _rhs.0[6]) as i64; let y6_19 = 19 * y[6]; // could be avoided by promoting to
let g7_19 = (19 * _rhs.0[7]) as i64; let y7_19 = 19 * y[7]; // u64 here instead of in m()
let g8_19 = (19 * _rhs.0[8]) as i64; let y8_19 = 19 * y[8];
let g9_19 = (19 * _rhs.0[9]) as i64; let y9_19 = 19 * y[9];
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; // What happens when we multiply x[i] with y[j] and place the
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; // result into the (i+j)-th limb?
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; // x[i] represents the value x[i]*2^ceil(i*51/2)
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; // y[j] represents the value y[j]*2^ceil(j*51/2)
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; // z[i+j] represents the value z[i+j]*2^ceil((i+j)*51/2)
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; // x[i]*y[j] represents the value x[i]*y[i]*2^(ceil(i*51/2)+ceil(j*51/2))
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; // Since the radix is already accounted for, the result placed
let h9 = f0*g9 + f1*g8 + f2*g7 + f3*g6 + f4*g5 + f5*g4 + f6*g3 + f7*g2 + f8*g1 + f9*g0; // into the (i+j)-th limb should be
//
// x[i]*y[i]*2^(ceil(i*51/2)+ceil(j*51/2) - ceil((i+j)*51/2)).
//
// The value of ceil(i*51/2)+ceil(j*51/2) - ceil((i+j)*51/2) is
// 1 when both i and j are odd, and 0 otherwise. So we add
//
// x[i]*y[j] if either i or j is even
// 2*x[i]*y[j] if i and j are both odd
//
// by using precomputed multiples of x[i] for odd i:
FieldElement32::reduce([h0, h1, h2, h3, h4, h5, h6, h7, h8, h9]) let x1_2 = 2 * x[1]; // This fits in a u32 iff 25 + b + 1 < 32
let x3_2 = 2 * x[3]; // iff b < 6
let x5_2 = 2 * x[5];
let x7_2 = 2 * x[7];
let x9_2 = 2 * x[9];
let z0 = m(x[0],y[0]) + m(x1_2,y9_19) + m(x[2],y8_19) + m(x3_2,y7_19) + m(x[4],y6_19) + m(x5_2,y5_19) + m(x[6],y4_19) + m(x7_2,y3_19) + m(x[8],y2_19) + m(x9_2,y1_19);
let z1 = m(x[0],y[1]) + m(x[1],y[0]) + m(x[2],y9_19) + m(x[3],y8_19) + m(x[4],y7_19) + m(x[5],y6_19) + m(x[6],y5_19) + m(x[7],y4_19) + m(x[8],y3_19) + m(x[9],y2_19);
let z2 = m(x[0],y[2]) + m(x1_2,y[1]) + m(x[2],y[0]) + m(x3_2,y9_19) + m(x[4],y8_19) + m(x5_2,y7_19) + m(x[6],y6_19) + m(x7_2,y5_19) + m(x[8],y4_19) + m(x9_2,y3_19);
let z3 = m(x[0],y[3]) + m(x[1],y[2]) + m(x[2],y[1]) + m(x[3],y[0]) + m(x[4],y9_19) + m(x[5],y8_19) + m(x[6],y7_19) + m(x[7],y6_19) + m(x[8],y5_19) + m(x[9],y4_19);
let z4 = m(x[0],y[4]) + m(x1_2,y[3]) + m(x[2],y[2]) + m(x3_2,y[1]) + m(x[4],y[0]) + m(x5_2,y9_19) + m(x[6],y8_19) + m(x7_2,y7_19) + m(x[8],y6_19) + m(x9_2,y5_19);
let z5 = m(x[0],y[5]) + m(x[1],y[4]) + m(x[2],y[3]) + m(x[3],y[2]) + m(x[4],y[1]) + m(x[5],y[0]) + m(x[6],y9_19) + m(x[7],y8_19) + m(x[8],y7_19) + m(x[9],y6_19);
let z6 = m(x[0],y[6]) + m(x1_2,y[5]) + m(x[2],y[4]) + m(x3_2,y[3]) + m(x[4],y[2]) + m(x5_2,y[1]) + m(x[6],y[0]) + m(x7_2,y9_19) + m(x[8],y8_19) + m(x9_2,y7_19);
let z7 = m(x[0],y[7]) + m(x[1],y[6]) + m(x[2],y[5]) + m(x[3],y[4]) + m(x[4],y[3]) + m(x[5],y[2]) + m(x[6],y[1]) + m(x[7],y[0]) + m(x[8],y9_19) + m(x[9],y8_19);
let z8 = m(x[0],y[8]) + m(x1_2,y[7]) + m(x[2],y[6]) + m(x3_2,y[5]) + m(x[4],y[4]) + m(x5_2,y[3]) + m(x[6],y[2]) + m(x7_2,y[1]) + m(x[8],y[0]) + m(x9_2,y9_19);
let z9 = m(x[0],y[9]) + m(x[1],y[8]) + m(x[2],y[7]) + m(x[3],y[6]) + m(x[4],y[5]) + m(x[5],y[4]) + m(x[6],y[3]) + m(x[7],y[2]) + m(x[8],y[1]) + m(x[9],y[0]);
// How big is the contribution to z[i+j] from x[i], y[j]?
//
// Using the bounds above, we get:
//
// i even, j even: x[i]*y[j] < 2^(26+b)*2^(26+b) = 2*2^(51+2*b)
// i odd, j even: x[i]*y[j] < 2^(25+b)*2^(26+b) = 1*2^(51+2*b)
// i even, j odd: x[i]*y[j] < 2^(26+b)*2^(25+b) = 1*2^(51+2*b)
// i odd, j odd: 2*x[i]*y[j] < 2*2^(25+b)*2^(25+b) = 1*2^(51+2*b)
//
// We perform inline reduction mod p by replacing 2^255 by 19
// (since 2^255 - 19 = 0 mod p). This adds a factor of 19, so
// we get the bounds (z0 is the biggest one, but calculated for
// posterity here in case finer estimation is needed later):
//
// z0 < ( 2 + 1*19 + 2*19 + 1*19 + 2*19 + 1*19 + 2*19 + 1*19 + 2*19 + 1*19 )*2^(51 + 2b) = 249*2^(51 + 2*b)
// z1 < ( 1 + 1 + 1*19 + 1*19 + 1*19 + 1*19 + 1*19 + 1*19 + 1*19 + 1*19 )*2^(51 + 2b) = 154*2^(51 + 2*b)
// z2 < ( 2 + 1 + 2 + 1*19 + 2*19 + 1*19 + 2*19 + 1*19 + 2*19 + 1*19 )*2^(51 + 2b) = 195*2^(51 + 2*b)
// z3 < ( 1 + 1 + 1 + 1 + 1*19 + 1*19 + 1*19 + 1*19 + 1*19 + 1*19 )*2^(51 + 2b) = 118*2^(51 + 2*b)
// z4 < ( 2 + 1 + 2 + 1 + 2 + 1*19 + 2*19 + 1*19 + 2*19 + 1*19 )*2^(51 + 2b) = 141*2^(51 + 2*b)
// z5 < ( 1 + 1 + 1 + 1 + 1 + 1 + 1*19 + 1*19 + 1*19 + 1*19 )*2^(51 + 2b) = 82*2^(51 + 2*b)
// z6 < ( 2 + 1 + 2 + 1 + 2 + 1 + 2 + 1*19 + 2*19 + 1*19 )*2^(51 + 2b) = 87*2^(51 + 2*b)
// z7 < ( 1 + 1 + 1 + 1 + 1 + 1 + 1 + 1 + 1*19 + 1*19 )*2^(51 + 2b) = 46*2^(51 + 2*b)
// z6 < ( 2 + 1 + 2 + 1 + 2 + 1 + 2 + 1 + 2 + 1*19 )*2^(51 + 2b) = 33*2^(51 + 2*b)
// z7 < ( 1 + 1 + 1 + 1 + 1 + 1 + 1 + 1 + 1 + 1 )*2^(51 + 2b) = 10*2^(51 + 2*b)
//
// So z[0] fits into a u64 if 51 + 2*b + lg(249) < 64
// if b < 2.5.
FieldElement32::reduce([z0, z1, z2, z3, z4, z5, z6, z7, z8, z9])
} }
} }
@ -202,7 +227,7 @@ impl<'a> Neg for &'a FieldElement32 {
impl ConditionallyAssignable for FieldElement32 { impl ConditionallyAssignable for FieldElement32 {
fn conditional_assign(&mut self, f: &FieldElement32, choice: u8) { fn conditional_assign(&mut self, f: &FieldElement32, choice: u8) {
let mask = -(choice as i32); let mask = (-(choice as i32)) as u32;
for i in 0..10 { for i in 0..10 {
self.0[i] ^= mask & (self.0[i] ^ f.0[i]); self.0[i] ^= mask & (self.0[i] ^ f.0[i]);
} }
@ -212,9 +237,20 @@ impl ConditionallyAssignable for FieldElement32 {
impl FieldElement32 { impl FieldElement32 {
/// Invert the sign of this field element /// Invert the sign of this field element
pub fn negate(&mut self) { pub fn negate(&mut self) {
for i in 0..10 { // Compute -b as ((2^4 * p) - b) to avoid underflow.
self.0[i] = -self.0[i]; let neg = FieldElement32::reduce([
} ((0x3ffffed << 4) - self.0[0]) as u64,
((0x1ffffff << 4) - self.0[1]) as u64,
((0x3ffffff << 4) - self.0[2]) as u64,
((0x1ffffff << 4) - self.0[3]) as u64,
((0x3ffffff << 4) - self.0[4]) as u64,
((0x1ffffff << 4) - self.0[5]) as u64,
((0x3ffffff << 4) - self.0[6]) as u64,
((0x1ffffff << 4) - self.0[7]) as u64,
((0x3ffffff << 4) - self.0[8]) as u64,
((0x1ffffff << 4) - self.0[9]) as u64,
]);
self.0 = neg.0;
} }
/// Construct zero. /// Construct zero.
@ -229,98 +265,64 @@ impl FieldElement32 {
/// Construct -1. /// Construct -1.
pub fn minus_one() -> FieldElement32 { pub fn minus_one() -> FieldElement32 {
FieldElement32([-1, 0, 0, 0, 0, 0, 0, 0, 0, 0 ]) FieldElement32([
0x3ffffec, 0x1ffffff, 0x3ffffff, 0x1ffffff, 0x3ffffff,
0x1ffffff, 0x3ffffff, 0x1ffffff, 0x3ffffff, 0x1ffffff,
])
} }
fn reduce(mut h: [i64; 10]) -> FieldElement32 { //FeCombine /// Given unreduced coefficients `z[0], ..., z[9]` of any size,
let mut c = [0i64; 10]; /// carry and reduce them mod p to obtain a `FieldElement32`
/// whose coefficients have excess `b < 0.007`.
///
/// In other words, each coefficient of the result is bounded by
/// either `2^(25 + 0.007)` or `2^(26 + 0.007)`, as appropriate.
fn reduce(mut z: [u64; 10]) -> FieldElement32 {
/* const LOW_25_BITS: u64 = (1 << 25) - 1;
|h[0]| <= (1.1*1.1*2^52*(1+19+19+19+19)+1.1*1.1*2^50*(38+38+38+38+38)) const LOW_26_BITS: u64 = (1 << 26) - 1;
i.e. |h[0]| <= 1.2*2^59; narrower ranges for h[2], h[4], h[6], h[8]
|h[1]| <= (1.1*1.1*2^51*(1+1+19+19+19+19+19+19+19+19))
i.e. |h[1]| <= 1.5*2^58; narrower ranges for h[3], h[5], h[7], h[9]
*/
c[0] = (h[0] + (1 << 25)) >> 26; /// Carry the value from limb i = 0..8 to limb i+1
h[1] += c[0]; #[inline(always)]
h[0] -= c[0] << 26; fn carry(z: &mut [u64; 10], i: usize) {
c[4] = (h[4] + (1 << 25)) >> 26; debug_assert!(i < 9);
h[5] += c[4]; if i % 2 == 0 {
h[4] -= c[4] << 26; // Even limbs have 26 bits
/* |h[0]| <= 2^25 */ z[i+1] += z[i] >> 26;
/* |h[4]| <= 2^25 */ z[i] &= LOW_26_BITS;
/* |h[1]| <= 1.51*2^58 */ } else {
/* |h[5]| <= 1.51*2^58 */ // Odd limbs have 25 bits
z[i+1] += z[i] >> 25;
z[i] &= LOW_25_BITS;
}
}
c[1] = (h[1] + (1 << 24)) >> 25; // Perform two halves of the carry chain in parallel.
h[2] += c[1]; carry(&mut z, 0); carry(&mut z, 4);
h[1] -= c[1] << 25; carry(&mut z, 1); carry(&mut z, 5);
c[5] = (h[5] + (1 << 24)) >> 25; carry(&mut z, 2); carry(&mut z, 6);
h[6] += c[5]; carry(&mut z, 3); carry(&mut z, 7);
h[5] -= c[5] << 25; // Since z[3] < 2^64, c < 2^(64-25) = 2^39,
/* |h[1]| <= 2^24; from now on fits into int32 */ // so z[4] < 2^26 + 2^39 < 2^39.0002
/* |h[5]| <= 2^24; from now on fits into int32 */ carry(&mut z, 4); carry(&mut z, 8);
/* |h[2]| <= 1.21*2^59 */ // Now z[4] < 2^26
/* |h[6]| <= 1.21*2^59 */ // and z[5] < 2^25 + 2^13.0002 < 2^25.0004 (good enough)
c[2] = (h[2] + (1 << 25)) >> 26; // Last carry has a multiplication by 19:
h[3] += c[2]; z[0] += 19*(z[9] >> 25);
h[2] -= c[2] << 26; z[9] &= LOW_25_BITS;
c[6] = (h[6] + (1 << 25)) >> 26;
h[7] += c[6];
h[6] -= c[6] << 26;
/* |h[2]| <= 2^25; from now on fits into int32 unchanged */
/* |h[6]| <= 2^25; from now on fits into int32 unchanged */
/* |h[3]| <= 1.51*2^58 */
/* |h[7]| <= 1.51*2^58 */
c[3] = (h[3] + (1 << 24)) >> 25; // Since z[9] < 2^64, c < 2^(64-25) = 2^39,
h[4] += c[3]; // so z[0] + 19*c < 2^26 + 2^43.248 < 2^43.249.
h[3] -= c[3] << 25; carry(&mut z, 0);
c[7] = (h[7] + (1 << 24)) >> 25; // Now z[1] < 2^25 - 2^(43.249 - 26)
h[8] += c[7]; // < 2^25.007 (good enough)
h[7] -= c[7] << 25; // and we're done.
/* |h[3]| <= 2^24; from now on fits into int32 unchanged */
/* |h[7]| <= 2^24; from now on fits into int32 unchanged */
/* |h[4]| <= 1.52*2^33 */
/* |h[8]| <= 1.52*2^33 */
c[4] = (h[4] + (1 << 25)) >> 26; FieldElement32([
h[5] += c[4]; z[0] as u32, z[1] as u32, z[2] as u32, z[3] as u32, z[4] as u32,
h[4] -= c[4] << 26; z[5] as u32, z[6] as u32, z[7] as u32, z[8] as u32, z[9] as u32,
c[8] = (h[8] + (1 << 25)) >> 26; ])
h[9] += c[8];
h[8] -= c[8] << 26;
/* |h[4]| <= 2^25; from now on fits into int32 unchanged */
/* |h[8]| <= 2^25; from now on fits into int32 unchanged */
/* |h[5]| <= 1.01*2^24 */
/* |h[9]| <= 1.51*2^58 */
c[9] = (h[9] + (1 << 24)) >> 25;
h[0] += c[9] * 19;
h[9] -= c[9] << 25;
/* |h[9]| <= 2^24; from now on fits into int32 unchanged */
/* |h[0]| <= 1.8*2^37 */
c[0] = (h[0] + (1 << 25)) >> 26;
h[1] += c[0];
h[0] -= c[0] << 26;
/* |h[0]| <= 2^25; from now on fits into int32 unchanged */
/* |h[1]| <= 1.01*2^24 */
let mut output = FieldElement32([0i32; 10]);
output.0[0] = h[0] as i32;
output.0[1] = h[1] as i32;
output.0[2] = h[2] as i32;
output.0[3] = h[3] as i32;
output.0[4] = h[4] as i32;
output.0[5] = h[5] as i32;
output.0[6] = h[6] as i32;
output.0[7] = h[7] as i32;
output.0[8] = h[8] as i32;
output.0[9] = h[9] as i32;
output
} }
/// Load a `FieldElement64` from the low 255 bits of a 256-bit /// Load a `FieldElement64` from the low 255 bits of a 256-bit
@ -334,11 +336,19 @@ impl FieldElement32 {
/// encoding of every field element should decode, re-encode to /// encoding of every field element should decode, re-encode to
/// the canonical encoding, and check that the input was /// the canonical encoding, and check that the input was
/// canonical. /// canonical.
///
/// XXX the above applies to the 64-bit implementation; check that
/// it applies here too.
pub fn from_bytes(data: &[u8; 32]) -> FieldElement32 { //FeFromBytes pub fn from_bytes(data: &[u8; 32]) -> FieldElement32 { //FeFromBytes
let mut h = [0i64;10]; #[inline]
fn load3(b: &[u8]) -> u64 {
(b[0] as u64) | ((b[1] as u64) << 8) | ((b[2] as u64) << 16)
}
#[inline]
fn load4(b: &[u8]) -> u64 {
(b[0] as u64) | ((b[1] as u64) << 8) | ((b[2] as u64) << 16) | ((b[3] as u64) << 24)
}
let mut h = [0u64;10];
const LOW_23_BITS: u64 = (1 << 23) - 1;
h[0] = load4(&data[ 0..]); h[0] = load4(&data[ 0..]);
h[1] = load3(&data[ 4..]) << 6; h[1] = load3(&data[ 4..]) << 6;
h[2] = load3(&data[ 7..]) << 5; h[2] = load3(&data[ 7..]) << 5;
@ -348,52 +358,33 @@ impl FieldElement32 {
h[6] = load3(&data[20..]) << 7; h[6] = load3(&data[20..]) << 7;
h[7] = load3(&data[23..]) << 5; h[7] = load3(&data[23..]) << 5;
h[8] = load3(&data[26..]) << 4; h[8] = load3(&data[26..]) << 4;
h[9] = (load3(&data[29..]) & 8388607) << 2; h[9] = (load3(&data[29..]) & LOW_23_BITS) << 2;
FieldElement32::reduce(h) FieldElement32::reduce(h)
} }
/// Serialize this `FieldElement64` to a 32-byte array. The /// Serialize this `FieldElement64` to a 32-byte array. The
/// encoding is canonical. /// encoding is canonical.
pub fn to_bytes(&self) -> [u8; 32] { //FeToBytes pub fn to_bytes(&self) -> [u8; 32] {
// 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: [i32; 10] = self.0;
let mut q:i32 = (19*h[9] + (1 << 24)) >> 25; let inp = &self.0;
q = (h[0] + q) >> 26; // Reduce the value represented by `in` to the range [0,2*p)
let mut h: [u32; 10] = FieldElement32::reduce([
// XXX this cast is annoying
inp[0] as u64, inp[1] as u64, inp[2] as u64, inp[3] as u64, inp[4] as u64,
inp[5] as u64, inp[6] as u64, inp[7] as u64, inp[8] as u64, inp[9] as u64,
]).0;
// Let h be the value to encode.
//
// Write h = pq + r with 0 <= r < p. We want to compute r = h mod p.
//
// Since h < 2*p, 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: u32 = (h[0] + 19) >> 26;
q = (h[1] + q) >> 25; q = (h[1] + q) >> 25;
q = (h[2] + q) >> 26; q = (h[2] + q) >> 26;
q = (h[3] + q) >> 25; q = (h[3] + q) >> 25;
@ -404,45 +395,40 @@ impl FieldElement32 {
q = (h[8] + q) >> 26; q = (h[8] + q) >> 26;
q = (h[9] + q) >> 25; q = (h[9] + q) >> 25;
// Goal: Output h-(2^255-19)q, which is between 0 and 2^255-20. debug_assert!( q == 0 || q == 1 );
h[0] += 19 * q;
// Goal: Output h-2^255 q, which is between 0 and 2^255-20.
carry[0] = h[0] >> 26; // Now we can compute r as r = h - pq = r - (2^255-19)q = r + 19q - 2^255q
h[1] += carry[0];
h[0] -= carry[0] << 26;
carry[1] = h[1] >> 25;
h[2] += carry[1];
h[1] -= carry[1] << 25;
carry[2] = h[2] >> 26;
h[3] += carry[2];
h[2] -= carry[2] << 26;
carry[3] = h[3] >> 25;
h[4] += carry[3];
h[3] -= carry[3] << 25;
carry[4] = h[4] >> 26;
h[5] += carry[4];
h[4] -= carry[4] << 26;
carry[5] = h[5] >> 25;
h[6] += carry[5];
h[5] -= carry[5] << 25;
carry[6] = h[6] >> 26;
h[7] += carry[6];
h[6] -= carry[6] << 26;
carry[7] = h[7] >> 25;
h[8] += carry[7];
h[7] -= carry[7] << 25;
carry[8] = h[8] >> 26;
h[9] += carry[8];
h[8] -= carry[8] << 26;
carry[9] = h[9] >> 25;
h[9] -= carry[9] << 25;
// h10 = carry9
// Goal: Output h[0]+...+2^255 h10-2^255 q, which is between 0 and 2^255-20. const LOW_25_BITS: u32 = (1 << 25) - 1;
// Have h[0]+...+2^230 h[9] between 0 and 2^255-1; const LOW_26_BITS: u32 = (1 << 26) - 1;
// evidently 2^255 h10-2^255 q = 0.
// Goal: Output h[0]+...+2^230 h[9]. h[0] += 19*q;
// Now carry the result to compute r + 19q...
h[1] += h[0] >> 26;
h[0] = h[0] & LOW_26_BITS;
h[2] += h[1] >> 25;
h[1] = h[1] & LOW_25_BITS;
h[3] += h[2] >> 26;
h[2] = h[2] & LOW_26_BITS;
h[4] += h[3] >> 25;
h[3] = h[3] & LOW_25_BITS;
h[5] += h[4] >> 26;
h[4] = h[4] & LOW_26_BITS;
h[6] += h[5] >> 25;
h[5] = h[5] & LOW_25_BITS;
h[7] += h[6] >> 26;
h[6] = h[6] & LOW_26_BITS;
h[8] += h[7] >> 25;
h[7] = h[7] & LOW_25_BITS;
h[9] += h[8] >> 26;
h[8] = h[8] & LOW_26_BITS;
// ... but instead of carrying the value
// (h[9] >> 25) = q*2^255 into another limb,
// discard it, subtracting the value from h.
debug_assert!( (h[9] >> 25) == 0 || (h[9] >> 25) == 1);
h[9] = h[9] & LOW_25_BITS;
let mut s = [0u8; 32]; let mut s = [0u8; 32];
s[0] = (h[0] >> 0) as u8; s[0] = (h[0] >> 0) as u8;
@ -484,77 +470,53 @@ impl FieldElement32 {
s s
} }
fn square_inner(&self) -> [i64; 10] { fn square_inner(&self) -> [u64; 10] {
let f0 = self.0[0] as i64; // Optimized version of multiplication for the case of squaring.
let f1 = self.0[1] as i64; // Pre- and post- conditions identical to multiplication function.
let f2 = self.0[2] as i64; let x = &self.0;
let f3 = self.0[3] as i64; let x0_2 = 2 * x[0];
let f4 = self.0[4] as i64; let x1_2 = 2 * x[1];
let f5 = self.0[5] as i64; let x2_2 = 2 * x[2];
let f6 = self.0[6] as i64; let x3_2 = 2 * x[3];
let f7 = self.0[7] as i64; let x4_2 = 2 * x[4];
let f8 = self.0[8] as i64; let x5_2 = 2 * x[5];
let f9 = self.0[9] as i64; let x6_2 = 2 * x[6];
let f0_2 = (2 * self.0[0]) as i64; let x7_2 = 2 * x[7];
let f1_2 = (2 * self.0[1]) as i64; let x5_19 = 19 * x[5];
let f2_2 = (2 * self.0[2]) as i64; let x6_19 = 19 * x[6];
let f3_2 = (2 * self.0[3]) as i64; let x7_19 = 19 * x[7];
let f4_2 = (2 * self.0[4]) as i64; let x8_19 = 19 * x[8];
let f5_2 = (2 * self.0[5]) as i64; let x9_19 = 19 * x[9];
let f6_2 = (2 * self.0[6]) as i64;
let f7_2 = (2 * self.0[7]) as i64;
let f5_38 = 38 * f5; // 1.31*2^30
let f6_19 = 19 * f6; // 1.31*2^30
let f7_38 = 38 * f7; // 1.31*2^30
let f8_19 = 19 * f8; // 1.31*2^30
let f9_38 = 38 * f9; // 1.31*2^30
let mut h = [0i64;10]; /// Helper function to multiply two 32-bit integers with 64 bits
h[0] = f0*f0 + f1_2*f9_38 + f2_2*f8_19 + f3_2*f7_38 + f4_2*f6_19 + f5*f5_38; /// of output.
h[1] = f0_2*f1 + f2*f9_38 + f3_2*f8_19 + f4*f7_38 + f5_2*f6_19; #[inline(always)]
h[2] = f0_2*f2 + f1_2*f1 + f3_2*f9_38 + f4_2*f8_19 + f5_2*f7_38 + f6*f6_19; fn m(x: u32, y: u32) -> u64 { (x as u64) * (y as u64) }
h[3] = f0_2*f3 + f1_2*f2 + f4*f9_38 + f5_2*f8_19 + f6*f7_38;
h[4] = f0_2*f4 + f1_2*f3_2 + f2*f2 + f5_2*f9_38 + f6_2*f8_19 + f7*f7_38;
h[5] = f0_2*f5 + f1_2*f4 + f2_2*f3 + f6*f9_38 + f7_2*f8_19;
h[6] = f0_2*f6 + f1_2*f5_2 + f2_2*f4 + f3_2*f3 + f7_2*f9_38 + f8*f8_19;
h[7] = f0_2*f7 + f1_2*f6 + f2_2*f5 + f3_2*f4 + f8*f9_38;
h[8] = f0_2*f8 + f1_2*f7_2 + f2_2*f6 + f3_2*f5_2 + f4*f4 + f9*f9_38;
h[9] = f0_2*f9 + f1_2*f8 + f2_2*f7 + f3_2*f6 + f4_2*f5;
h // This block is rearranged so that instead of doing a 32-bit multiplication by 38, we do a
// 64-bit multiplication by 2 on the results. This is because lg(38) is too big: we would
// have less than 1 bit of headroom left, which is too little.
let mut z = [0u64;10];
z[0] = m(x[0],x[0]) + m(x2_2,x8_19) + m(x4_2,x6_19) + (m(x1_2,x9_19) + m(x3_2,x7_19) + m(x[5],x5_19))*2;
z[1] = m(x0_2,x[1]) + m(x3_2,x8_19) + m(x5_2,x6_19) + (m(x[2],x9_19) + m(x[4],x7_19))*2;
z[2] = m(x0_2,x[2]) + m(x1_2,x[1]) + m(x4_2,x8_19) + m(x[6],x6_19) + (m(x3_2,x9_19) + m(x5_2,x7_19))*2;
z[3] = m(x0_2,x[3]) + m(x1_2,x[2]) + m(x5_2,x8_19) + (m(x[4],x9_19) + m(x[6],x7_19))*2;
z[4] = m(x0_2,x[4]) + m(x1_2,x3_2) + m(x[2],x[2]) + m(x6_2,x8_19) + (m(x5_2,x9_19) + m(x[7],x7_19))*2;
z[5] = m(x0_2,x[5]) + m(x1_2,x[4]) + m(x2_2,x[3]) + m(x7_2,x8_19) + m(x[6],x9_19)*2;
z[6] = m(x0_2,x[6]) + m(x1_2,x5_2) + m(x2_2,x[4]) + m(x3_2,x[3]) + m(x[8],x8_19) + m(x7_2,x9_19)*2;
z[7] = m(x0_2,x[7]) + m(x1_2,x[6]) + m(x2_2,x[5]) + m(x3_2,x[4]) + m(x[8],x9_19)*2;
z[8] = m(x0_2,x[8]) + m(x1_2,x7_2) + m(x2_2,x[6]) + m(x3_2,x5_2) + m(x[4],x[4]) + m(x[9],x9_19)*2;
z[9] = m(x0_2,x[9]) + m(x1_2,x[8]) + m(x2_2,x[7]) + m(x3_2,x[6]) + m(x4_2,x[5]) ;
z
} }
/// Calculates h = f*f. Can overlap h with f. /// Compute `self^2`.
///
/// 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.
///
/// # Postconditions
///
/// * |h[i]| bounded by 1.1*2^25, 1.1*2^24, 1.1*2^25, 1.1*2^24, etc.
pub fn square(&self) -> FieldElement32 { pub fn square(&self) -> FieldElement32 {
FieldElement32::reduce(self.square_inner()) FieldElement32::reduce(self.square_inner())
} }
/// Square this field element and multiply the result by 2. /// Compute `2*self^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.
///
/// # Postconditions
///
/// * |h[i]| bounded by 1.01*2^25, 1.01*2^24, 1.01*2^25, 1.01*2^24, etc.
///
/// # Notes
///
/// See fe_mul.c in ref10 implementation for discussion of implementation
/// strategy.
pub fn square2(&self) -> FieldElement32 { pub fn square2(&self) -> FieldElement32 {
let mut coeffs = self.square_inner(); let mut coeffs = self.square_inner();
for i in 0..self.0.len() { for i in 0..self.0.len() {

15
src/field_64bit.rs
View file

@ -25,8 +25,6 @@ use core::ops::Neg;
use subtle::ConditionallyAssignable; use subtle::ConditionallyAssignable;
use utils::load8;
/// In the 64-bit implementation, field elements are represented in /// In the 64-bit implementation, field elements are represented in
/// radix 2^51 as five `u64`s. /// radix 2^51 as five `u64`s.
pub type Limb = u64; pub type Limb = u64;
@ -55,7 +53,7 @@ pub struct FieldElement64(pub (crate) [u64; 5]);
impl Debug for FieldElement64 { impl Debug for FieldElement64 {
fn fmt(&self, f: &mut ::core::fmt::Formatter) -> ::core::fmt::Result { fn fmt(&self, f: &mut ::core::fmt::Formatter) -> ::core::fmt::Result {
write!(f, "FieldElement64: {:?}", &self.0[..]) write!(f, "FieldElement64({:?})", &self.0[..])
} }
} }
@ -247,6 +245,17 @@ impl FieldElement64 {
/// canonical. /// canonical.
/// ///
pub fn from_bytes(bytes: &[u8; 32]) -> FieldElement64 { pub fn from_bytes(bytes: &[u8; 32]) -> FieldElement64 {
let 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)
};
let low_51_bit_mask = (1u64 << 51) - 1; let low_51_bit_mask = (1u64 << 51) - 1;
FieldElement64( FieldElement64(
// load bits [ 0, 64), no shift // load bits [ 0, 64), no shift

4
src/lib.rs
View file

@ -80,10 +80,6 @@ pub mod montgomery;
pub mod ristretto; pub mod ristretto;
// Other miscelaneous utilities.
pub mod utils;
// Low-level curve and point constants, as well as pre-computed curve group elements. // Low-level curve and point constants, as well as pre-computed curve group elements.
pub mod constants; pub mod constants;

13
src/montgomery.rs
View file

@ -426,7 +426,6 @@ impl<'a, 'b> Mul<&'b MontgomeryPoint> for &'a Scalar {
#[cfg(test)] #[cfg(test)]
mod test { mod test {
use constants::ED25519_BASEPOINT_TABLE;
use constants::BASE_COMPRESSED_MONTGOMERY; use constants::BASE_COMPRESSED_MONTGOMERY;
use edwards::Identity; use edwards::Identity;
use super::*; use super::*;
@ -484,26 +483,29 @@ mod test {
} }
#[test] #[test]
#[cfg(feature="precomputed_tables")]
fn montgomery_ct_eq_ne() { fn montgomery_ct_eq_ne() {
let mut csprng: OsRng = OsRng::new().unwrap(); let mut csprng: OsRng = OsRng::new().unwrap();
let s1: Scalar = Scalar::random(&mut csprng); let s1: Scalar = Scalar::random(&mut csprng);
let s2: Scalar = Scalar::random(&mut csprng); let s2: Scalar = Scalar::random(&mut csprng);
let p1: MontgomeryPoint = (&s1 * &ED25519_BASEPOINT_TABLE).to_montgomery(); let p1: MontgomeryPoint = (&s1 * &constants::ED25519_BASEPOINT_TABLE).to_montgomery();
let p2: MontgomeryPoint = (&s2 * &ED25519_BASEPOINT_TABLE).to_montgomery(); let p2: MontgomeryPoint = (&s2 * &constants::ED25519_BASEPOINT_TABLE).to_montgomery();
assert_eq!(p1.ct_eq(&p2), 0); assert_eq!(p1.ct_eq(&p2), 0);
} }
#[test] #[test]
#[cfg(feature="precomputed_tables")]
fn montgomery_ct_eq_eq() { fn montgomery_ct_eq_eq() {
let mut csprng: OsRng = OsRng::new().unwrap(); let mut csprng: OsRng = OsRng::new().unwrap();
let s1: Scalar = Scalar::random(&mut csprng); let s1: Scalar = Scalar::random(&mut csprng);
let p1: MontgomeryPoint = (&s1 * &ED25519_BASEPOINT_TABLE).to_montgomery(); let p1: MontgomeryPoint = (&s1 * &constants::ED25519_BASEPOINT_TABLE).to_montgomery();
assert_eq!(p1.ct_eq(&p1), 1); assert_eq!(p1.ct_eq(&p1), 1);
} }
#[test] #[test]
#[cfg(feature="precomputed_tables")]
fn differential_add_matches_edwards_model() { fn differential_add_matches_edwards_model() {
let mut csprng: OsRng = OsRng::new().unwrap(); let mut csprng: OsRng = OsRng::new().unwrap();
@ -523,6 +525,7 @@ mod test {
} }
#[test] #[test]
#[cfg(feature="precomputed_tables")]
fn ladder_matches_scalarmult() { fn ladder_matches_scalarmult() {
let mut csprng: OsRng = OsRng::new().unwrap(); let mut csprng: OsRng = OsRng::new().unwrap();
@ -547,6 +550,7 @@ mod test {
#[test] #[test]
#[should_panic(expected = "assertion failed: self[31] <= 127")] #[should_panic(expected = "assertion failed: self[31] <= 127")]
#[cfg(feature="precomputed_tables")]
fn ladder_matches_scalarmult_with_scalar_high_bit_set() { fn ladder_matches_scalarmult_with_scalar_high_bit_set() {
let mut s: Scalar = Scalar::one(); let mut s: Scalar = Scalar::one();
@ -560,6 +564,7 @@ mod test {
} }
#[cfg(all(test, feature = "bench"))] #[cfg(all(test, feature = "bench"))]
#[cfg(feature="precomputed_tables")]
mod bench { mod bench {
use rand::OsRng; use rand::OsRng;
use constants::ED25519_BASEPOINT_TABLE; use constants::ED25519_BASEPOINT_TABLE;

4
src/ristretto.rs
View file

@ -1133,6 +1133,7 @@ mod test {
} }
#[test] #[test]
#[cfg(feature="precomputed_tables")]
fn four_torsion_random() { fn four_torsion_random() {
let mut rng = OsRng::new().unwrap(); let mut rng = OsRng::new().unwrap();
let B = &constants::RISTRETTO_BASEPOINT_TABLE; let B = &constants::RISTRETTO_BASEPOINT_TABLE;
@ -1195,6 +1196,7 @@ mod test {
} }
#[test] #[test]
#[cfg(feature="precomputed_tables")]
fn random_roundtrip() { fn random_roundtrip() {
let mut rng = OsRng::new().unwrap(); let mut rng = OsRng::new().unwrap();
let B = &constants::RISTRETTO_BASEPOINT_TABLE; let B = &constants::RISTRETTO_BASEPOINT_TABLE;
@ -1227,6 +1229,7 @@ mod bench {
use super::*; use super::*;
#[bench] #[bench]
#[cfg(feature="precomputed_tables")]
fn decompression(b: &mut Bencher) { fn decompression(b: &mut Bencher) {
let mut rng = OsRng::new().unwrap(); let mut rng = OsRng::new().unwrap();
let B = &constants::RISTRETTO_BASEPOINT_TABLE; let B = &constants::RISTRETTO_BASEPOINT_TABLE;
@ -1236,6 +1239,7 @@ mod bench {
} }
#[bench] #[bench]
#[cfg(feature="precomputed_tables")]
fn compression(b: &mut Bencher) { fn compression(b: &mut Bencher) {
let mut rng = OsRng::new().unwrap(); let mut rng = OsRng::new().unwrap();
let B = &constants::RISTRETTO_BASEPOINT_TABLE; let B = &constants::RISTRETTO_BASEPOINT_TABLE;

44
src/utils.rs
View file

@ -1,44 +0,0 @@
// -*- mode: rust; -*-
//
// This file is part of curve25519-dalek.
// Copyright (c) 2016-2017 Isis Lovecruft, Henry de Valence
// See LICENSE for licensing information.
//
// Authors:
// - Isis Agora Lovecruft <isis@patternsinthevoid.net>
// - Henry de Valence <hdevalence@hdevalence.ca>
//! Miscellaneous common utility functions.
/// Convert an array of (at least) three bytes into an i64.
#[inline]
//#[allow(dead_code)]
pub fn load3(input: &[u8]) -> i64 {
(input[0] as i64)
| ((input[1] as i64) << 8)
| ((input[2] as i64) << 16)
}
/// Convert an array of (at least) four bytes into an i64.
#[inline]
//#[allow(dead_code)]
pub fn load4(input: &[u8]) -> i64 {
(input[0] as i64)
| ((input[1] as i64) << 8)
| ((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)
}