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

Merge branch 'release/0.13.0'

Browse source
This commit is contained in:
Isis Lovecruft 2017-11-16 00:51:14 +00:00
commit 80ea4a480a
Failed to extract signature
20 changed files with 3264 additions and 1522 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

58
Cargo.toml
View file

@ -1,6 +1,6 @@
[package]
name = "curve25519-dalek"
version = "0.12.1"
version = "0.13.0"
authors = ["Isis Lovecruft <isis@patternsinthevoid.net>",
"Henry de Valence <hdevalence@hdevalence.ca>"]
readme = "README.md"
@ -15,6 +15,9 @@ exclude = [
".gitignore"
]
[package.metadata.docs.rs]
rustdoc-args = ["--html-in-header rustdoc-include-katex-header.html"]
[badges]
travis-ci = { repository = "isislovecruft/curve25519-dalek", branch = "master"}
@ -51,59 +54,8 @@ nightly = ["radix_51", "subtle/nightly"]
default = ["std"]
std = ["rand", "subtle/std"]
alloc = []
# This isn't used at the moment, but keep it around for future yolocrypto features.
yolocrypto = []
bench = []
# Radix-51 arithmetic using u128
radix_51 = []
# The development profile, used for `cargo build`.
[profile.dev]
opt-level = 0 # controls the `--opt-level` the compiler builds with
debug = true # controls whether the compiler passes `-g`
rpath = false # controls whether the compiler passes `-C rpath`
lto = false # controls `-C lto` for binaries and staticlibs
debug-assertions = true # controls whether debug assertions are enabled
codegen-units = 1 # controls whether the compiler passes `-C codegen-units`
# `codegen-units` is ignored when `lto = true`
panic = 'unwind' # panic strategy (`-C panic=...`), can also be 'abort'
# The release profile, used for `cargo build --release`.
[profile.release]
opt-level = 3
debug = false
rpath = false
lto = false
debug-assertions = false
codegen-units = 1
panic = 'unwind'
# The testing profile, used for `cargo test`.
[profile.test]
opt-level = 0
debug = true
rpath = false
lto = false
debug-assertions = true
codegen-units = 1
panic = 'unwind'
required-features = ['yolocrypto']
# The benchmarking profile, used for `cargo bench`.
[profile.bench]
opt-level = 3
debug = false
rpath = false
lto = false
debug-assertions = false
codegen-units = 1
panic = 'unwind'
# The documentation profile, used for `cargo doc`.
[profile.doc]
opt-level = 0
debug = true
rpath = false
lto = false
debug-assertions = true
codegen-units = 1
panic = 'unwind'

3
Makefile Normal file
View file

@ -0,0 +1,3 @@
doc:
cargo rustdoc --features "nightly yolocrypto" -- --html-in-header katex-header.html

10
README.md
View file

@ -44,7 +44,7 @@ Extensive documentation is available [here](https://docs.rs/curve25519-dalek).
To install, add the following to the dependencies section of your project's
`Cargo.toml`:
curve25519-dalek = "^0.12"
curve25519-dalek = "^0.13"
Then, in your library or executable source, add:
@ -60,3 +60,11 @@ We intend to stabilise the following before curve25519-dalek-1.0.0:
* Implement hashing to a point on the curve (Elligator).
* Finish Ristretto (Decaf for curve25519) implementation.
## Contributing
Please see
[CONTRIBUTING.md](https://github.com/isislovecruft/curve25519-dalek/blob/master/CONTRIBUTING.md).
Patches and pull requests should be make against the `develop`
branch, **not** `master`.

12
fuzz/Cargo.toml
View file

@ -10,7 +10,11 @@ cargo-fuzz = true
[dependencies.curve25519-dalek]
path = ".."
features = ["yolocrypto"]
[features]
yolocrypto = ["curve25519-dalek/yolocrypto"]
nightly = ["curve25519-dalek/nightly"]
radix_51 = ["curve25519-dalek/radix_51"]
[dependencies.libfuzzer-sys]
git = "https://github.com/rust-fuzz/libfuzzer-sys.git"
@ -21,4 +25,8 @@ members = ["."]
[[bin]]
name = "decaf"
path = "fuzzers/decaf.rs"
path = "fuzz_targets/decaf.rs"
[[bin]]
name = "scalar_constructor_accepts_256bit_values"
path = "fuzz_targets/scalar_constructor_accepts_256bit_values.rs"

0
fuzz/fuzzers/decaf.rs → fuzz/fuzz_targets/decaf.rs
View file

36
fuzz/fuzz_targets/scalar_constructor_accepts_256bit_values.rs Normal file
View file

@ -0,0 +1,36 @@
#![no_main]
#[macro_use] extern crate libfuzzer_sys;
extern crate curve25519_dalek;
use curve25519_dalek::scalar::Scalar;
/// Check that the Scalar constructor accepts 256-bit input values and
/// behaves correctly on them.
///
/// Specifically, we take 256-bit values `a` and `b` from the fuzzer
/// input data and check that `(a mod l) * (b mod l) == (a * b) mod l`.
fuzz_target!(|data: &[u8]| {
if data.len() != 64 {
return;
}
let mut a_bytes = [0u8; 32];
let mut b_bytes = [0u8; 32];
// Set a, b to be random 256-bit integers
a_bytes.copy_from_slice(&data[ 0..32]);
b_bytes.copy_from_slice(&data[32..64]);
// Compute c = a*b (mod l)
let c1 = &Scalar(a_bytes) * &Scalar(b_bytes);
// Compute c = (a mod l) * (b mod l)
let mut tmp = [0u8; 64];
tmp[0..32].copy_from_slice(&a_bytes[..]);
let a_mod_l = Scalar::reduce(&tmp);
tmp[0..32].copy_from_slice(&b_bytes[..]);
let b_mod_l = Scalar::reduce(&tmp);
let c2 = &a_mod_l * &b_mod_l;
assert_eq!(c1, c2);
});

10
rustdoc-include-katex-header.html Normal file
View file

@ -0,0 +1,10 @@
<link rel="stylesheet" href="https://cdnjs.cloudflare.com/ajax/libs/KaTeX/0.8.3/katex.min.css" integrity="sha384-B41nY7vEWuDrE9Mr+J2nBL0Liu+nl/rBXTdpQal730oTHdlrlXHzYMOhDU60cwde" crossorigin="anonymous">
<script src="https://cdnjs.cloudflare.com/ajax/libs/KaTeX/0.8.3/katex.min.js" integrity="sha384-L9gv4ooDLrYwW0QCM6zY3EKSSPrsuUncpx26+erN0pJX4wv1B1FzVW1SvpcJPx/8" crossorigin="anonymous"></script>
<script src="https://cdnjs.cloudflare.com/ajax/libs/KaTeX/0.8.3/contrib/auto-render.min.js" integrity="sha384-RkgGHBDdR8eyBOoWeZ/vpGg1cOvSAJRflCUDACusAAIVwkwPrOUYykglPeqWakZu" crossorigin="anonymous"></script>
<script>
document.addEventListener("DOMContentLoaded", function() { renderMathInElement(document.body); });
</script>
<style>
.katex { font-size: 1em !important; }
.docblock code, .dockblock-short code { font-size: 0.85em !important; }
</style>

131
src/constants.rs
View file

@ -11,16 +11,26 @@
//! This module contains various constants (such as curve parameters
//! and useful field elements like `sqrt(-1)`), as well as
//! lookup tables of pre-computed points.
//!
//! Most of the constants are given with
//! `LONG_DESCRIPTIVE_UPPER_CASE_NAMES`, but they can be brought into
//! scope using a `let` binding:
//!
//! ```
//! use curve25519_dalek::constants;
//! use curve25519_dalek::edwards::IsIdentity;
//!
//! let B = &constants::RISTRETTO_BASEPOINT_TABLE;
//! let l = &constants::BASEPOINT_ORDER;
//!
//! let A = l * B;
//! assert!(A.is_identity());
//! ```
#![allow(dead_code)]
#![allow(non_snake_case)]
#![allow(non_upper_case_globals)]
#![allow(missing_docs)]
#![allow(non_snake_case)]
use edwards::CompressedEdwardsY;
#[cfg(feature = "yolocrypto")]
use decaf::{DecafPoint, DecafBasepointTable};
use ristretto::{RistrettoPoint, RistrettoBasepointTable};
use montgomery::CompressedMontgomeryU;
use scalar::Scalar;
@ -29,20 +39,6 @@ pub use constants_64bit::*;
#[cfg(not(feature="radix_51"))]
pub use constants_32bit::*;
/// (p-1)/2, in little-endian bytes.
pub const HALF_P_MINUS_1_BYTES: [u8; 32] =
[0xf6, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff,
0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff,
0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff,
0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0x3f];
/// `HALF_Q_MINUS_1_BYTES` is (2^255-20)/2 expressed in little endian form.
pub const HALF_Q_MINUS_1_BYTES: [u8; 32] = [ // halfQMinus1Bytes
0xf6, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff,
0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff,
0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff,
0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0x3f, ];
/// Basepoint has y = 4/5.
///
/// Generated with Sage: these are the bytes of 4/5 in 𝔽_p. The
@ -61,36 +57,38 @@ pub const BASE_COMPRESSED_MONTGOMERY: CompressedMontgomeryU =
0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00]);
/// The Ed25519 basepoint, as a `DecafPoint`. This is called `_POINT` to distinguish it from
/// The Ed25519 basepoint, as a `RistrettoPoint`. This is called `_POINT` to distinguish it from
/// `_TABLE`, which provides fast scalar multiplication.
#[cfg(feature = "yolocrypto")] pub const DECAF_ED25519_BASEPOINT_POINT: DecafPoint =
DecafPoint(ED25519_BASEPOINT_POINT);
pub const RISTRETTO_BASEPOINT_POINT: RistrettoPoint = RistrettoPoint(ED25519_BASEPOINT_POINT);
/// `l` is the order of base point, i.e. 2^252 +
/// 27742317777372353535851937790883648493, in little-endian form
pub const l: Scalar = Scalar([ 0xed, 0xd3, 0xf5, 0x5c, 0x1a, 0x63, 0x12, 0x58,
0xd6, 0x9c, 0xf7, 0xa2, 0xde, 0xf9, 0xde, 0x14,
0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00,
0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x10 ]);
/// `BASEPOINT_ORDER` is the order of base point, i.e. `l = 2^252 +
/// 27742317777372353535851937790883648493`, in little-endian bytes.
pub const BASEPOINT_ORDER: Scalar = Scalar([
0xed, 0xd3, 0xf5, 0x5c, 0x1a, 0x63, 0x12, 0x58,
0xd6, 0x9c, 0xf7, 0xa2, 0xde, 0xf9, 0xde, 0x14,
0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00,
0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x10,
]);
/// `l_minus_1` is the order of base point minus one, i.e. 2^252 +
/// 27742317777372353535851937790883648493 - 1, in little-endian form
pub const l_minus_1: Scalar = Scalar([ 0xec, 0xd3, 0xf5, 0x5c, 0x1a, 0x63, 0x12, 0x58,
0xd6, 0x9c, 0xf7, 0xa2, 0xde, 0xf9, 0xde, 0x14,
0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00,
0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x10 ]);
/// `BASEPOINT_ORDER_MINUS_1` is the order of base point minus one, i.e. `l-1`, in little-endian bytes.
pub const BASEPOINT_ORDER_MINUS_1: Scalar = Scalar([
0xec, 0xd3, 0xf5, 0x5c, 0x1a, 0x63, 0x12, 0x58,
0xd6, 0x9c, 0xf7, 0xa2, 0xde, 0xf9, 0xde, 0x14,
0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00,
0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x10,
]);
/// `lminus1` is the order of base point minus two, i.e. 2^252 +
/// 27742317777372353535851937790883648493 - 2, in little-endian form
pub const l_minus_2: Scalar = Scalar([ 0xeb, 0xd3, 0xf5, 0x5c, 0x1a, 0x63, 0x12, 0x58,
0xd6, 0x9c, 0xf7, 0xa2, 0xde, 0xf9, 0xde, 0x14,
0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00,
0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x10 ]);
/// `BASEPOINT_ORDER_MINUS_2` is the order of base point minus two, i.e. `l-2`, in little-endian bytes.
pub const BASEPOINT_ORDER_MINUS_2: Scalar = Scalar([
0xeb, 0xd3, 0xf5, 0x5c, 0x1a, 0x63, 0x12, 0x58,
0xd6, 0x9c, 0xf7, 0xa2, 0xde, 0xf9, 0xde, 0x14,
0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00,
0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x10,
]);
#[cfg(feature = "yolocrypto")]
/// The Ed25519 basepoint
pub const DECAF_ED25519_BASEPOINT_TABLE: DecafBasepointTable
= DecafBasepointTable(ED25519_BASEPOINT_TABLE);
/// The Ed25519 basepoint, as a RistrettoPoint
pub const RISTRETTO_BASEPOINT_TABLE: RistrettoBasepointTable
= RistrettoBasepointTable(ED25519_BASEPOINT_TABLE);
#[cfg(test)]
mod test {
@ -126,13 +124,6 @@ mod test {
}
}
#[test]
fn test_half() {
let one = FieldElement::one();
let two = &one + &one;
assert_eq!(one, &two * &constants::HALF);
}
/// Test that the constant for sqrt(-486664) really is a square
/// root of -486664.
#[test]
@ -158,27 +149,22 @@ mod test {
}
#[test]
/// Test that SQRT_M1 and MSQRT_M1 are square roots of -1
/// Test that SQRT_M1 is a square root of -1
fn test_sqrt_minus_one() {
let minus_one = FieldElement::minus_one();
let sqrt_m1_sq = &constants::SQRT_M1 * &constants::SQRT_M1;
let msqrt_m1_sq = &constants::MSQRT_M1 * &constants::MSQRT_M1;
assert_eq!(minus_one, sqrt_m1_sq);
assert_eq!(minus_one, msqrt_m1_sq);
}
#[test]
fn test_sqrt_constants_sign() {
let one = FieldElement::one();
let minus_one = FieldElement::minus_one();
let (was_nonzero_square, invsqrt_m1) = minus_one.invsqrt();
assert_eq!(was_nonzero_square, 1u8);
let sign_test_sqrt = &invsqrt_m1 * &constants::SQRT_M1;
let sign_test_msqrt = &invsqrt_m1 * &constants::MSQRT_M1;
// XXX it seems we have flipped the sign relative to
// the invsqrt function?
assert_eq!(sign_test_sqrt, minus_one);
assert_eq!(sign_test_msqrt, one);
}
/// Test that d = -121665/121666
@ -190,8 +176,8 @@ mod test {
let b = FieldElement32([ 121666,0,0,0,0,0,0,0,0,0]);
let d = &a * &b.invert();
let d2 = &d + &d;
assert_eq!(d, constants::d);
assert_eq!(d2, constants::d2);
assert_eq!(d, constants::EDWARDS_D);
assert_eq!(d2, constants::EDWARDS_D2);
}
/// Test that d = -121665/121666
@ -203,27 +189,16 @@ mod test {
let b = FieldElement64([121666,0,0,0,0]);
let d = &a * &b.invert();
let d2 = &d + &d;
assert_eq!(d, constants::d);
assert_eq!(d2, constants::d2);
assert_eq!(d, constants::EDWARDS_D);
assert_eq!(d2, constants::EDWARDS_D2);
}
#[test]
fn test_d4() {
let mut four = FieldElement::zero();
// XXX should have a way to create small field elements
four.0[0] = 4;
assert_eq!(&constants::d * &four, constants::d4);
}
#[test]
fn test_a_minus_d() {
fn test_sqrt_ad_minus_one() {
let a = FieldElement::minus_one();
let a_minus_d = &a - &constants::d;
assert_eq!(a_minus_d, constants::a_minus_d);
let (_, invsqrt_a_minus_d) = constants::a_minus_d.invsqrt();
assert_eq!(invsqrt_a_minus_d, constants::invsqrt_a_minus_d);
let inv_a_minus_d = invsqrt_a_minus_d.square();
assert_eq!(inv_a_minus_d, constants::inv_a_minus_d);
assert_eq!(&inv_a_minus_d * &a_minus_d, FieldElement::one());
let ad_minus_one = &(&a * &constants::EDWARDS_D) + &a;
let should_be_ad_minus_one = constants::SQRT_AD_MINUS_ONE.square();
assert_eq!(should_be_ad_minus_one, ad_minus_one);
}
}

90
src/constants_32bit.rs
View file

@ -12,82 +12,49 @@
//! and useful field elements like `sqrt(-1)`), as well as
//! lookup tables of pre-computed points.
#![allow(dead_code)]
#![allow(non_snake_case)]
#![allow(non_upper_case_globals)]
#![allow(missing_docs)]
#![allow(non_snake_case)]
use field_32bit::FieldElement32;
use scalar_32bit::Scalar32;
use edwards::ExtendedPoint;
use edwards::AffineNielsPoint;
use edwards::EdwardsBasepointTable;
pub const d: FieldElement32 = FieldElement32([
/// Edwards `d` value, equal to `-121665/121666 mod p`.
pub(crate) const EDWARDS_D: FieldElement32 = FieldElement32([
-10913610, 13857413, -15372611, 6949391, 114729,
-8787816, -6275908, -3247719, -18696448, -12055116, ]);
pub const d2: FieldElement32 = FieldElement32([
/// Edwards `2*d` value, equal to `2*(-121665/121666) mod p`.
pub(crate) const EDWARDS_D2: FieldElement32 = FieldElement32([
-21827239, -5839606, -30745221, 13898782, 229458,
15978800, -12551817, -6495438, 29715968, 9444199, ]);
pub const d4: FieldElement32 = FieldElement32([
23454405, -11679213, 5618422, -5756869, 458917,
-1596832, -25103633, -12990876, -7676928, -14666033 ]);
/// `= sqrt(a*d - 1)`, where `a = -1 (mod p)`, `d` are the Edwards curve parameters.
pub(crate) const SQRT_AD_MINUS_ONE: FieldElement32 = FieldElement32([
24849947, -153582, -23613485, 6347715, -21072328, -667138, -25271143, -15367704, -870347, 14525639
]);
pub const a_minus_d: FieldElement32 = FieldElement32([
10913609, -13857413, 15372611, -6949391, -114729,
8787816, 6275908, 3247719, 18696448, 12055116, ]);
pub const invsqrt_a_minus_d: FieldElement32 = FieldElement32([
/// `= 1/sqrt(a-d)`, where `a = -1 (mod p)`, `d` are the Edwards curve parameters.
pub(crate) const INVSQRT_A_MINUS_D: FieldElement32 = FieldElement32([
6111485, 4156064, -27798727, 12243468, -25904040,
120897, 20826367, -7060776, 6093568, -1986012
]);
#[cfg(not(feature="radix_51"))]
pub const inv_a_minus_d: FieldElement32 = FieldElement32([
-121666, 0, 0, 0, 0, 0, 0, 0, 0, 0
]);
/// (p-1)/2, in little-endian bytes.
pub const HALF_P_MINUS_1_BYTES: [u8; 32] =
[0xf6, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff,
0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff,
0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff,
0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0x3f];
/// Precomputed value of one of the square roots of -1 (mod p)
pub const SQRT_M1: FieldElement32 = FieldElement32([
pub(crate) const SQRT_M1: FieldElement32 = FieldElement32([
-32595792, -7943725, 9377950, 3500415, 12389472,
-272473, -25146209, -2005654, 326686, 11406482, ]);
/// Precomputed value of the other square root of -1 (mod p),
/// i.e., `MSQRT_M1 = -SQRT_M1`.
pub const MSQRT_M1: FieldElement32 = FieldElement32([
32595792, 7943725, -9377950, -3500415, -12389472,
272473, 25146209, 2005654, -326686, -11406482, ]);
/// Precomputed value of 1/2 (mod p).
pub const HALF: FieldElement32 = FieldElement32([
10, 0, 0, 0, 0, 0, 0, 0, 0, -16777216, ]);
/// In Montgomery form y² = x³+Ax²+x, Curve25519 has A=486662.
pub const A: FieldElement32 = FieldElement32([
pub(crate) const MONTGOMERY_A: FieldElement32 = FieldElement32([
486662, 0, 0, 0, 0, 0, 0, 0, 0, 0, ]);
/// `APLUS2_OVER_FOUR` is (A+2)/4. (This is used internally within the Montgomery ladder.)
pub const APLUS2_OVER_FOUR: FieldElement32 = FieldElement32([121666, 0, 0, 0, 0, 0, 0, 0, 0, 0]);
/// `SQRT_MINUS_A` is sqrt(-486662)
// XXX I think that this was used in Adam's code for his elligator
// implementation, but that should maybe be using sqrt(-486664)
// instead...? - hdevalence
pub const SQRT_MINUS_A: FieldElement32 = FieldElement32([ // sqrtMinusA
12222970, 8312128, 11511410, -9067497, 15300785,
241793, -25456130, -14121551, 12187136, -3972024, ]);
pub(crate) const APLUS2_OVER_FOUR: FieldElement32 = FieldElement32([121666, 0, 0, 0, 0, 0, 0, 0, 0, 0]);
/// `SQRT_MINUS_APLUS2` is sqrt(-486664)
pub const SQRT_MINUS_APLUS2: FieldElement32 = FieldElement32([
pub(crate) const SQRT_MINUS_APLUS2: FieldElement32 = FieldElement32([
-12222970, -8312128, -11511410, 9067497, -15300785,
-241793, 25456130, 14121551, -12187136, 3972024]);
@ -96,8 +63,28 @@ pub const SQRT_MINUS_HALF: FieldElement32 = FieldElement32([ // sqrtMinusHalf
-17256545, 3971863, 28865457, -1750208, 27359696,
-16640980, 12573105, 1002827, -163343, 11073975, ]);
/// 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).
/// `L` is the order of base point, i.e. 2^252 +
/// 27742317777372353535851937790883648493
pub(crate) const L: Scalar32 = Scalar32([ 0x1cf5d3ed, 0x009318d2, 0x1de73596, 0x1df3bd45,
0x0000014d, 0x00000000, 0x00000000, 0x00000000,
0x00100000 ]);
/// `L` * `LFACTOR` = -1 (mod 2^29)
pub(crate) const LFACTOR: u32 = 0x12547e1b;
/// `R` = R % L where R = 2^261
pub(crate) const R: Scalar32 = Scalar32([ 0x114df9ed, 0x1a617303, 0x0f7c098c, 0x16793167,
0x1ffd656e, 0x1fffffff, 0x1fffffff, 0x1fffffff,
0x000fffff ]);
/// `RR` = (R^2) % L where R = 2^261
pub(crate) const RR: Scalar32 = Scalar32([ 0x0b5f9d12, 0x1e141b17, 0x158d7f3d, 0x143f3757,
0x1972d781, 0x042feb7c, 0x1ceec73d, 0x1e184d1e,
0x0005046d ]);
/// The Ed25519 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{
X: FieldElement32([-14297830, -7645148, 16144683, -16471763, 27570974, -2696100, -26142465, 8378389, 20764389, 8758491]),
Y: FieldElement32([-26843541, -6710886, 13421773, -13421773, 26843546, 6710886, -13421773, 13421773, -26843546, -6710886]),
@ -163,7 +150,8 @@ pub const EIGHT_TORSION: [ExtendedPoint; 8] = [
},
];
pub const bi: [AffineNielsPoint; 8] = [
/// Odd multiples of the basepoint `[B, 3B, 5B, 7B, 9B, 11B, 13B, 15B]`.
pub(crate) const AFFINE_ODD_MULTIPLES_OF_BASEPOINT: [AffineNielsPoint; 8] = [
AffineNielsPoint{
y_plus_x: FieldElement32([25967493, -14356035, 29566456, 3660896, -12694345, 4014787, 27544626, -11754271, -6079156, 2047605]),
y_minus_x: FieldElement32([-12545711, 934262, -2722910, 3049990, -727428, 9406986, 12720692, 5043384, 19500929, -15469378]),

69
src/constants_64bit.rs
View file

@ -12,65 +12,57 @@
//! and useful field elements like `sqrt(-1)`), as well as
//! lookup tables of pre-computed points.
#![allow(dead_code)]
#![allow(non_snake_case)]
#![allow(non_upper_case_globals)]
#![allow(missing_docs)]
#![allow(non_snake_case)]
use field_64bit::FieldElement64;
use scalar_64bit::Scalar64;
use edwards::ExtendedPoint;
use edwards::AffineNielsPoint;
use edwards::EdwardsBasepointTable;
pub const p: FieldElement64 = FieldElement64([2251799813685229, 2251799813685247, 2251799813685247, 2251799813685247, 2251799813685247]);
/// Edwards `d` value, equal to `-121665/121666 mod p`.
pub(crate) const EDWARDS_D: FieldElement64 = FieldElement64([929955233495203, 466365720129213, 1662059464998953, 2033849074728123, 1442794654840575]);
pub const d: FieldElement64 = FieldElement64([929955233495203, 466365720129213, 1662059464998953, 2033849074728123, 1442794654840575]);
/// Edwards `2*d` value, equal to `2*(-121665/121666) mod p`.
pub(crate) const EDWARDS_D2: FieldElement64 = FieldElement64([1859910466990425, 932731440258426, 1072319116312658, 1815898335770999, 633789495995903]);
pub const d2: FieldElement64 = FieldElement64([1859910466990425, 932731440258426, 1072319116312658, 1815898335770999, 633789495995903]);
/// `= sqrt(a*d - 1)`, where `a = -1 (mod p)`, `d` are the Edwards curve parameters.
pub(crate) const SQRT_AD_MINUS_ONE: FieldElement64 = FieldElement64([
2241493124984347, 425987919032274, 2207028919301688, 1220490630685848, 974799131293748
]);
pub const d4: FieldElement64 = FieldElement64([1468021120295602, 1865462880516853, 2144638232625316, 1379996857856750, 1267578991991807]);
pub const a_minus_d: FieldElement64 = FieldElement64([1321844580190025, 1785434093556034, 589740348686294, 217950738957124, 809005158844672]);
pub const invsqrt_a_minus_d: FieldElement64 = FieldElement64([
/// `= 1/sqrt(a-d)`, where `a = -1 (mod p)`, `d` are the Edwards curve parameters.
pub(crate) const INVSQRT_A_MINUS_D: FieldElement64 = FieldElement64([
278908739862762, 821645201101625, 8113234426968, 1777959178193151, 2118520810568447
]);
pub const inv_a_minus_d: FieldElement64 = FieldElement64([
2251799813563563, 2251799813685247, 2251799813685247, 2251799813685247, 2251799813685247
]);
/// Precomputed value of one of the square roots of -1 (mod p)
pub const SQRT_M1: FieldElement64 = FieldElement64([1718705420411056, 234908883556509, 2233514472574048, 2117202627021982, 765476049583133]);
/// Precomputed value of the other square root of -1 (mod p),
/// i.e., `MSQRT_M1 = -SQRT_M1`.
pub const MSQRT_M1: FieldElement64 = FieldElement64([533094393274173, 2016890930128738, 18285341111199, 134597186663265, 1486323764102114]);
/// Precomputed value of 1/2 (mod p).
pub const HALF: FieldElement64 = FieldElement64([2251799813685239, 2251799813685247, 2251799813685247, 2251799813685247, 1125899906842623]);
pub(crate) const SQRT_M1: FieldElement64 = FieldElement64([1718705420411056, 234908883556509, 2233514472574048, 2117202627021982, 765476049583133]);
/// In Montgomery form y² = x³+Ax²+x, Curve25519 has A=486662.
pub const A: FieldElement64 = FieldElement64([486662, 0, 0, 0, 0]);
pub(crate) const MONTGOMERY_A: FieldElement64 = FieldElement64([486662, 0, 0, 0, 0]);
/// `APLUS2_OVER_FOUR` is (A+2)/4. (This is used internally within the Montgomery ladder.)
pub const APLUS2_OVER_FOUR: FieldElement64 = FieldElement64([121666, 0, 0, 0, 0]);
/// `SQRT_MINUS_A` is sqrt(-486662)
// XXX I think that this was used in Adam's code for his elligator
// implementation, but that should maybe be using sqrt(-486664)
// instead...? - hdevalence
pub const SQRT_MINUS_A: FieldElement64 = FieldElement64([557817479725543, 1643290402203250, 16226468853936, 1304118542701054, 1985241807451647]);
pub(crate) const APLUS2_OVER_FOUR: FieldElement64 = FieldElement64([121666, 0, 0, 0, 0]);
/// `SQRT_MINUS_APLUS2` is sqrt(-486664)
pub const SQRT_MINUS_APLUS2: FieldElement64 = FieldElement64([1693982333959686, 608509411481997, 2235573344831311, 947681270984193, 266558006233600]);
pub(crate) const SQRT_MINUS_APLUS2: FieldElement64 = FieldElement64([1693982333959686, 608509411481997, 2235573344831311, 947681270984193, 266558006233600]);
/// `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(crate) const L: Scalar64 = Scalar64([ 0x0002631a5cf5d3ed, 0x000dea2f79cd6581, 0x000000000014def9, 0x0000000000000000, 0x0000100000000000 ]);
/// 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).
/// `L` * `LFACTOR` = -1 (mod 2^51)
pub(crate) const LFACTOR: u64 = 0x51da312547e1b;
/// `R` = R % L where R = 2^260
pub(crate) const R: Scalar64 = Scalar64([ 0x000f48bd6721e6ed, 0x0003bab5ac67e45a, 0x000fffffeb35e51b, 0x000fffffffffffff, 0x00000fffffffffff ]);
/// `RR` = (R^2) % L where R = 2^260
pub(crate) const RR: Scalar64 = Scalar64([ 0x0009d265e952d13b, 0x000d63c715bea69f, 0x0005be65cb687604, 0x0003dceec73d217f, 0x000009411b7c309a ]);
/// The Ed25519 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{
X: FieldElement64([1738742601995546, 1146398526822698, 2070867633025821, 562264141797630, 587772402128613]),
Y: FieldElement64([1801439850948184, 1351079888211148, 450359962737049, 900719925474099, 1801439850948198]),
@ -143,7 +135,8 @@ pub const EIGHT_TORSION: [ExtendedPoint; 8] = [
}
];
pub const bi: [AffineNielsPoint; 8] = [
/// Odd multiples of the basepoint `[B, 3B, 5B, 7B, 9B, 11B, 13B, 15B]`.
pub(crate) const AFFINE_ODD_MULTIPLES_OF_BASEPOINT: [AffineNielsPoint; 8] = [
AffineNielsPoint {
y_plus_x: FieldElement64([1288382639258501, 245678601348599, 269427782077623, 1462984067271730, 137412439391563]),
y_minus_x: FieldElement64([62697248952638, 204681361388450, 631292143396476, 338455783676468, 1213667448819585]),

871
src/decaf.rs
View file

@ -1,871 +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>
//! An implementation of Mike Hamburg's Decaf cofactor-eliminating
//! point-compression scheme, providing a prime-order group on top of
//! a non-prime-order elliptic curve.
//!
//! Note: this code is currently feature-gated with the `yolocrypto`
//! feature flag, because our implementation is still unfinished.
// We allow non snake_case names because coordinates in projective space are
// traditionally denoted by the capitalisation of their respective
// counterparts in affine space. Yeah, you heard me, rustc, I'm gonna have my
// affine and projective cakes and eat both of them too.
#![allow(non_snake_case)]
use core::fmt::Debug;
#[cfg(feature = "std")]
use rand::Rng;
use digest::Digest;
use generic_array::typenum::U32;
use constants;
use field::FieldElement;
use core::ops::{Add, Sub, Neg};
use core::ops::{AddAssign, SubAssign};
use core::ops::{Mul, MulAssign};
use edwards;
use edwards::ExtendedPoint;
use edwards::CompletedPoint;
use edwards::EdwardsBasepointTable;
use edwards::Identity;
use scalar::Scalar;
use subtle::ConditionallyAssignable;
use subtle::ConditionallyNegatable;
// ------------------------------------------------------------------------
// Compressed points
// ------------------------------------------------------------------------
/// A point serialized using Mike Hamburg's Decaf scheme.
///
/// XXX think about how this API should work
#[derive(Copy, Clone, Eq, PartialEq)]
pub struct CompressedDecaf(pub [u8; 32]);
/// The result of compressing a `DecafPoint`.
impl CompressedDecaf {
/// Convert this `CompressedDecaf` to an array of bytes.
pub fn to_bytes(&self) -> [u8; 32] {
self.0
}
/// View this `CompressedDecaf` as an array of bytes.
pub fn as_bytes<'a>(&'a self) -> &'a [u8; 32] {
&self.0
}
/// Attempt to decompress to an `DecafPoint`.
pub fn decompress(&self) -> Option<DecafPoint> {
// XXX should decoding be CT ?
// XXX need to check that xy is nonnegative and reject otherwise
let s = FieldElement::from_bytes(self.as_bytes());
// Check that s = |s| and reject otherwise.
let mut abs_s = s;
let neg = abs_s.is_negative_decaf();
abs_s.conditional_negate(neg);
if abs_s != s { return None; }
let ss = s.square();
let X = &s + &s; // X = 2s
let Z = &FieldElement::one() - &ss; // Z = 1+as^2
let ZZ = Z.square();
let u = &ZZ- &(&constants::d4 * &ss); // u = Z^2 - 4ds^2
let uss = &u * &ss;
let ussZZ = &uss * &ZZ;
if Z.is_zero() == 1u8 { return None; }
// Batch inversion: set b = 1/sqrt(us^2 Z^2)
let (ussZZ_is_nonzero_square, b) = ussZZ.invsqrt();
if (ussZZ_is_nonzero_square | uss.is_zero()) == 0u8 {
return None; // us^2 is nonzero nonsquare
}
let mut v = &b * &Z; // now v = 1/sqrt(us^2)
let Zinv = &b * &(&v * &uss); // now Zinv = b^2 Z us^2 = 1/Z
// Now v = 1/sqrt(us^2) if us^2 is a nonzero square, 0 if us^2 is zero.
let uv = &v * &u;
if uv.is_negative_decaf() == 1u8 {
v.negate();
}
let mut two_minus_Z = -&Z; two_minus_Z.0[0] += 2;
let mut w = &v * &(&s * &two_minus_Z);
w.conditional_assign(&FieldElement::one(), s.is_zero());
let Y = &w * &Z;
let T = &w * &X;
// "To decode the point, one must decode it to affine form
// instead of projective, and check that xy is non-negative."
// Use the value of 1/Z previously computed in the batch inversion
let xy = &T * &Zinv;
if (Y.is_nonzero() & xy.is_nonnegative_decaf()) == 1u8 {
Some(DecafPoint(ExtendedPoint{ X: X, Y: Y, Z: Z, T: T }))
} else {
None
}
}
}
impl Identity for CompressedDecaf {
fn identity() -> CompressedDecaf {
CompressedDecaf([0u8; 32])
}
}
// ------------------------------------------------------------------------
// Serde support
// ------------------------------------------------------------------------
// Serializes to and from `DecafPoint` directly, doing compression
// and decompression internally. This means that users can create
// structs containing `DecafPoint`s and use Serde's derived
// serializers to serialize those structures.
#[cfg(feature = "serde")]
use serde::{self, Serialize, Deserialize, Serializer, Deserializer};
#[cfg(feature = "serde")]
use serde::de::Visitor;
#[cfg(feature = "serde")]
impl Serialize for DecafPoint {
fn serialize<S>(&self, serializer: S) -> Result<S::Ok, S::Error>
where S: Serializer
{
serializer.serialize_bytes(self.compress().as_bytes())
}
}
#[cfg(feature = "serde")]
impl<'de> Deserialize<'de> for DecafPoint {
fn deserialize<D>(deserializer: D) -> Result<Self, D::Error>
where D: Deserializer<'de>
{
struct DecafPointVisitor;
impl<'de> Visitor<'de> for DecafPointVisitor {
type Value = DecafPoint;
fn expecting(&self, formatter: &mut ::core::fmt::Formatter) -> ::core::fmt::Result {
formatter.write_str("a valid point in Decaf format")
}
fn visit_bytes<E>(self, v: &[u8]) -> Result<DecafPoint, E>
where E: serde::de::Error
{
if v.len() == 32 {
let arr32 = array_ref!(v, 0, 32); // &[u8;32] from &[u8]
CompressedDecaf(*arr32)
.decompress()
.ok_or(serde::de::Error::custom("decompression failed"))
} else {
Err(serde::de::Error::invalid_length(v.len(), &self))
}
}
}
deserializer.deserialize_bytes(DecafPointVisitor)
}
}
// ------------------------------------------------------------------------
// Internal point representations
// ------------------------------------------------------------------------
/// A point in a prime-order group.
///
// XXX think about how this API should work
#[derive(Copy, Clone)]
pub struct DecafPoint(pub ExtendedPoint);
impl DecafPoint {
/// Compress in Decaf format.
pub fn compress(&self) -> CompressedDecaf {
// Q: Do we want to encode twisted or untwisted?
//
// Notes:
// Recall that the twisted Edwards curve E_{a,d} is of the form
//
// ax^2 + y^2 = 1 + dx^2y^2.
//
// Internally, we operate on the curve with a = -1, d =
// -121665/121666, a.k.a., the twist. But maybe we would like
// to use Decaf on the untwisted curve with a = 1, d =
// 121665/121666. (why? interop?)
//
// Fix i, a square root of -1 (mod p).
//
// The map x -> ix is an isomorphism from E_{a,d} to E_{-a,-d}.
// Its inverse is x -> -ix.
// let untwisted_X = &self.X * &constants::MSQRT_M1;
// etc.
//
// Step 0: pre-rotation, needed for Decaf with E[8] = Z/8.
//
// We want to select a point (x,y) in the coset P + E[4] with
// y nonzero and xy nonnegative. The naive approach is as
// follows. First, compute xy = T/Z and check that Y is
// nonzero and xy is nonnegative. If not, then "rotate" the
// original point by adding (x,y) + (i,0) = (iy,ix) = (x',y'):
// this rotated point has x'y' = -xy. Then perform the normal
// Decaf encoding, as described in Appendix A.1 of the Decaf
// paper, using the rotated point (x',y').
//
// This is straightforward but requires an extra inversion.
// We would like to batch the inversion in xy = T/Z with the
// inverse square root in the computation of
//
// r = invsqrt((a-d)*(Z+Y)*(Z-Y))
// = invsqrt(a-d)*invsqrt(Z^2-Y^2),
//
// but the X and Y we are trying to decode depend on whether
// we rotated the coset representative!
//
// However, it is possible to batch these inversions. Credit:
// the following explanation (and trick) is adapted from an
// email from Mike Hamburg, but of course any errors are ours.
//
// Let the initial point be ( X_0 : Y_0 : Z_0 : T_0).
// The rotated point is then (iY_0 : iX_0 : Z_0 : -T_0).
//
// We want to relate the computation of:
//
// invsqrt(Z^2 - Y^2) = invsqrt(Z_0^2 - Y_0^2) [non-rotated]
// invsqrt(Z^2 - Y^2) = invsqrt(Z_0^2 + X_0^2) [rotated]
//
// The curve equation in extended coordinates is
//
// 0 = (-X^2 + Y^2)*Z^2 - Z^4 - d*X^2*Y^2,
//
// so
// 0 = (-X^2 + Y^2)*Z^2 - Z^4 - d*T^2*Z^2 since XY=TZ
// = (-X^2 + Y^2 - Z^2 - d*T^2)*Z^2
// = ( X^2 - Y^2 + Z^2 + d*T^2)*Z^2 mult by -1
// -T^2*Z^2 = (X^2 - Y^2 + Z^2 + d*T^2)*Z^2 - T^2*Z^2 sub T^2*Z^2
// -T^2*Z^2 = (X^2 - Y^2 + Z^2 - T^2)*Z^2 + d*T^2*Z^2
// (-1-d)*T^2*Z^2 = (X^2 - Y^2 + Z^2 - T^2)*Z^2
//
// for any point (X:Y:Z:T) in extended coordinates. Therefore,
//
// (Z^2 - Y^2)*(Z^2 + X^2) = Z^4 + Z^2*X^2 - Y^2*Z^2 - Y^2*X^2
// = Z^4 + Z^2*X^2 - Y^2*Z^2 - T^2*Z^2 since XY=TZ
// = Z^2*(X^2 - Y^2 + Z^2 - T^2)
// = (-1-d)*T^2*Z^2.
//
// Taking square roots of both sides and rearranging, we get
//
// invsqrt(Z^2 - Y^2) = invsqrt(-1-d)*(Z^2+X^2)*(1/TZ)*invsqrt(Z^2+X^2)
// `-----------' `---------------------'
// curve constant batchable
//
// for any point (X:Y:Z:T) in extended coordinates.
//
// Therefore, we can do the computation with only one inverse
// square root like so:
//
// W <--- invsqrt((T_0 * Z_0)^2 * (Z_0^2+X_0^2))
// = 1/(T_0 * Z_0 * sqrt(Z_0^2 + X_0^2))
//
// xy <--- T_0^2 * W^2 * (T_0 * Z_0) * (Z_0^2 + X_0^2)
// = T_0 / Z_0 = xy
//
// if Y_0 nonzero and xy nonnegative:
// (X : Y : Z : T) <--- (X_0 : Y_0 : Z_0 : T_0)
// r <--- (1/(-1-d)) * (Z_0^2 + X_0^2) * W
// = invsqrt(a-d) * invsqrt(Z_0^2 - Y_0^2) since a = -1
// = invsqrt(a-d) * invsqrt(Z^2 - Y^2)
// otherwise:
// (X : Y : Z : T) <--- (i*Y_0 : i*X_0 : Z_0 : -T_0)
// r <--- invsqrt(a-d) * (T_0 * Z_0) * W
// = invsqrt(a-d) * invsqrt(Z_0^2 + X_0^2)
// = invsqrt(a-d) * invsqrt(Z^2 - Y^2)
//
// The rest of the compression follows the steps in the
// appendix of the Decaf paper.
let mut X = self.0.X;
let mut Y = self.0.Y;
let mut T = self.0.T;
let Z = &self.0.Z;
let TZ = &T * Z;
let ZZ_plus_XX = &Z.square() + &X.square();
let tmp = &TZ.square() * &ZZ_plus_XX;
let (tmp_is_nonzero_square, W) = tmp.invsqrt();
// tmp should always be a square (why? related to being in the
// image of the isogeny?)
debug_assert_eq!(tmp_is_nonzero_square | tmp.is_zero(), 1u8);
let xy = &T.square() * &(&W.square() * &(&TZ * &ZZ_plus_XX));
let rotate = 1u8 & !(Y.is_nonzero() & xy.is_nonnegative_decaf());
let mut r = &W * &(&ZZ_plus_XX * &constants::inv_a_minus_d);
let r_rot = &W * &(&TZ * &constants::invsqrt_a_minus_d);
let iX = &X * &constants::SQRT_M1;
let iY = &Y * &constants::SQRT_M1;
r.conditional_assign(&r_rot, rotate);
X.conditional_assign(&iY, rotate);
Y.conditional_assign(&iX, rotate);
T.conditional_negate(rotate);
// Step 2: Compute u = (a-d)r
let u = &constants::a_minus_d * &r;
// Step 3: Negate r if -2uZ is negative.
let uZ = &u * Z;
let m2uZ = -&(&uZ + &uZ);
r.conditional_negate(m2uZ.is_negative_decaf());
// Step 4: Compute s = | u(r(aZX - dYT)+Y)/a|
// = |u(r(-ZX - dYT)+Y)| since a = -1
let minus_ZX = -&(Z * &X);
let dYT = &constants::d * &(&Y * &T);
// Compute s = u(r(aZX - dYT)+Y) and cnegate for abs
let mut s = &u * &(&(&r * &(&minus_ZX - &dYT)) + &Y);
let neg = s.is_negative_decaf();
s.conditional_negate(neg);
CompressedDecaf(s.to_bytes())
}
/// Return the coset self + E[4], for debugging.
fn coset4(&self) -> [ExtendedPoint; 4] {
[ self.0
, &self.0 + &constants::EIGHT_TORSION[2]
, &self.0 + &constants::EIGHT_TORSION[4]
, &self.0 + &constants::EIGHT_TORSION[6]
]
}
/// Computes the Elligator map as described in the Decaf paper.
///
/// # Note
///
/// This method is not public because it's just used for hashing
/// to a point -- proper elligator support is deferred for now.
pub fn elligator_decaf_flavour(r_0: &FieldElement) -> DecafPoint {
// Follows Appendix C of the Decaf paper.
// Use n = 2 as the quadratic nonresidue so that n*x = x + x.
let minus_one = -&FieldElement::one();
// 1. Compute r <--- nr_0^2.
let r_0_squared = r_0.square();
let r = &r_0_squared + &r_0_squared;
// 2. Compute D <--- (dr + (a-d)) * (dr - (d + ar))
let dr = &constants::d * &r;
// D = (dr + (a-d)) * (dr - (d + ar))
// = (dr + (a-d)) * (dr - (d-r)) since a=-1
// writing as
// = (dr + (a-d)) * dr - (dr + (a-d)) * (d - r)
// avoids two consecutive additions (could cause overflow)
let dr_plus_amd = &dr + &constants::a_minus_d;
let D = &(&dr_plus_amd * &dr) - &(&dr_plus_amd * &(&constants::d - &r));
// 3. Compute N <--- (r+1) * (a-2d)
let N = &(&r + &FieldElement::one()) * &(&minus_one - &constants::d2);
// 4. Compute
// / +1, 1 / sqrt(ND) if ND is square
// c, e <--- | +1, 0 if N or D = 0
// \ -1, nr_0 / sqrt(nND) otherwise
let ND = &N * &D;
let nND = &ND + &ND;
let mut c = FieldElement::one();
let mut e = FieldElement::zero();
let (ND_is_nonzero_square, ND_invsqrt) = ND.invsqrt();
e.conditional_assign(&ND_invsqrt, ND_is_nonzero_square);
let (nND_is_nonzero_square, nND_invsqrt) = nND.invsqrt();
let nr_0_nND_invsqrt = &nND_invsqrt * &(r_0 + r_0);
c.conditional_assign(&minus_one, nND_is_nonzero_square);
e.conditional_assign(&nr_0_nND_invsqrt, nND_is_nonzero_square);
// 5. Compute s <--- c*|N*e|
let mut s = &N * &e;
let neg = s.is_negative_decaf();
s.conditional_negate(neg);
s *= &c;
// 6. Compute t <--- -c*N*(r-1)* ((a-2d)*e)^2 -1
let a_minus_2d_e_sq = (&(&minus_one - &constants::d2) * &e).square();
let c_N_r_minus_1 = &c * &(&N * &(&r + &minus_one));
let t = &minus_one - &(&c_N_r_minus_1 * &a_minus_2d_e_sq);
// 7. Apply the isogeny:
// (x,y) = ((2s)/(1+as^2), (1-as^2)/(t))
let as_sq = &minus_one * &s.square();
let P = CompletedPoint{
X: &s + &s,
Z: &FieldElement::one() + &as_sq,
Y: &FieldElement::one() - &as_sq,
T: t,
};
// Convert to extended and return.
DecafPoint(P.to_extended())
}
/// Return a `DecafPoint` chosen uniformly at random using a user-provided RNG.
///
/// # Inputs
///
/// * `rng`: any RNG which implements the `rand::Rng` interface.
///
/// # Returns
///
/// A random element of the Decaf group.
///
/// # Implementation
///
/// Uses the Decaf-flavoured Elligator 2 map, so that the discrete log of the
/// output point with respect to any other point should be unknown.
#[cfg(feature = "std")]
pub fn random<T: Rng>(rng: &mut T) -> Self {
let mut field_bytes = [0u8; 32];
rng.fill_bytes(&mut field_bytes);
let r_0 = FieldElement::from_bytes(&field_bytes);
DecafPoint::elligator_decaf_flavour(&r_0)
}
/// Hash a slice of bytes into a `DecafPoint`.
///
/// Takes a type parameter `D`, which is any `Digest` producing 32
/// bytes (256 bits) of output.
///
/// Convenience wrapper around `from_hash`.
///
/// # Implementation
///
/// Uses the Decaf-flavoured Elligator 2 map, so that the discrete log of the
/// output point with respect to any other point should be unknown.
///
/// # Example
///
/// ```
/// # extern crate curve25519_dalek;
/// # use curve25519_dalek::decaf::DecafPoint;
/// extern crate sha2;
/// use sha2::Sha256;
///
/// # // Need fn main() here in comment so the doctest compiles
/// # // See https://doc.rust-lang.org/book/documentation.html#documentation-as-tests
/// # fn main() {
/// let msg = "To really appreciate architecture, you may even need to commit a murder";
/// let P = DecafPoint::hash_from_bytes::<Sha256>(msg.as_bytes());
/// # }
/// ```
///
pub fn hash_from_bytes<D>(input: &[u8]) -> DecafPoint
where D: Digest<OutputSize = U32> + Default
{
let mut hash = D::default();
hash.input(input);
DecafPoint::from_hash(hash)
}
/// Construct a `DecafPoint` from an existing `Digest` instance.
///
/// Use this instead of `hash_from_bytes` if it is more convenient
/// to stream data into the `Digest` than to pass a single byte
/// slice.
pub fn from_hash<D>(hash: D) -> DecafPoint
where D: Digest<OutputSize = U32> + Default
{
// XXX this seems clumsy
let mut output = [0u8; 32];
output.copy_from_slice(hash.result().as_slice());
let r_0 = FieldElement::from_bytes(&output);
DecafPoint::elligator_decaf_flavour(&r_0)
}
}
impl Identity for DecafPoint {
fn identity() -> DecafPoint {
DecafPoint(ExtendedPoint::identity())
}
}
// ------------------------------------------------------------------------
// Equality
// ------------------------------------------------------------------------
/// XXX check whether there's a simple way to do equality checking
/// with cofactor 8, not just cofactor 4, and add a CT equality function?
impl PartialEq for DecafPoint {
fn eq(&self, other: &DecafPoint) -> bool {
let self_compressed = self.compress();
let other_compressed = other.compress();
self_compressed == other_compressed
}
}
impl Eq for DecafPoint {}
// ------------------------------------------------------------------------
// Arithmetic
// ------------------------------------------------------------------------
impl<'a, 'b> Add<&'b DecafPoint> for &'a DecafPoint {
type Output = DecafPoint;
fn add(self, other: &'b DecafPoint) -> DecafPoint {
DecafPoint(&self.0 + &other.0)
}
}
impl<'b> AddAssign<&'b DecafPoint> for DecafPoint {
fn add_assign(&mut self, _rhs: &DecafPoint) {
*self = (self as &DecafPoint) + _rhs;
}
}
impl<'a, 'b> Sub<&'b DecafPoint> for &'a DecafPoint {
type Output = DecafPoint;
fn sub(self, other: &'b DecafPoint) -> DecafPoint {
DecafPoint(&self.0 - &other.0)
}
}
impl<'b> SubAssign<&'b DecafPoint> for DecafPoint {
fn sub_assign(&mut self, _rhs: &DecafPoint) {
*self = (self as &DecafPoint) - _rhs;
}
}
impl<'a> Neg for &'a DecafPoint {
type Output = DecafPoint;
fn neg(self) -> DecafPoint {
DecafPoint(-&self.0)
}
}
impl<'b> MulAssign<&'b Scalar> for DecafPoint {
fn mul_assign(&mut self, scalar: &'b Scalar) {
let result = (self as &DecafPoint) * scalar;
*self = result;
}
}
impl<'a, 'b> Mul<&'b Scalar> for &'a DecafPoint {
type Output = DecafPoint;
/// Scalar multiplication: compute `scalar * self`.
fn mul(self, scalar: &'b Scalar) -> DecafPoint {
DecafPoint(&self.0 * scalar)
}
}
impl<'a, 'b> Mul<&'b DecafPoint> for &'a Scalar {
type Output = DecafPoint;
/// Scalar multiplication: compute `self * scalar`.
fn mul(self, point: &'b DecafPoint) -> DecafPoint {
DecafPoint(self * &point.0)
}
}
/// Given a vector of (possibly secret) scalars and a vector of
/// (possibly secret) points, compute `c_1 P_1 + ... + c_n P_n`.
///
/// This function has the same behaviour as
/// `vartime::multiscalar_mult` but is constant-time.
///
/// # Input
///
/// A vector of `Scalar`s and a vector of `DecafPoints`. It is an
/// error to call this function with two vectors of different lengths.
#[cfg(any(feature = "alloc", feature = "std"))]
pub fn multiscalar_mult<'a, 'b, I, J>(scalars: I, points: J) -> DecafPoint
where I: IntoIterator<Item = &'a Scalar>,
J: IntoIterator<Item = &'b DecafPoint>,
{
let extended_points = points.into_iter().map(|P| &P.0);
DecafPoint(edwards::multiscalar_mult(scalars, extended_points))
}
/// Precomputation
#[derive(Clone)]
pub struct DecafBasepointTable(pub EdwardsBasepointTable);
impl<'a, 'b> Mul<&'b Scalar> for &'a DecafBasepointTable {
type Output = DecafPoint;
fn mul(self, scalar: &'b Scalar) -> DecafPoint {
DecafPoint(&self.0 * scalar)
}
}
impl<'a, 'b> Mul<&'a DecafBasepointTable> for &'b Scalar {
type Output = DecafPoint;
fn mul(self, basepoint_table: &'a DecafBasepointTable) -> DecafPoint {
DecafPoint(self * &basepoint_table.0)
}
}
impl DecafBasepointTable {
/// Create a precomputed table of multiples of the given `basepoint`.
pub fn create(basepoint: &DecafPoint) -> DecafBasepointTable {
DecafBasepointTable(EdwardsBasepointTable::create(&basepoint.0))
}
/// Get the basepoint for this table as a `DecafPoint`.
pub fn basepoint(&self) -> DecafPoint {
DecafPoint(self.0.basepoint())
}
}
// ------------------------------------------------------------------------
// Constant-time conditional assignment
// ------------------------------------------------------------------------
impl ConditionallyAssignable for DecafPoint {
/// Conditionally assign `other` to `self`, if `choice == 1u8`.
///
/// # Example
///
/// ```
/// # extern crate subtle;
/// # extern crate curve25519_dalek;
/// #
/// # use subtle::ConditionallyAssignable;
/// #
/// # use curve25519_dalek::edwards::Identity;
/// # use curve25519_dalek::decaf::DecafPoint;
/// # use curve25519_dalek::constants;
/// # fn main() {
/// let A = DecafPoint::identity();
/// let B = constants::DECAF_ED25519_BASEPOINT_POINT;
///
/// let mut P = A;
///
/// P.conditional_assign(&B, 0u8);
/// assert!(P == A);
/// P.conditional_assign(&B, 1u8);
/// assert!(P == B);
/// # }
/// ```
fn conditional_assign(&mut self, other: &DecafPoint, choice: u8) {
self.0.X.conditional_assign(&other.0.X, choice);
self.0.Y.conditional_assign(&other.0.Y, choice);
self.0.Z.conditional_assign(&other.0.Z, choice);
self.0.T.conditional_assign(&other.0.T, choice);
}
}
// ------------------------------------------------------------------------
// Debug traits
// ------------------------------------------------------------------------
impl Debug for CompressedDecaf {
fn fmt(&self, f: &mut ::core::fmt::Formatter) -> ::core::fmt::Result {
write!(f, "CompressedDecaf: {:?}", self.as_bytes())
}
}
impl Debug for DecafPoint {
fn fmt(&self, f: &mut ::core::fmt::Formatter) -> ::core::fmt::Result {
let coset = self.coset4();
write!(f, "DecafPoint: coset \n{:?}\n{:?}\n{:?}\n{:?}",
coset[0], coset[1], coset[2], coset[3])
}
}
// ------------------------------------------------------------------------
// Variable-time functions
// ------------------------------------------------------------------------
pub mod vartime {
//! Variable-time operations on decaf points, useful for non-secret data.
use super::*;
/// Given a vector of public scalars and a vector of (possibly secret)
/// points, compute
///
/// c_1 P_1 + ... + c_n P_n.
///
/// # Input
///
/// A vector of `Scalar`s and a vector of `ExtendedPoints`. It is an
/// error to call this function with two vectors of different lengths.
pub fn multiscalar_mult<'a, 'b, I, J>(scalars: I, points: J) -> DecafPoint
where I: IntoIterator<Item = &'a Scalar>,
J: IntoIterator<Item = &'b DecafPoint>
{
let extended_points = points.into_iter().map(|P| &P.0);
DecafPoint(edwards::vartime::multiscalar_mult(scalars, extended_points))
}
}
// ------------------------------------------------------------------------
// Tests
// ------------------------------------------------------------------------
#[cfg(test)]
mod test {
use rand::OsRng;
use scalar::Scalar;
use constants;
use edwards::CompressedEdwardsY;
use edwards::Identity;
use edwards::ValidityCheck;
use super::*;
#[cfg(feature = "serde")]
use serde_cbor;
#[test]
#[cfg(feature = "serde")]
fn serde_cbor_basepoint_roundtrip() {
let output = serde_cbor::to_vec(&constants::DECAF_ED25519_BASEPOINT_POINT).unwrap();
let parsed: DecafPoint = serde_cbor::from_slice(&output).unwrap();
assert_eq!(parsed, constants::DECAF_ED25519_BASEPOINT_POINT);
}
#[test]
fn decaf_decompress_negative_s_fails() {
// constants::d is neg, so decompression should fail as |d| != d.
let bad_compressed = CompressedDecaf(constants::d.to_bytes());
assert!(bad_compressed.decompress().is_none());
}
#[test]
fn decaf_decompress_id() {
let compressed_id = CompressedDecaf::identity();
let id = compressed_id.decompress().unwrap();
assert_eq!(id.0.compress(), CompressedEdwardsY::identity());
}
#[test]
fn decaf_compress_id() {
let id = DecafPoint::identity();
assert_eq!(id.compress(), CompressedDecaf::identity());
}
#[test]
fn decaf_basepoint_roundtrip() {
let bp_compressed_decaf = constants::DECAF_ED25519_BASEPOINT_POINT.compress();
let bp_recaf = bp_compressed_decaf.decompress().unwrap().0;
// Check that bp_recaf differs from bp by a point of order 4
let diff = &constants::ED25519_BASEPOINT_POINT - &bp_recaf;
let diff4 = diff.mult_by_pow_2(4); // XXX this is wrong
assert_eq!(diff4.compress(), CompressedEdwardsY::identity());
}
#[test]
fn encodings_of_small_multiples_of_basepoint() {
// Table of encodings of (1+i)*basepoint
// Generated using the previous naive implementation.
let compressed = [
CompressedDecaf([141, 190, 226, 107, 177, 201, 35, 118, 14, 55, 160, 165, 242, 207, 121, 161, 177, 80, 8, 132, 205, 254, 101, 169, 233, 65, 124, 96, 255, 182, 249, 40]),
CompressedDecaf([131, 57, 148, 16, 8, 196, 141, 82, 144, 220, 105, 112, 66, 33, 48, 16, 182, 198, 173, 35, 248, 181, 92, 231, 222, 35, 85, 56, 5, 252, 91, 40]),
CompressedDecaf([199, 132, 32, 144, 156, 143, 81, 170, 240, 56, 232, 6, 178, 37, 118, 190, 110, 201, 26, 173, 156, 97, 59, 162, 240, 247, 226, 107, 197, 111, 107, 26]),
CompressedDecaf([210, 120, 34, 214, 175, 27, 61, 6, 229, 181, 216, 36, 11, 245, 146, 232, 130, 215, 77, 29, 210, 30, 54, 155, 191, 81, 59, 124, 174, 3, 135, 36]),
CompressedDecaf([155, 52, 159, 52, 189, 27, 181, 0, 245, 131, 0, 197, 79, 208, 252, 122, 104, 161, 245, 143, 67, 94, 13, 129, 153, 173, 129, 179, 118, 231, 90, 52]),
CompressedDecaf([42, 117, 252, 118, 8, 1, 72, 25, 111, 246, 247, 103, 236, 86, 235, 29, 100, 156, 186, 209, 159, 21, 61, 26, 249, 25, 137, 228, 84, 23, 10, 27]),
CompressedDecaf([21, 126, 181, 117, 58, 90, 216, 28, 184, 57, 9, 23, 158, 68, 159, 171, 109, 150, 232, 140, 144, 73, 139, 122, 124, 105, 125, 160, 94, 185, 150, 52]),
CompressedDecaf([232, 167, 112, 233, 126, 33, 105, 63, 151, 6, 88, 225, 181, 17, 223, 12, 116, 138, 203, 47, 243, 225, 50, 171, 21, 220, 186, 179, 132, 20, 48, 6]),
CompressedDecaf([99, 44, 97, 48, 242, 174, 78, 198, 112, 154, 146, 36, 239, 34, 94, 4, 0, 244, 175, 34, 46, 0, 83, 187, 5, 163, 225, 63, 51, 237, 234, 22]),
CompressedDecaf([2, 33, 89, 176, 178, 123, 159, 75, 235, 172, 251, 11, 137, 177, 90, 122, 149, 186, 52, 243, 153, 190, 185, 202, 59, 137, 204, 160, 150, 152, 148, 55]),
CompressedDecaf([245, 79, 78, 226, 114, 69, 247, 112, 18, 54, 90, 225, 176, 77, 231, 235, 196, 123, 49, 221, 34, 205, 151, 228, 244, 112, 82, 58, 30, 31, 58, 12]),
CompressedDecaf([135, 53, 175, 167, 13, 94, 62, 31, 29, 248, 13, 132, 29, 69, 7, 188, 145, 49, 62, 55, 181, 109, 214, 11, 248, 162, 70, 15, 236, 126, 100, 60]),
CompressedDecaf([98, 150, 69, 229, 144, 122, 237, 107, 127, 177, 33, 64, 59, 173, 210, 102, 74, 34, 23, 16, 252, 117, 14, 97, 231, 178, 63, 193, 157, 28, 178, 17]),
CompressedDecaf([222, 104, 6, 1, 72, 12, 72, 178, 204, 238, 128, 70, 41, 150, 235, 96, 153, 150, 18, 4, 141, 206, 0, 38, 122, 112, 249, 51, 94, 251, 20, 57]),
CompressedDecaf([7, 221, 140, 57, 13, 146, 248, 27, 56, 4, 128, 23, 145, 120, 126, 4, 158, 173, 52, 213, 164, 250, 26, 55, 89, 96, 187, 111, 211, 18, 63, 19]),
CompressedDecaf([91, 213, 193, 10, 102, 92, 199, 124, 61, 176, 1, 47, 111, 59, 183, 91, 79, 56, 208, 109, 172, 209, 17, 167, 229, 216, 3, 236, 200, 208, 15, 20]),
];
let mut bp = constants::DECAF_ED25519_BASEPOINT_POINT;
for i in 0..16 {
assert_eq!(bp.compress(), compressed[i]);
bp = &bp + &constants::DECAF_ED25519_BASEPOINT_POINT;
}
}
#[test]
fn decaf_four_torsion_basepoint() {
let bp = constants::DECAF_ED25519_BASEPOINT_POINT;
let bp_coset = bp.coset4();
for i in 0..4 {
assert_eq!(bp, DecafPoint(bp_coset[i]));
}
}
#[test]
fn decaf_four_torsion_random() {
let mut rng = OsRng::new().unwrap();
let B = &constants::DECAF_ED25519_BASEPOINT_TABLE;
let P = B * &Scalar::random(&mut rng);
let P_coset = P.coset4();
for i in 0..4 {
assert_eq!(P, DecafPoint(P_coset[i]));
}
}
#[test]
fn decaf_random_roundtrip() {
let mut rng = OsRng::new().unwrap();
let B = &constants::DECAF_ED25519_BASEPOINT_TABLE;
for _ in 0..100 {
let P = B * &Scalar::random(&mut rng);
let compressed_P = P.compress();
let Q = compressed_P.decompress().unwrap();
assert_eq!(P, Q);
}
}
#[test]
fn decaf_random_is_valid() {
let mut rng = OsRng::new().unwrap();
for _ in 0..100 {
let P = DecafPoint::random(&mut rng);
// Check that P is on the curve
assert!(P.0.is_valid());
// Check that P is in the image of the decaf map
P.compress();
}
}
}
#[cfg(all(test, feature = "bench"))]
mod bench {
use rand::OsRng;
use test::Bencher;
use super::*;
#[bench]
fn decompression(b: &mut Bencher) {
let mut rng = OsRng::new().unwrap();
let B = &constants::DECAF_ED25519_BASEPOINT_TABLE;
let P = B * &Scalar::random(&mut rng);
let P_compressed = P.compress();
b.iter(|| P_compressed.decompress().unwrap());
}
#[bench]
fn compression(b: &mut Bencher) {
let mut rng = OsRng::new().unwrap();
let B = &constants::DECAF_ED25519_BASEPOINT_TABLE;
let P = B * &Scalar::random(&mut rng);
b.iter(|| P.compress());
}
}

35
src/edwards.rs
View file

@ -135,15 +135,15 @@ impl CompressedEdwardsY {
let Y = FieldElement::from_bytes(self.as_bytes());
let Z = FieldElement::one();
let YY = Y.square();
let u = &YY - &Z; // u = y²-1
let v = &(&YY * &constants::d) + &Z; // v = dy²+1
let u = &YY - &Z; // u = y²-1
let v = &(&YY * &constants::EDWARDS_D) + &Z; // v = dy²+1
let (is_nonzero_square, mut X) = FieldElement::sqrt_ratio(&u, &v);
if is_nonzero_square != 1u8 { return None; }
// Flip the sign of X if it's not correct
let compressed_sign_bit = self.as_bytes()[31] >> 7;
let current_sign_bit = X.is_negative_ed25519();
let current_sign_bit = X.is_negative();
X.conditional_negate(current_sign_bit ^ compressed_sign_bit);
Some(ExtendedPoint{ X: X, Y: Y, Z: Z, T: &X * &Y })
@ -342,7 +342,7 @@ impl ValidityCheck for ProjectivePoint {
let ZZ = self.Z.square();
let ZZZZ = ZZ.square();
let lhs = &(&YY - &XX) * &ZZ;
let rhs = &ZZZZ + &(&constants::d * &(&XX * &YY));
let rhs = &ZZZZ + &(&constants::EDWARDS_D * &(&XX * &YY));
lhs == rhs
}
@ -442,7 +442,7 @@ impl ProjectivePoint {
let mut s: [u8; 32];
s = y.to_bytes();
s[31] ^= (x.is_negative_ed25519() << 7) as u8;
s[31] ^= (x.is_negative() << 7) as u8;
CompressedEdwardsY(s)
}
@ -518,7 +518,7 @@ impl ExtendedPoint {
Y_plus_X: &self.Y + &self.X,
Y_minus_X: &self.Y - &self.X,
Z: self.Z,
T2d: &self.T * &constants::d2,
T2d: &self.T * &constants::EDWARDS_D2,
}
}
@ -541,7 +541,7 @@ impl ExtendedPoint {
let recip = self.Z.invert();
let x = &self.X * &recip;
let y = &self.Y * &recip;
let xy2d = &(&x * &y) * &constants::d2;
let xy2d = &(&x * &y) * &constants::EDWARDS_D2;
AffineNielsPoint{
y_plus_x: &y + &x,
y_minus_x: &y - &x,
@ -1233,6 +1233,7 @@ pub mod vartime {
}
let odd_multiples_of_A = OddMultiples::create(A);
let odd_multiples_of_B = &constants::AFFINE_ODD_MULTIPLES_OF_BASEPOINT;
let mut r = ProjectivePoint::identity();
loop {
@ -1245,9 +1246,9 @@ pub mod vartime {
}
if b_naf[i] > 0 {
t = &t.to_extended() + &constants::bi[( b_naf[i]/2) as usize];
t = &t.to_extended() + &odd_multiples_of_B[( b_naf[i]/2) as usize];
} else if b_naf[i] < 0 {
t = &t.to_extended() - &constants::bi[(-b_naf[i]/2) as usize];
t = &t.to_extended() - &odd_multiples_of_B[(-b_naf[i]/2) as usize];
}
r = t.to_projective();
@ -1269,8 +1270,6 @@ pub mod vartime {
#[cfg(test)]
mod test {
#[cfg(feature = "yolocrypto")]
use decaf::DecafPoint;
use field::FieldElement;
use scalar::Scalar;
use subtle::ConditionallyAssignable;
@ -1433,7 +1432,7 @@ mod test {
#[test]
fn basepoint_mult_by_basepoint_order() {
let B = &constants::ED25519_BASEPOINT_TABLE;
let should_be_id = B * &constants::l;
let should_be_id = B * &constants::BASEPOINT_ORDER;
assert!(should_be_id.is_identity());
}
@ -1550,18 +1549,6 @@ mod test {
assert!(P1.compress().to_bytes() == P2.compress().to_bytes());
}
#[test]
#[cfg(feature = "yolocrypto")]
fn scalarmult_decafpoint_works_both_ways() {
let P: DecafPoint = DecafPoint(constants::ED25519_BASEPOINT_POINT);
let s: Scalar = A_SCALAR;
let P1 = &P * &s;
let P2 = &s * &P;
assert!(P1.compress().as_bytes() == P2.compress().as_bytes());
}
mod vartime {
use super::super::*;
use super::{A_SCALAR, B_SCALAR, A_TIMES_BASEPOINT, DOUBLE_SCALAR_MULT_RESULT};

56
src/field.rs
View file

@ -23,6 +23,7 @@ use core::cmp::{Eq, PartialEq};
use subtle::slices_equal;
use subtle::byte_is_nonzero;
use subtle::ConditionallyAssignable;
use subtle::ConditionallyNegatable;
use subtle::Equal;
use constants;
@ -81,44 +82,11 @@ impl FieldElement {
/// # Return
///
/// If negative, return `1u8`. Otherwise, return `0u8`.
pub fn is_negative_ed25519(&self) -> u8 { //FeIsNegative
pub fn is_negative(&self) -> u8 {
let bytes = self.to_bytes();
(bytes[0] & 1) as u8
}
/// Determine if this `FieldElement` is negative, in the
/// sense used by Decaf: `x` is nonnegative if the least
/// absolute residue for `x` lies in `[0, (p-1)/2]`, and
/// is negative otherwise.
///
/// # Return
///
/// Returns `1u8` if negative, `0u8` if nonnegative.
///
/// # Implementation
///
/// Uses a trick borrowed from Mike Hamburg's code. Let `x \in
/// F_p` and let `y \in Z` be the least absolute residue for `x`.
/// Suppose `y ≤ (p-1)/2`. Then `2y < p` so `2y = 2y mod p` and
/// `2y mod p` is even. On the other hand, if `y > (p-1)/2` then
/// `2y ≥ p`; since `y < p`, `2y \in [p, 2p)`, so `2y mod p =
/// 2y-p`, which is odd.
///
/// Thus we can test whether `y ≤ (p-1)/2` by checking whether `2y
/// mod p` is even.
pub fn is_negative_decaf(&self) -> u8 {
let y = self + self;
(y.to_bytes()[0] & 1) as u8
}
/// Determine if this `FieldElement` is nonnegative, in the
/// sense used by Decaf: `x` is nonnegative if the least
/// absolute residue for `x` lies in `[0, (p-1)/2]`, and
/// is negative otherwise.
pub fn is_nonnegative_decaf(&self) -> u8 {
1u8 & (!self.is_negative_decaf())
}
/// Determine if this `FieldElement` is zero.
///
/// # Return
@ -232,6 +200,8 @@ impl FieldElement {
/// Given `FieldElements` `u` and `v`, attempt to compute
/// `sqrt(u/v)` in constant time.
///
/// This function always returns the nonnegative square root, if it exists.
///
/// It would be much better to use an `Option` type here, but
/// doing so forces the caller to branch, which we don't want to
/// do. This seems like the least bad solution.
@ -278,6 +248,10 @@ impl FieldElement {
let r_prime = &constants::SQRT_M1 * &r;
r.conditional_assign(&r_prime, flipped_sign_sqrt);
// Choose the nonnegative square root.
let r_is_negative = r.is_negative();
r.conditional_negate(r_is_negative);
let was_nonzero_square = correct_sign_sqrt | flipped_sign_sqrt;
(was_nonzero_square, r)
@ -465,6 +439,20 @@ mod test {
x.conditional_negate(1u8);
assert_eq!(x, one);
}
#[test]
fn encoding_is_canonical() {
// Encode 1 wrongly as 1 + (2^255 - 19) = 2^255 - 18
let one_encoded_wrongly_bytes: [u8;32] = [0xee, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0x7f];
// Decode to a field element
let one = FieldElement::from_bytes(&one_encoded_wrongly_bytes);
// .. then check that the encoding is correct
let one_bytes = one.to_bytes();
assert_eq!(one_bytes[0], 1);
for i in 1..32 {
assert_eq!(one_bytes[i], 0);
}
}
}
#[cfg(all(test, feature = "bench"))]

14
src/lib.rs
View file

@ -38,9 +38,8 @@
#[cfg(all(test, feature = "bench"))]
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,12 +70,15 @@ 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;
// Feature gate decaf while our implementation is unfinished and probably incorrect.
#[cfg(feature = "yolocrypto")]
pub mod decaf;
pub mod ristretto;
// Other miscelaneous utilities.

5
src/montgomery.rs
View file

@ -146,8 +146,9 @@ impl CompressedMontgomeryU {
/// actually a square and `0` if otherwise, along with a `FieldElement`: the
/// Montgomery `v` corresponding to this `u`.
pub fn to_montgomery_v(u: &FieldElement) -> (u8, FieldElement) {
let A = &constants::MONTGOMERY_A;
let one: FieldElement = FieldElement::one();
let v_squared: FieldElement = u * &(&u.square() + &(&(&constants::A * u) + &one));
let v_squared: FieldElement = u * &(&u.square() + &(&(A * u) + &one));
let (okay, v_inv) = v_squared.invsqrt();
let v = &v_inv * &v_squared;
@ -183,7 +184,7 @@ impl CompressedMontgomeryU {
pub fn to_edwards_x(u: &FieldElement, v: &FieldElement, sign: &u8) -> FieldElement {
let mut x: FieldElement = &(u * &v.invert()) * &constants::SQRT_MINUS_APLUS2;
let neg_x: FieldElement = -(&x);
let current_sign: u8 = x.is_negative_ed25519();
let current_sign: u8 = x.is_negative();
// Negate x to match the sign:
x.conditional_assign(&neg_x, current_sign ^ sign);

1245
src/ristretto.rs Normal file
View file

File diff suppressed because it is too large Load diff

502
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,
@ -857,6 +582,36 @@ mod test {
0,0,0,11,0,0,0,0,0,15,0,0,0,0,0,-9,0,0,0,0,0,0,0,-1,0,0,0,0,0,0,0,7,
0,0,0,0,0,-15,0,0,0,0,0,15,0,0,0,0,15,0,0,0,0,15,0,0,0,0,0,1,0,0,0,0];
#[test]
fn fuzzer_testcase_reduction() {
// LE bytes of 24519928653854221733733552434404946937899825954937634815
let a_bytes = [255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 0, 0, 0, 0, 0, 0, 0, 0, 0];
// LE bytes of 4975441334397345751130612518500927154628011511324180036903450236863266160640
let b_bytes = [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 255, 210, 210, 210, 255, 255, 255, 255, 10];
// LE bytes of 6432735165214683820902750800207468552549813371247423777071615116673864412038
let c_bytes = [134, 171, 119, 216, 180, 128, 178, 62, 171, 132, 32, 62, 34, 119, 104, 193, 47, 215, 181, 250, 14, 207, 172, 93, 75, 207, 211, 103, 144, 204, 56, 14];
let a = Scalar(a_bytes);
let b = Scalar(b_bytes);
let c = Scalar(c_bytes);
let mut tmp = [0u8; 64];
// also_a = (a mod l)
tmp[0..32].copy_from_slice(&a_bytes[..]);
let also_a = Scalar::reduce(&tmp);
// also_b = (b mod l)
tmp[0..32].copy_from_slice(&b_bytes[..]);
let also_b = Scalar::reduce(&tmp);
let expected_c = &a * &b;
let also_expected_c = &also_a * &also_b;
assert_eq!(c, expected_c);
assert_eq!(c, also_expected_c);
}
#[test]
fn non_adjacent_form() {
let naf = A_SCALAR.non_adjacent_form();
@ -900,7 +655,7 @@ mod test {
#[test]
fn impl_sub() {
let should_be_one = &constants::l - &constants::l_minus_1;
let should_be_one = &constants::BASEPOINT_ORDER - &constants::BASEPOINT_ORDER_MINUS_1;
assert_eq!(should_be_one, Scalar::one());
}
@ -952,6 +707,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());
}
@ -966,6 +722,47 @@ mod test {
assert_eq!(should_be_X, X);
}
#[test]
fn to_bytes_from_bytes_roundtrips() {
let unpacked = X.unpack();
let bytes = unpacked.to_bytes();
let should_be_unpacked = UnpackedScalar::from_bytes(&bytes);
assert_eq!(should_be_unpacked.0, unpacked.0);
}
#[test]
fn montgomery_reduce_matches_reduce() {
let mut bignum = [0u8; 64];
// set bignum = x + 2^256x
for i in 0..32 {
bignum[ i] = X[i];
bignum[32+i] = X[i];
}
// x + 2^256x (mod l)
// = 3958878930004874126169954872055634648693766179881526445624823978500314864344
let expected = Scalar([216, 154, 179, 139, 210, 121, 2, 71,
69, 99, 158, 216, 23, 173, 63, 100,
204, 0, 91, 50, 219, 153, 57, 249,
28, 82, 31, 197, 100, 165, 192, 8]);
let reduced = Scalar::reduce(&bignum);
// The reduced scalar should match the expected
assert_eq!(reduced.0, expected.0);
// (x + 2^256x) * R
let interim = UnpackedScalar::mul_internal(&UnpackedScalar::from_bytes_wide(&bignum),
&constants::R);
// ((x + 2^256x) * R) / R (mod l)
let montgomery_reduced = UnpackedScalar::montgomery_reduce(&interim);
// The Montgomery reduced scalar should match the reduced one, as well as the expected
assert_eq!(montgomery_reduced.0, reduced.unpack().0);
assert_eq!(montgomery_reduced.0, expected.unpack().0)
}
#[cfg(feature = "serde")]
use serde_cbor;
@ -984,7 +781,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 +790,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));
}
}

555
src/scalar_32bit.rs Normal file
View file

@ -0,0 +1,555 @@
//! 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.
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 / 256 bit scalar into 9 29-bit limbs.
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 << 24) - 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 ;
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)]
pub (crate) 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)]
pub (crate) 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));
}
}

473
src/scalar_64bit.rs Normal file
View file

@ -0,0 +1,473 @@
//! 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).
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.
pub fn from_bytes(bytes: &[u8; 32]) -> Scalar64 {
let mut words = [0u64; 4];
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 << 48) - 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 = [0u64; 8];
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)]
pub (crate) fn mul_internal(a: &Scalar64, b: &Scalar64) -> [u128; 9] {
let mut z = [0u128; 9];
z[0] = m(a[0],b[0]);
z[1] = m(a[0],b[1]) + m(a[1],b[0]);
z[2] = m(a[0],b[2]) + m(a[1],b[1]) + m(a[2],b[0]);
z[3] = m(a[0],b[3]) + m(a[1],b[2]) + m(a[2],b[1]) + m(a[3],b[0]);
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]);
z[5] = m(a[1],b[4]) + m(a[2],b[3]) + m(a[3],b[2]) + m(a[4],b[1]);
z[6] = m(a[2],b[4]) + m(a[3],b[3]) + m(a[4],b[2]);
z[7] = m(a[3],b[4]) + m(a[4],b[3]);
z[8] = m(a[4],b[4]);
z
}
/// 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)]
pub (crate) 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));
}
}

611
vendor/ristretto.sage vendored Normal file
View file

@ -0,0 +1,611 @@
import binascii
class InvalidEncodingException(Exception): pass
class NotOnCurveException(Exception): pass
class SpecException(Exception): pass
def lobit(x): return int(x) & 1
def hibit(x): return lobit(2*x)
def negative(x): return lobit(x)
def enc_le(x,n): return bytearray([int(x)>>(8*i) & 0xFF for i in xrange(n)])
def dec_le(x): return sum(b<<(8*i) for i,b in enumerate(x))
def randombytes(n): return bytearray([randint(0,255) for _ in range(n)])
def optimized_version_of(spec):
"""Decorator: This function is an optimized version of some specification"""
def decorator(f):
def wrapper(self,*args,**kwargs):
def pr(x):
if isinstance(x,bytearray): return binascii.hexlify(x)
else: return str(x)
try: spec_ans = getattr(self,spec,spec)(*args,**kwargs),None
except Exception as e: spec_ans = None,e
try: opt_ans = f(self,*args,**kwargs),None
except Exception as e: opt_ans = None,e
if spec_ans[1] is None and opt_ans[1] is not None:
raise
#raise SpecException("Mismatch in %s: spec returned %s but opt threw %s"
# % (f.__name__,str(spec_ans[0]),str(opt_ans[1])))
if spec_ans[1] is not None and opt_ans[1] is None:
raise
#raise SpecException("Mismatch in %s: spec threw %s but opt returned %s"
# % (f.__name__,str(spec_ans[1]),str(opt_ans[0])))
if spec_ans[0] != opt_ans[0]:
raise SpecException("Mismatch in %s: %s != %s"
% (f.__name__,pr(spec_ans[0]),pr(opt_ans[0])))
if opt_ans[1] is not None: raise
else: return opt_ans[0]
wrapper.__name__ = f.__name__
return wrapper
return decorator
def xsqrt(x,exn=InvalidEncodingException("Not on curve")):
"""Return sqrt(x)"""
if not is_square(x): raise exn
s = sqrt(x)
if negative(s): s=-s
return s
def isqrt(x,exn=InvalidEncodingException("Not on curve")):
"""Return 1/sqrt(x)"""
if x==0: return 0
if not is_square(x): raise exn
return 1/sqrt(x)
def isqrt_i(x):
"""Return 1/sqrt(x) or 1/sqrt(zeta * x)"""
if x==0: return True,0
gen = x.parent(-1)
while is_square(gen): gen = sqrt(gen)
if is_square(x): return True,1/sqrt(x)
else: return False,1/sqrt(x*gen)
class QuotientEdwardsPoint(object):
"""Abstract class for point an a quotiented Edwards curve; needs F,a,d,cofactor to work"""
def __init__(self,x=0,y=1):
x = self.x = self.F(x)
y = self.y = self.F(y)
if y^2 + self.a*x^2 != 1 + self.d*x^2*y^2:
raise NotOnCurveException(str(self))
def __repr__(self):
return "%s(0x%x,0x%x)" % (self.__class__.__name__, self.x, self.y)
def __iter__(self):
yield self.x
yield self.y
def __add__(self,other):
x,y = self
X,Y = other
a,d = self.a,self.d
return self.__class__(
(x*Y+y*X)/(1+d*x*y*X*Y),
(y*Y-a*x*X)/(1-d*x*y*X*Y)
)
def __neg__(self): return self.__class__(-self.x,self.y)
def __sub__(self,other): return self + (-other)
def __rmul__(self,other): return self*other
def __eq__(self,other):
"""NB: this is the only method that is different from the usual one"""
x,y = self
X,Y = other
return x*Y == X*y or (self.cofactor==8 and -self.a*x*X == y*Y)
def __ne__(self,other): return not (self==other)
def __mul__(self,exp):
exp = int(exp)
if exp < 0: exp,self = -exp,-self
total = self.__class__()
work = self
while exp != 0:
if exp & 1: total += work
work += work
exp >>= 1
return total
def xyzt(self):
x,y = self
z = self.F.random_element()
return x*z,y*z,z,x*y*z
def torque(self):
"""Apply cofactor group, except keeping the point even"""
if self.cofactor == 8:
if self.a == -1: return self.__class__(self.y*self.i, self.x*self.i)
if self.a == 1: return self.__class__(-self.y, self.x)
else:
return self.__class__(-self.x, -self.y)
# Utility functions
@classmethod
def bytesToGf(cls,bytes,mustBeProper=True,mustBePositive=False):
"""Convert little-endian bytes to field element, sanity check length"""
if len(bytes) != cls.encLen:
raise InvalidEncodingException("wrong length %d" % len(bytes))
s = dec_le(bytes)
if mustBeProper and s >= cls.F.modulus():
raise InvalidEncodingException("%d out of range!" % s)
s = cls.F(s)
if mustBePositive and negative(s):
raise InvalidEncodingException("%d is negative!" % s)
return s
@classmethod
def gfToBytes(cls,x,mustBePositive=False):
"""Convert little-endian bytes to field element, sanity check length"""
if negative(x) and mustBePositive: x = -x
return enc_le(x,cls.encLen)
class RistrettoPoint(QuotientEdwardsPoint):
"""The new Ristretto group"""
def encodeSpec(self):
"""Unoptimized specification for encoding"""
x,y = self
if self.cofactor==8 and (negative(x*y) or y==0): (x,y) = self.torque()
if y == -1: y = 1 # Avoid divide by 0; doesn't affect impl
if negative(x): x,y = -x,-y
s = xsqrt(self.mneg*(1-y)/(1+y),exn=Exception("Unimplemented: point is odd: " + str(self)))
return self.gfToBytes(s)
@classmethod
def decodeSpec(cls,s):
"""Unoptimized specification for decoding"""
s = cls.bytesToGf(s,mustBePositive=True)
a,d = cls.a,cls.d
x = xsqrt(4*s^2 / (a*d*(1+a*s^2)^2 - (1-a*s^2)^2))
y = (1+a*s^2) / (1-a*s^2)
if cls.cofactor==8 and (negative(x*y) or y==0):
raise InvalidEncodingException("x*y has high bit")
return cls(x,y)
@optimized_version_of("encodeSpec")
def encode(self):
"""Encode, optimized version"""
a,d,mneg = self.a,self.d,self.mneg
x,y,z,t = self.xyzt()
if self.cofactor==8:
u1 = mneg*(z+y)*(z-y)
u2 = x*y # = t*z
isr = isqrt(u1*u2^2)
i1 = isr*u1 # sqrt(mneg*(z+y)*(z-y))/(x*y)
i2 = isr*u2 # 1/sqrt(a*(y+z)*(y-z))
z_inv = i1*i2*t # 1/z
if negative(t*z_inv):
if a==-1:
x,y = y*self.i,x*self.i
den_inv = self.magic * i1
else:
x,y = -y,x
den_inv = self.i * self.magic * i1
else:
den_inv = i2
if negative(x*z_inv): y = -y
s = (z-y) * den_inv
else:
num = mneg*(z+y)*(z-y)
isr = isqrt(num*y^2)
if negative(isr^2*num*y*t): y = -y
s = isr*y*(z-y)
return self.gfToBytes(s,mustBePositive=True)
@classmethod
@optimized_version_of("decodeSpec")
def decode(cls,s):
"""Decode, optimized version"""
s = cls.bytesToGf(s,mustBePositive=True)
a,d = cls.a,cls.d
yden = 1-a*s^2
ynum = 1+a*s^2
yden_sqr = yden^2
xden_sqr = a*d*ynum^2 - yden_sqr
isr = isqrt(xden_sqr * yden_sqr)
xden_inv = isr * yden
yden_inv = xden_inv * isr * xden_sqr
x = 2*s*xden_inv
if negative(x): x = -x
y = ynum * yden_inv
if cls.cofactor==8 and (negative(x*y) or y==0):
raise InvalidEncodingException("x*y is invalid: %d, %d" % (x,y))
return cls(x,y)
@classmethod
def fromJacobiQuartic(cls,s,t,sgn=1):
"""Convert point from its Jacobi Quartic representation"""
a,d = cls.a,cls.d
assert s^4 - 2*cls.a*(1-2*d/(d-a))*s^2 + 1 == t^2
x = 2*s*cls.magic / t
y = (1+a*s^2) / (1-a*s^2)
return cls(sgn*x,y)
@classmethod
def elligatorSpec(cls,r0):
a,d = cls.a,cls.d
r = cls.qnr * cls.bytesToGf(r0)^2
den = (d*r-a)*(a*r-d)
n1 = cls.a*(r+1)*(a+d)*(d-a)/den
n2 = r*n1
if is_square(n1):
sgn,s,t = 1, xsqrt(n1), -(r-1)*(a+d)^2 / den - 1
else:
sgn,s,t = -1,-xsqrt(n2), r*(r-1)*(a+d)^2 / den - 1
return cls.fromJacobiQuartic(s,t)
@classmethod
@optimized_version_of("elligatorSpec")
def elligator(cls,r0):
a,d = cls.a,cls.d
r0 = cls.bytesToGf(r0)
r = cls.qnr * r0^2
den = (d*r-a)*(a*r-d)
num = cls.a*(r+1)*(a+d)*(d-a)
iss,isri = isqrt_i(num*den)
if iss: sgn,twiddle = 1,1
else: sgn,twiddle = -1,r0*cls.qnr
isri *= twiddle
s = isri*num
t = -sgn*isri*s*(r-1)*(d+a)^2 - 1
if negative(s) == iss: s = -s
return cls.fromJacobiQuartic(s,t)
class Decaf_1_1_Point(QuotientEdwardsPoint):
"""Like current decaf but tweaked for simplicity"""
def encodeSpec(self):
"""Unoptimized specification for encoding"""
a,d = self.a,self.d
x,y = self
if x==0 or y==0: return(self.gfToBytes(0))
if self.cofactor==8 and negative(x*y*self.isoMagic):
x,y = self.torque()
isr2 = isqrt(a*(y^2-1)) * sqrt(a*d-1)
sr = xsqrt(1-a*x^2)
assert sr in [isr2*x*y,-isr2*x*y]
altx = 1/isr2*self.isoMagic
if negative(altx): s = (1+x*y*isr2)/(a*x)
else: s = (1-x*y*isr2)/(a*x)
return self.gfToBytes(s,mustBePositive=True)
@classmethod
def decodeSpec(cls,s):
"""Unoptimized specification for decoding"""
a,d = cls.a,cls.d
s = cls.bytesToGf(s,mustBePositive=True)
if s==0: return cls()
isr = isqrt(s^4 + 2*(a-2*d)*s^2 + 1)
altx = 2*s*isr*cls.isoMagic
if negative(altx): isr = -isr
x = 2*s / (1+a*s^2)
y = (1-a*s^2) * isr
if cls.cofactor==8 and (negative(x*y*cls.isoMagic) or y==0):
raise InvalidEncodingException("x*y is invalid: %d, %d" % (x,y))
return cls(x,y)
@optimized_version_of("encodeSpec")
def encode(self):
"""Encode, optimized version"""
a,d = self.a,self.d
x,y,z,t = self.xyzt()
if self.cofactor == 8:
# Cofactor 8 version
num = (z+y)*(z-y)
den = x*y
tmp = isqrt(num*(a-d)*den^2)
if negative(tmp^2*den*num*(a-d)*t^2*self.isoMagic):
den,num = num,den
tmp *= sqrt(a-d) # witness that cofactor is 8
yisr = x*sqrt(a)
toggle = (a==1)
else:
yisr = y*(a*d-1)
toggle = False
tiisr = tmp*num
altx = tiisr*t*self.isoMagic
if negative(altx) != toggle: tiisr =- tiisr
s = tmp*den*yisr*(tiisr*z - 1)
else:
# Much simpler cofactor 4 version
num = (x+t)*(x-t)
isr = isqrt(num*(a-d)*x^2)
ratio = isr*num
if negative(ratio*self.isoMagic): ratio=-ratio
s = (a-d)*isr*x*(ratio*z - t)
return self.gfToBytes(s,mustBePositive=True)
@classmethod
@optimized_version_of("decodeSpec")
def decode(cls,s):
"""Decode, optimized version"""
a,d = cls.a,cls.d
s = cls.bytesToGf(s,mustBePositive=True)
if s==0: return cls()
s2 = s^2
den = 1+a*s2
num = den^2 - 4*d*s2
isr = isqrt(num*den^2)
altx = 2*s*isr*den*cls.isoMagic
if negative(altx): isr = -isr
x = 2*s *isr^2*den*num
y = (1-a*s^2) * isr*den
if cls.cofactor==8 and (negative(x*y*cls.isoMagic) or y==0):
raise InvalidEncodingException("x*y is invalid: %d, %d" % (x,y))
return cls(x,y)
@classmethod
def fromJacobiQuartic(cls,s,t,sgn=1):
"""Convert point from its Jacobi Quartic representation"""
a,d = cls.a,cls.d
if s==0: return cls()
x = 2*s / (1+a*s^2)
y = (1-a*s^2) / t
return cls(x,sgn*y)
@classmethod
def elligatorSpec(cls,r0):
a,d = cls.a,cls.d
r = cls.qnr * cls.bytesToGf(r0)^2
den = (d*r-(d-a))*((d-a)*r-d)
n1 = (r+1)*(a-2*d)/den
n2 = r*n1
if is_square(n1):
sgn,s,t = 1, xsqrt(n1), -(r-1)*(a-2*d)^2 / den - 1
else:
sgn,s,t = -1, -xsqrt(n2), r*(r-1)*(a-2*d)^2 / den - 1
return cls.fromJacobiQuartic(s,t)
@classmethod
@optimized_version_of("elligatorSpec")
def elligator(cls,r0):
a,d = cls.a,cls.d
r0 = cls.bytesToGf(r0)
r = cls.qnr * r0^2
den = (d*r-(d-a))*((d-a)*r-d)
num = (r+1)*(a-2*d)
iss,isri = isqrt_i(num*den)
if iss: sgn,twiddle = 1,1
else: sgn,twiddle = -1,r0*cls.qnr
isri *= twiddle
s = isri*num
t = -sgn*isri*s*(r-1)*(a-2*d)^2 - 1
if negative(s) == iss: s = -s
return cls.fromJacobiQuartic(s,t)
class Ed25519Point(RistrettoPoint):
F = GF(2^255-19)
d = F(-121665/121666)
a = F(-1)
i = sqrt(F(-1))
mneg = F(1)
qnr = i
magic = isqrt(a*d-1)
cofactor = 8
encLen = 32
@classmethod
def base(cls):
return cls( 15112221349535400772501151409588531511454012693041857206046113283949847762202, 46316835694926478169428394003475163141307993866256225615783033603165251855960
)
class NegEd25519Point(RistrettoPoint):
F = GF(2^255-19)
d = F(121665/121666)
a = F(1)
i = sqrt(F(-1))
mneg = F(-1) # TODO checkme vs 1-ad or whatever
qnr = i
magic = isqrt(a*d-1)
cofactor = 8
encLen = 32
@classmethod
def base(cls):
y = cls.F(4/5)
x = sqrt((y^2-1)/(cls.d*y^2-cls.a))
if negative(x): x = -x
return cls(x,y)
class IsoEd448Point(RistrettoPoint):
F = GF(2^448-2^224-1)
d = F(39082/39081)
a = F(1)
mneg = F(-1)
qnr = -1
magic = isqrt(a*d-1)
cofactor = 4
encLen = 56
@classmethod
def base(cls):
return cls( # RFC has it wrong
-345397493039729516374008604150537410266655260075183290216406970281645695073672344430481787759340633221708391583424041788924124567700732,
-363419362147803445274661903944002267176820680343659030140745099590306164083365386343198191849338272965044442230921818680526749009182718
)
class TwistedEd448GoldilocksPoint(Decaf_1_1_Point):
F = GF(2^448-2^224-1)
d = F(-39082)
a = F(-1)
qnr = -1
magic = isqrt(a*d-1)
cofactor = 4
encLen = 56
isoMagic = IsoEd448Point.magic
@classmethod
def base(cls):
return cls.decodeSpec(Ed448GoldilocksPoint.base().encodeSpec())
class Ed448GoldilocksPoint(Decaf_1_1_Point):
F = GF(2^448-2^224-1)
d = F(-39081)
a = F(1)
qnr = -1
magic = isqrt(a*d-1)
cofactor = 4
encLen = 56
isoMagic = IsoEd448Point.magic
@classmethod
def base(cls):
return -2*cls( # FIXME: make not negative
224580040295924300187604334099896036246789641632564134246125461686950415467406032909029192869357953282578032075146446173674602635247710, 298819210078481492676017930443930673437544040154080242095928241372331506189835876003536878655418784733982303233503462500531545062832660
)
class IsoEd25519Point(Decaf_1_1_Point):
# TODO: twisted iso too!
# TODO: twisted iso might have to IMAGINE_TWIST or whatever
F = GF(2^255-19)
d = F(-121665)
a = F(1)
i = sqrt(F(-1))
qnr = i
magic = isqrt(a*d-1)
cofactor = 8
encLen = 32
isoMagic = Ed25519Point.magic
isoA = Ed25519Point.a
@classmethod
def base(cls):
return cls.decodeSpec(Ed25519Point.base().encode())
class TestFailedException(Exception): pass
def test(cls,n):
print "Testing curve %s" % cls.__name__
specials = [1]
ii = cls.F(-1)
while is_square(ii):
specials.append(ii)
ii = sqrt(ii)
specials.append(ii)
for i in specials:
if negative(cls.F(i)): i = -i
i = enc_le(i,cls.encLen)
try:
Q = cls.decode(i)
QE = Q.encode()
if QE != i:
raise TestFailedException("Round trip special %s != %s" %
(binascii.hexlify(QE),binascii.hexlify(i)))
except NotOnCurveException: pass
except InvalidEncodingException: pass
P = cls.base()
print "base", list(P.encode())
for i in xrange(16):
Q = P*i
print i, list(Q.encode())
Q = cls()
for i in xrange(n):
#print i, binascii.hexlify(Q.encode())
QQ = cls.decode(Q.encode())
if QQ != Q: raise TestFailedException("Round trip %s != %s" % (str(QQ),str(Q)))
QT = Q
QE = Q.encode()
for h in xrange(cls.cofactor):
QT = QT.torque()
if QT.encode() != QE:
raise TestFailedException("Can't torque %s,%d" % (str(Q),h+1))
Q0 = Q + P
if Q0 == Q: raise TestFailedException("Addition doesn't work")
if Q0-P != Q: raise TestFailedException("Subtraction doesn't work")
r = randint(1,1000)
Q1 = Q0*r
Q2 = Q0*(r+1)
if Q1 + Q0 != Q2: raise TestFailedException("Scalarmul doesn't work")
Q = Q1
test(Ed25519Point,100)
#test(NegEd25519Point,100)
#test(IsoEd25519Point,100)
#test(IsoEd448Point,100)
#test(TwistedEd448GoldilocksPoint,100)
#test(Ed448GoldilocksPoint,100)
def testElligator(cls,n):
print "Testing elligator on %s" % cls.__name__
for i in xrange(n):
r = randombytes(cls.encLen)
Q = cls.elligator(r)
print list(r), list(Q.encode())
testElligator(Ed25519Point,100)
#testElligator(NegEd25519Point,100)
#testElligator(IsoEd448Point,100)
#testElligator(Ed448GoldilocksPoint,100)
#testElligator(TwistedEd448GoldilocksPoint,100)
def gangtest(classes,n):
print "Gang test",[cls.__name__ for cls in classes]
specials = [1]
ii = classes[0].F(-1)
while is_square(ii):
specials.append(ii)
ii = sqrt(ii)
specials.append(ii)
for i in xrange(n):
rets = [bytes((cls.base()*i).encode()) for cls in classes]
if len(set(rets)) != 1:
print "Divergence in encode at %d" % i
for c,ret in zip(classes,rets):
print c,binascii.hexlify(ret)
print
if i < len(specials): r0 = enc_le(specials[i],classes[0].encLen)
else: r0 = randombytes(classes[0].encLen)
rets = [bytes((cls.elligator(r0)*i).encode()) for cls in classes]
if len(set(rets)) != 1:
print "Divergence in elligator at %d" % i
for c,ret in zip(classes,rets):
print c,binascii.hexlify(ret)
print
gangtest([IsoEd448Point,TwistedEd448GoldilocksPoint,Ed448GoldilocksPoint],100)
gangtest([Ed25519Point,IsoEd25519Point],100)