arnaucube
arnaucube/testudoPublic Self-hosted repository
testudo/src/scalar/ristretto255.rs
src/scalar/ristretto255.rsRaw
   1//! This module provides an implementation of the Curve25519's scalar field $\mathbb{F}_q$
   2//! where `q = 2^252 + 27742317777372353535851937790883648493 = 0x1000000000000000 0000000000000000 14def9dea2f79cd6 5812631a5cf5d3ed`
   3//! The entire file is an adaptation from bls12-381 crate. We modify various constants (MODULUS, R, R2, etc.) to appropriate values for Curve25519 and update tests
   4//! We borrow the `invert` method from curve25519-dalek crate
   5//! See NOTICE.md for more details
   6#![allow(clippy::all)]
   7use core::borrow::Borrow;
   8use core::convert::TryFrom;
   9use core::fmt;
  10use core::iter::{Product, Sum};
  11use core::ops::{Add, AddAssign, Mul, MulAssign, Neg, Sub, SubAssign};
  12use rand_core::{CryptoRng, RngCore};
  13use serde::{Deserialize, Serialize};
  14use subtle::{Choice, ConditionallySelectable, ConstantTimeEq, CtOption};
  15use zeroize::Zeroize;
  16
  17// use crate::util::{adc, mac, sbb};
  18/// Compute a + b + carry, returning the result and the new carry over.
  19#[inline(always)]
  20pub const fn adc(a: u64, b: u64, carry: u64) -> (u64, u64) {
  21  let ret = (a as u128) + (b as u128) + (carry as u128);
  22  (ret as u64, (ret >> 64) as u64)
  23}
  24
  25/// Compute a - (b + borrow), returning the result and the new borrow.
  26#[inline(always)]
  27pub const fn sbb(a: u64, b: u64, borrow: u64) -> (u64, u64) {
  28  let ret = (a as u128).wrapping_sub((b as u128) + ((borrow >> 63) as u128));
  29  (ret as u64, (ret >> 64) as u64)
  30}
  31
  32/// Compute a + (b * c) + carry, returning the result and the new carry over.
  33#[inline(always)]
  34pub const fn mac(a: u64, b: u64, c: u64, carry: u64) -> (u64, u64) {
  35  let ret = (a as u128) + ((b as u128) * (c as u128)) + (carry as u128);
  36  (ret as u64, (ret >> 64) as u64)
  37}
  38
  39macro_rules! impl_add_binop_specify_output {
  40  ($lhs:ident, $rhs:ident, $output:ident) => {
  41    impl<'b> Add<&'b $rhs> for $lhs {
  42      type Output = $output;
  43
  44      #[inline]
  45      fn add(self, rhs: &'b $rhs) -> $output {
  46        &self + rhs
  47      }
  48    }
  49
  50    impl<'a> Add<$rhs> for &'a $lhs {
  51      type Output = $output;
  52
  53      #[inline]
  54      fn add(self, rhs: $rhs) -> $output {
  55        self + &rhs
  56      }
  57    }
  58
  59    impl Add<$rhs> for $lhs {
  60      type Output = $output;
  61
  62      #[inline]
  63      fn add(self, rhs: $rhs) -> $output {
  64        &self + &rhs
  65      }
  66    }
  67  };
  68}
  69
  70macro_rules! impl_sub_binop_specify_output {
  71  ($lhs:ident, $rhs:ident, $output:ident) => {
  72    impl<'b> Sub<&'b $rhs> for $lhs {
  73      type Output = $output;
  74
  75      #[inline]
  76      fn sub(self, rhs: &'b $rhs) -> $output {
  77        &self - rhs
  78      }
  79    }
  80
  81    impl<'a> Sub<$rhs> for &'a $lhs {
  82      type Output = $output;
  83
  84      #[inline]
  85      fn sub(self, rhs: $rhs) -> $output {
  86        self - &rhs
  87      }
  88    }
  89
  90    impl Sub<$rhs> for $lhs {
  91      type Output = $output;
  92
  93      #[inline]
  94      fn sub(self, rhs: $rhs) -> $output {
  95        &self - &rhs
  96      }
  97    }
  98  };
  99}
 100
 101macro_rules! impl_binops_additive_specify_output {
 102  ($lhs:ident, $rhs:ident, $output:ident) => {
 103    impl_add_binop_specify_output!($lhs, $rhs, $output);
 104    impl_sub_binop_specify_output!($lhs, $rhs, $output);
 105  };
 106}
 107
 108macro_rules! impl_binops_multiplicative_mixed {
 109  ($lhs:ident, $rhs:ident, $output:ident) => {
 110    impl<'b> Mul<&'b $rhs> for $lhs {
 111      type Output = $output;
 112
 113      #[inline]
 114      fn mul(self, rhs: &'b $rhs) -> $output {
 115        &self * rhs
 116      }
 117    }
 118
 119    impl<'a> Mul<$rhs> for &'a $lhs {
 120      type Output = $output;
 121
 122      #[inline]
 123      fn mul(self, rhs: $rhs) -> $output {
 124        self * &rhs
 125      }
 126    }
 127
 128    impl Mul<$rhs> for $lhs {
 129      type Output = $output;
 130
 131      #[inline]
 132      fn mul(self, rhs: $rhs) -> $output {
 133        &self * &rhs
 134      }
 135    }
 136  };
 137}
 138
 139macro_rules! impl_binops_additive {
 140  ($lhs:ident, $rhs:ident) => {
 141    impl_binops_additive_specify_output!($lhs, $rhs, $lhs);
 142
 143    impl SubAssign<$rhs> for $lhs {
 144      #[inline]
 145      fn sub_assign(&mut self, rhs: $rhs) {
 146        *self = &*self - &rhs;
 147      }
 148    }
 149
 150    impl AddAssign<$rhs> for $lhs {
 151      #[inline]
 152      fn add_assign(&mut self, rhs: $rhs) {
 153        *self = &*self + &rhs;
 154      }
 155    }
 156
 157    impl<'b> SubAssign<&'b $rhs> for $lhs {
 158      #[inline]
 159      fn sub_assign(&mut self, rhs: &'b $rhs) {
 160        *self = &*self - rhs;
 161      }
 162    }
 163
 164    impl<'b> AddAssign<&'b $rhs> for $lhs {
 165      #[inline]
 166      fn add_assign(&mut self, rhs: &'b $rhs) {
 167        *self = &*self + rhs;
 168      }
 169    }
 170  };
 171}
 172
 173macro_rules! impl_binops_multiplicative {
 174  ($lhs:ident, $rhs:ident) => {
 175    impl_binops_multiplicative_mixed!($lhs, $rhs, $lhs);
 176
 177    impl MulAssign<$rhs> for $lhs {
 178      #[inline]
 179      fn mul_assign(&mut self, rhs: $rhs) {
 180        *self = &*self * &rhs;
 181      }
 182    }
 183
 184    impl<'b> MulAssign<&'b $rhs> for $lhs {
 185      #[inline]
 186      fn mul_assign(&mut self, rhs: &'b $rhs) {
 187        *self = &*self * rhs;
 188      }
 189    }
 190  };
 191}
 192
 193/// Represents an element of the scalar field $\mathbb{F}_q$ of the Curve25519 elliptic
 194/// curve construction.
 195// The internal representation of this type is four 64-bit unsigned
 196// integers in little-endian order. `Scalar` values are always in
 197// Montgomery form; i.e., Scalar(a) = aR mod q, with R = 2^256.
 198#[derive(Clone, Copy, Eq, Serialize, Deserialize)]
 199pub struct Scalar(pub(crate) [u64; 4]);
 200
 201impl fmt::Debug for Scalar {
 202  fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result {
 203    let tmp = self.to_bytes();
 204    write!(f, "0x")?;
 205    for &b in tmp.iter().rev() {
 206      write!(f, "{:02x}", b)?;
 207    }
 208    Ok(())
 209  }
 210}
 211
 212impl From<u64> for Scalar {
 213  fn from(val: u64) -> Scalar {
 214    Scalar([val, 0, 0, 0]) * R2
 215  }
 216}
 217
 218impl ConstantTimeEq for Scalar {
 219  fn ct_eq(&self, other: &Self) -> Choice {
 220    self.0[0].ct_eq(&other.0[0])
 221      & self.0[1].ct_eq(&other.0[1])
 222      & self.0[2].ct_eq(&other.0[2])
 223      & self.0[3].ct_eq(&other.0[3])
 224  }
 225}
 226
 227impl PartialEq for Scalar {
 228  #[inline]
 229  fn eq(&self, other: &Self) -> bool {
 230    self.ct_eq(other).unwrap_u8() == 1
 231  }
 232}
 233
 234impl ConditionallySelectable for Scalar {
 235  fn conditional_select(a: &Self, b: &Self, choice: Choice) -> Self {
 236    Scalar([
 237      u64::conditional_select(&a.0[0], &b.0[0], choice),
 238      u64::conditional_select(&a.0[1], &b.0[1], choice),
 239      u64::conditional_select(&a.0[2], &b.0[2], choice),
 240      u64::conditional_select(&a.0[3], &b.0[3], choice),
 241    ])
 242  }
 243}
 244
 245/// Constant representing the modulus
 246/// q = 2^252 + 27742317777372353535851937790883648493
 247/// 0x1000000000000000 0000000000000000 14def9dea2f79cd6 5812631a5cf5d3ed
 248const MODULUS: Scalar = Scalar([
 249  0x5812_631a_5cf5_d3ed,
 250  0x14de_f9de_a2f7_9cd6,
 251  0x0000_0000_0000_0000,
 252  0x1000_0000_0000_0000,
 253]);
 254
 255impl<'a> Neg for &'a Scalar {
 256  type Output = Scalar;
 257
 258  #[inline]
 259  fn neg(self) -> Scalar {
 260    self.neg()
 261  }
 262}
 263
 264impl Neg for Scalar {
 265  type Output = Scalar;
 266
 267  #[inline]
 268  fn neg(self) -> Scalar {
 269    -&self
 270  }
 271}
 272
 273impl<'a, 'b> Sub<&'b Scalar> for &'a Scalar {
 274  type Output = Scalar;
 275
 276  #[inline]
 277  fn sub(self, rhs: &'b Scalar) -> Scalar {
 278    self.sub(rhs)
 279  }
 280}
 281
 282impl<'a, 'b> Add<&'b Scalar> for &'a Scalar {
 283  type Output = Scalar;
 284
 285  #[inline]
 286  fn add(self, rhs: &'b Scalar) -> Scalar {
 287    self.add(rhs)
 288  }
 289}
 290
 291impl<'a, 'b> Mul<&'b Scalar> for &'a Scalar {
 292  type Output = Scalar;
 293
 294  #[inline]
 295  fn mul(self, rhs: &'b Scalar) -> Scalar {
 296    self.mul(rhs)
 297  }
 298}
 299
 300impl_binops_additive!(Scalar, Scalar);
 301impl_binops_multiplicative!(Scalar, Scalar);
 302
 303/// INV = -(q^{-1} mod 2^64) mod 2^64
 304const INV: u64 = 0xd2b5_1da3_1254_7e1b;
 305
 306/// R = 2^256 mod q
 307const R: Scalar = Scalar([
 308  0xd6ec_3174_8d98_951d,
 309  0xc6ef_5bf4_737d_cf70,
 310  0xffff_ffff_ffff_fffe,
 311  0x0fff_ffff_ffff_ffff,
 312]);
 313
 314/// R^2 = 2^512 mod q
 315const R2: Scalar = Scalar([
 316  0xa406_11e3_449c_0f01,
 317  0xd00e_1ba7_6885_9347,
 318  0xceec_73d2_17f5_be65,
 319  0x0399_411b_7c30_9a3d,
 320]);
 321
 322/// R^3 = 2^768 mod q
 323const R3: Scalar = Scalar([
 324  0x2a9e_4968_7b83_a2db,
 325  0x2783_24e6_aef7_f3ec,
 326  0x8065_dc6c_04ec_5b65,
 327  0x0e53_0b77_3599_cec7,
 328]);
 329
 330impl Default for Scalar {
 331  #[inline]
 332  fn default() -> Self {
 333    Self::zero()
 334  }
 335}
 336
 337impl<T> Product<T> for Scalar
 338where
 339  T: Borrow<Scalar>,
 340{
 341  fn product<I>(iter: I) -> Self
 342  where
 343    I: Iterator<Item = T>,
 344  {
 345    iter.fold(Scalar::one(), |acc, item| acc * item.borrow())
 346  }
 347}
 348
 349impl<T> Sum<T> for Scalar
 350where
 351  T: Borrow<Scalar>,
 352{
 353  fn sum<I>(iter: I) -> Self
 354  where
 355    I: Iterator<Item = T>,
 356  {
 357    iter.fold(Scalar::zero(), |acc, item| acc + item.borrow())
 358  }
 359}
 360
 361impl Zeroize for Scalar {
 362  fn zeroize(&mut self) {
 363    self.0 = [0u64; 4];
 364  }
 365}
 366
 367impl Scalar {
 368  /// Returns zero, the additive identity.
 369  #[inline]
 370  pub const fn zero() -> Scalar {
 371    Scalar([0, 0, 0, 0])
 372  }
 373
 374  /// Returns one, the multiplicative identity.
 375  #[inline]
 376  pub const fn one() -> Scalar {
 377    R
 378  }
 379
 380  pub fn random<Rng: RngCore + CryptoRng>(rng: &mut Rng) -> Self {
 381    let mut limbs = [0u64; 8];
 382    for i in 0..8 {
 383      limbs[i] = rng.next_u64();
 384    }
 385    Scalar::from_u512(limbs)
 386  }
 387
 388  /// Doubles this field element.
 389  #[inline]
 390  pub const fn double(&self) -> Scalar {
 391    // TODO: This can be achieved more efficiently with a bitshift.
 392    self.add(self)
 393  }
 394
 395  /// Attempts to convert a little-endian byte representation of
 396  /// a scalar into a `Scalar`, failing if the input is not canonical.
 397  pub fn from_bytes(bytes: &[u8; 32]) -> CtOption<Scalar> {
 398    let mut tmp = Scalar([0, 0, 0, 0]);
 399
 400    tmp.0[0] = u64::from_le_bytes(<[u8; 8]>::try_from(&bytes[0..8]).unwrap());
 401    tmp.0[1] = u64::from_le_bytes(<[u8; 8]>::try_from(&bytes[8..16]).unwrap());
 402    tmp.0[2] = u64::from_le_bytes(<[u8; 8]>::try_from(&bytes[16..24]).unwrap());
 403    tmp.0[3] = u64::from_le_bytes(<[u8; 8]>::try_from(&bytes[24..32]).unwrap());
 404
 405    // Try to subtract the modulus
 406    let (_, borrow) = sbb(tmp.0[0], MODULUS.0[0], 0);
 407    let (_, borrow) = sbb(tmp.0[1], MODULUS.0[1], borrow);
 408    let (_, borrow) = sbb(tmp.0[2], MODULUS.0[2], borrow);
 409    let (_, borrow) = sbb(tmp.0[3], MODULUS.0[3], borrow);
 410
 411    // If the element is smaller than MODULUS then the
 412    // subtraction will underflow, producing a borrow value
 413    // of 0xffff...ffff. Otherwise, it'll be zero.
 414    let is_some = (borrow as u8) & 1;
 415
 416    // Convert to Montgomery form by computing
 417    // (a.R^0 * R^2) / R = a.R
 418    tmp *= &R2;
 419
 420    CtOption::new(tmp, Choice::from(is_some))
 421  }
 422
 423  /// Converts an element of `Scalar` into a byte representation in
 424  /// little-endian byte order.
 425  pub fn to_bytes(&self) -> [u8; 32] {
 426    // Turn into canonical form by computing
 427    // (a.R) / R = a
 428    let tmp = Scalar::montgomery_reduce(self.0[0], self.0[1], self.0[2], self.0[3], 0, 0, 0, 0);
 429
 430    let mut res = [0; 32];
 431    res[0..8].copy_from_slice(&tmp.0[0].to_le_bytes());
 432    res[8..16].copy_from_slice(&tmp.0[1].to_le_bytes());
 433    res[16..24].copy_from_slice(&tmp.0[2].to_le_bytes());
 434    res[24..32].copy_from_slice(&tmp.0[3].to_le_bytes());
 435
 436    res
 437  }
 438
 439  /// Converts a 512-bit little endian integer into
 440  /// a `Scalar` by reducing by the modulus.
 441  pub fn from_bytes_wide(bytes: &[u8; 64]) -> Scalar {
 442    Scalar::from_u512([
 443      u64::from_le_bytes(<[u8; 8]>::try_from(&bytes[0..8]).unwrap()),
 444      u64::from_le_bytes(<[u8; 8]>::try_from(&bytes[8..16]).unwrap()),
 445      u64::from_le_bytes(<[u8; 8]>::try_from(&bytes[16..24]).unwrap()),
 446      u64::from_le_bytes(<[u8; 8]>::try_from(&bytes[24..32]).unwrap()),
 447      u64::from_le_bytes(<[u8; 8]>::try_from(&bytes[32..40]).unwrap()),
 448      u64::from_le_bytes(<[u8; 8]>::try_from(&bytes[40..48]).unwrap()),
 449      u64::from_le_bytes(<[u8; 8]>::try_from(&bytes[48..56]).unwrap()),
 450      u64::from_le_bytes(<[u8; 8]>::try_from(&bytes[56..64]).unwrap()),
 451    ])
 452  }
 453
 454  fn from_u512(limbs: [u64; 8]) -> Scalar {
 455    // We reduce an arbitrary 512-bit number by decomposing it into two 256-bit digits
 456    // with the higher bits multiplied by 2^256. Thus, we perform two reductions
 457    //
 458    // 1. the lower bits are multiplied by R^2, as normal
 459    // 2. the upper bits are multiplied by R^2 * 2^256 = R^3
 460    //
 461    // and computing their sum in the field. It remains to see that arbitrary 256-bit
 462    // numbers can be placed into Montgomery form safely using the reduction. The
 463    // reduction works so long as the product is less than R=2^256 multipled by
 464    // the modulus. This holds because for any `c` smaller than the modulus, we have
 465    // that (2^256 - 1)*c is an acceptable product for the reduction. Therefore, the
 466    // reduction always works so long as `c` is in the field; in this case it is either the
 467    // constant `R2` or `R3`.
 468    let d0 = Scalar([limbs[0], limbs[1], limbs[2], limbs[3]]);
 469    let d1 = Scalar([limbs[4], limbs[5], limbs[6], limbs[7]]);
 470    // Convert to Montgomery form
 471    d0 * R2 + d1 * R3
 472  }
 473
 474  /// Converts from an integer represented in little endian
 475  /// into its (congruent) `Scalar` representation.
 476  pub const fn from_raw(val: [u64; 4]) -> Self {
 477    (&Scalar(val)).mul(&R2)
 478  }
 479
 480  /// Squares this element.
 481  #[inline]
 482  pub const fn square(&self) -> Scalar {
 483    let (r1, carry) = mac(0, self.0[0], self.0[1], 0);
 484    let (r2, carry) = mac(0, self.0[0], self.0[2], carry);
 485    let (r3, r4) = mac(0, self.0[0], self.0[3], carry);
 486
 487    let (r3, carry) = mac(r3, self.0[1], self.0[2], 0);
 488    let (r4, r5) = mac(r4, self.0[1], self.0[3], carry);
 489
 490    let (r5, r6) = mac(r5, self.0[2], self.0[3], 0);
 491
 492    let r7 = r6 >> 63;
 493    let r6 = (r6 << 1) | (r5 >> 63);
 494    let r5 = (r5 << 1) | (r4 >> 63);
 495    let r4 = (r4 << 1) | (r3 >> 63);
 496    let r3 = (r3 << 1) | (r2 >> 63);
 497    let r2 = (r2 << 1) | (r1 >> 63);
 498    let r1 = r1 << 1;
 499
 500    let (r0, carry) = mac(0, self.0[0], self.0[0], 0);
 501    let (r1, carry) = adc(0, r1, carry);
 502    let (r2, carry) = mac(r2, self.0[1], self.0[1], carry);
 503    let (r3, carry) = adc(0, r3, carry);
 504    let (r4, carry) = mac(r4, self.0[2], self.0[2], carry);
 505    let (r5, carry) = adc(0, r5, carry);
 506    let (r6, carry) = mac(r6, self.0[3], self.0[3], carry);
 507    let (r7, _) = adc(0, r7, carry);
 508
 509    Scalar::montgomery_reduce(r0, r1, r2, r3, r4, r5, r6, r7)
 510  }
 511
 512  /// Exponentiates `self` by `by`, where `by` is a
 513  /// little-endian order integer exponent.
 514  pub fn pow(&self, by: &[u64; 4]) -> Self {
 515    let mut res = Self::one();
 516    for e in by.iter().rev() {
 517      for i in (0..64).rev() {
 518        res = res.square();
 519        let mut tmp = res;
 520        tmp *= self;
 521        res.conditional_assign(&tmp, (((*e >> i) & 0x1) as u8).into());
 522      }
 523    }
 524    res
 525  }
 526
 527  /// Exponentiates `self` by `by`, where `by` is a
 528  /// little-endian order integer exponent.
 529  ///
 530  /// **This operation is variable time with respect
 531  /// to the exponent.** If the exponent is fixed,
 532  /// this operation is effectively constant time.
 533  pub fn pow_vartime(&self, by: &[u64; 4]) -> Self {
 534    let mut res = Self::one();
 535    for e in by.iter().rev() {
 536      for i in (0..64).rev() {
 537        res = res.square();
 538
 539        if ((*e >> i) & 1) == 1 {
 540          res.mul_assign(self);
 541        }
 542      }
 543    }
 544    res
 545  }
 546
 547  pub fn invert(&self) -> CtOption<Self> {
 548    // Uses the addition chain from
 549    // https://briansmith.org/ecc-inversion-addition-chains-01#curve25519_scalar_inversion
 550    // implementation adapted from curve25519-dalek
 551    let _1 = self;
 552    let _10 = _1.square();
 553    let _100 = _10.square();
 554    let _11 = &_10 * _1;
 555    let _101 = &_10 * &_11;
 556    let _111 = &_10 * &_101;
 557    let _1001 = &_10 * &_111;
 558    let _1011 = &_10 * &_1001;
 559    let _1111 = &_100 * &_1011;
 560
 561    // _10000
 562    let mut y = &_1111 * _1;
 563
 564    #[inline]
 565    fn square_multiply(y: &mut Scalar, squarings: usize, x: &Scalar) {
 566      for _ in 0..squarings {
 567        *y = y.square();
 568      }
 569      *y = y.mul(x);
 570    }
 571
 572    square_multiply(&mut y, 123 + 3, &_101);
 573    square_multiply(&mut y, 2 + 2, &_11);
 574    square_multiply(&mut y, 1 + 4, &_1111);
 575    square_multiply(&mut y, 1 + 4, &_1111);
 576    square_multiply(&mut y, 4, &_1001);
 577    square_multiply(&mut y, 2, &_11);
 578    square_multiply(&mut y, 1 + 4, &_1111);
 579    square_multiply(&mut y, 1 + 3, &_101);
 580    square_multiply(&mut y, 3 + 3, &_101);
 581    square_multiply(&mut y, 3, &_111);
 582    square_multiply(&mut y, 1 + 4, &_1111);
 583    square_multiply(&mut y, 2 + 3, &_111);
 584    square_multiply(&mut y, 2 + 2, &_11);
 585    square_multiply(&mut y, 1 + 4, &_1011);
 586    square_multiply(&mut y, 2 + 4, &_1011);
 587    square_multiply(&mut y, 6 + 4, &_1001);
 588    square_multiply(&mut y, 2 + 2, &_11);
 589    square_multiply(&mut y, 3 + 2, &_11);
 590    square_multiply(&mut y, 3 + 2, &_11);
 591    square_multiply(&mut y, 1 + 4, &_1001);
 592    square_multiply(&mut y, 1 + 3, &_111);
 593    square_multiply(&mut y, 2 + 4, &_1111);
 594    square_multiply(&mut y, 1 + 4, &_1011);
 595    square_multiply(&mut y, 3, &_101);
 596    square_multiply(&mut y, 2 + 4, &_1111);
 597    square_multiply(&mut y, 3, &_101);
 598    square_multiply(&mut y, 1 + 2, &_11);
 599
 600    CtOption::new(y, !self.ct_eq(&Self::zero()))
 601  }
 602
 603  pub fn batch_invert(inputs: &mut [Scalar]) -> Scalar {
 604    // This code is essentially identical to the FieldElement
 605    // implementation, and is documented there.  Unfortunately,
 606    // it's not easy to write it generically, since here we want
 607    // to use `UnpackedScalar`s internally, and `Scalar`s
 608    // externally, but there's no corresponding distinction for
 609    // field elements.
 610
 611    use zeroize::Zeroizing;
 612
 613    let n = inputs.len();
 614    let one = Scalar::one();
 615
 616    // Place scratch storage in a Zeroizing wrapper to wipe it when
 617    // we pass out of scope.
 618    let scratch_vec = vec![one; n];
 619    let mut scratch = Zeroizing::new(scratch_vec);
 620
 621    // Keep an accumulator of all of the previous products
 622    let mut acc = Scalar::one();
 623
 624    // Pass through the input vector, recording the previous
 625    // products in the scratch space
 626    for (input, scratch) in inputs.iter().zip(scratch.iter_mut()) {
 627      *scratch = acc;
 628
 629      acc = acc * input;
 630    }
 631
 632    // acc is nonzero iff all inputs are nonzero
 633    debug_assert!(acc != Scalar::zero());
 634
 635    // Compute the inverse of all products
 636    acc = acc.invert().unwrap();
 637
 638    // We need to return the product of all inverses later
 639    let ret = acc;
 640
 641    // Pass through the vector backwards to compute the inverses
 642    // in place
 643    for (input, scratch) in inputs.iter_mut().rev().zip(scratch.iter().rev()) {
 644      let tmp = &acc * input.clone();
 645      *input = &acc * scratch;
 646      acc = tmp;
 647    }
 648
 649    ret
 650  }
 651
 652  #[inline(always)]
 653  const fn montgomery_reduce(
 654    r0: u64,
 655    r1: u64,
 656    r2: u64,
 657    r3: u64,
 658    r4: u64,
 659    r5: u64,
 660    r6: u64,
 661    r7: u64,
 662  ) -> Self {
 663    // The Montgomery reduction here is based on Algorithm 14.32 in
 664    // Handbook of Applied Cryptography
 665    // <http://cacr.uwaterloo.ca/hac/about/chap14.pdf>.
 666
 667    let k = r0.wrapping_mul(INV);
 668    let (_, carry) = mac(r0, k, MODULUS.0[0], 0);
 669    let (r1, carry) = mac(r1, k, MODULUS.0[1], carry);
 670    let (r2, carry) = mac(r2, k, MODULUS.0[2], carry);
 671    let (r3, carry) = mac(r3, k, MODULUS.0[3], carry);
 672    let (r4, carry2) = adc(r4, 0, carry);
 673
 674    let k = r1.wrapping_mul(INV);
 675    let (_, carry) = mac(r1, k, MODULUS.0[0], 0);
 676    let (r2, carry) = mac(r2, k, MODULUS.0[1], carry);
 677    let (r3, carry) = mac(r3, k, MODULUS.0[2], carry);
 678    let (r4, carry) = mac(r4, k, MODULUS.0[3], carry);
 679    let (r5, carry2) = adc(r5, carry2, carry);
 680
 681    let k = r2.wrapping_mul(INV);
 682    let (_, carry) = mac(r2, k, MODULUS.0[0], 0);
 683    let (r3, carry) = mac(r3, k, MODULUS.0[1], carry);
 684    let (r4, carry) = mac(r4, k, MODULUS.0[2], carry);
 685    let (r5, carry) = mac(r5, k, MODULUS.0[3], carry);
 686    let (r6, carry2) = adc(r6, carry2, carry);
 687
 688    let k = r3.wrapping_mul(INV);
 689    let (_, carry) = mac(r3, k, MODULUS.0[0], 0);
 690    let (r4, carry) = mac(r4, k, MODULUS.0[1], carry);
 691    let (r5, carry) = mac(r5, k, MODULUS.0[2], carry);
 692    let (r6, carry) = mac(r6, k, MODULUS.0[3], carry);
 693    let (r7, _) = adc(r7, carry2, carry);
 694
 695    // Result may be within MODULUS of the correct value
 696    (&Scalar([r4, r5, r6, r7])).sub(&MODULUS)
 697  }
 698
 699  /// Multiplies `rhs` by `self`, returning the result.
 700  #[inline]
 701  pub const fn mul(&self, rhs: &Self) -> Self {
 702    // Schoolbook multiplication
 703
 704    let (r0, carry) = mac(0, self.0[0], rhs.0[0], 0);
 705    let (r1, carry) = mac(0, self.0[0], rhs.0[1], carry);
 706    let (r2, carry) = mac(0, self.0[0], rhs.0[2], carry);
 707    let (r3, r4) = mac(0, self.0[0], rhs.0[3], carry);
 708
 709    let (r1, carry) = mac(r1, self.0[1], rhs.0[0], 0);
 710    let (r2, carry) = mac(r2, self.0[1], rhs.0[1], carry);
 711    let (r3, carry) = mac(r3, self.0[1], rhs.0[2], carry);
 712    let (r4, r5) = mac(r4, self.0[1], rhs.0[3], carry);
 713
 714    let (r2, carry) = mac(r2, self.0[2], rhs.0[0], 0);
 715    let (r3, carry) = mac(r3, self.0[2], rhs.0[1], carry);
 716    let (r4, carry) = mac(r4, self.0[2], rhs.0[2], carry);
 717    let (r5, r6) = mac(r5, self.0[2], rhs.0[3], carry);
 718
 719    let (r3, carry) = mac(r3, self.0[3], rhs.0[0], 0);
 720    let (r4, carry) = mac(r4, self.0[3], rhs.0[1], carry);
 721    let (r5, carry) = mac(r5, self.0[3], rhs.0[2], carry);
 722    let (r6, r7) = mac(r6, self.0[3], rhs.0[3], carry);
 723
 724    Scalar::montgomery_reduce(r0, r1, r2, r3, r4, r5, r6, r7)
 725  }
 726
 727  /// Subtracts `rhs` from `self`, returning the result.
 728  #[inline]
 729  pub const fn sub(&self, rhs: &Self) -> Self {
 730    let (d0, borrow) = sbb(self.0[0], rhs.0[0], 0);
 731    let (d1, borrow) = sbb(self.0[1], rhs.0[1], borrow);
 732    let (d2, borrow) = sbb(self.0[2], rhs.0[2], borrow);
 733    let (d3, borrow) = sbb(self.0[3], rhs.0[3], borrow);
 734
 735    // If underflow occurred on the final limb, borrow = 0xfff...fff, otherwise
 736    // borrow = 0x000...000. Thus, we use it as a mask to conditionally add the modulus.
 737    let (d0, carry) = adc(d0, MODULUS.0[0] & borrow, 0);
 738    let (d1, carry) = adc(d1, MODULUS.0[1] & borrow, carry);
 739    let (d2, carry) = adc(d2, MODULUS.0[2] & borrow, carry);
 740    let (d3, _) = adc(d3, MODULUS.0[3] & borrow, carry);
 741
 742    Scalar([d0, d1, d2, d3])
 743  }
 744
 745  /// Adds `rhs` to `self`, returning the result.
 746  #[inline]
 747  pub const fn add(&self, rhs: &Self) -> Self {
 748    let (d0, carry) = adc(self.0[0], rhs.0[0], 0);
 749    let (d1, carry) = adc(self.0[1], rhs.0[1], carry);
 750    let (d2, carry) = adc(self.0[2], rhs.0[2], carry);
 751    let (d3, _) = adc(self.0[3], rhs.0[3], carry);
 752
 753    // Attempt to subtract the modulus, to ensure the value
 754    // is smaller than the modulus.
 755    (&Scalar([d0, d1, d2, d3])).sub(&MODULUS)
 756  }
 757
 758  /// Negates `self`.
 759  #[inline]
 760  pub const fn neg(&self) -> Self {
 761    // Subtract `self` from `MODULUS` to negate. Ignore the final
 762    // borrow because it cannot underflow; self is guaranteed to
 763    // be in the field.
 764    let (d0, borrow) = sbb(MODULUS.0[0], self.0[0], 0);
 765    let (d1, borrow) = sbb(MODULUS.0[1], self.0[1], borrow);
 766    let (d2, borrow) = sbb(MODULUS.0[2], self.0[2], borrow);
 767    let (d3, _) = sbb(MODULUS.0[3], self.0[3], borrow);
 768
 769    // `tmp` could be `MODULUS` if `self` was zero. Create a mask that is
 770    // zero if `self` was zero, and `u64::max_value()` if self was nonzero.
 771    let mask = (((self.0[0] | self.0[1] | self.0[2] | self.0[3]) == 0) as u64).wrapping_sub(1);
 772
 773    Scalar([d0 & mask, d1 & mask, d2 & mask, d3 & mask])
 774  }
 775}
 776
 777impl<'a> From<&'a Scalar> for [u8; 32] {
 778  fn from(value: &'a Scalar) -> [u8; 32] {
 779    value.to_bytes()
 780  }
 781}
 782
 783#[cfg(test)]
 784mod tests {
 785  use super::*;
 786
 787  #[test]
 788  fn test_inv() {
 789    // Compute -(q^{-1} mod 2^64) mod 2^64 by exponentiating
 790    // by totient(2**64) - 1
 791
 792    let mut inv = 1u64;
 793    for _ in 0..63 {
 794      inv = inv.wrapping_mul(inv);
 795      inv = inv.wrapping_mul(MODULUS.0[0]);
 796    }
 797    inv = inv.wrapping_neg();
 798
 799    assert_eq!(inv, INV);
 800  }
 801
 802  #[cfg(feature = "std")]
 803  #[test]
 804  fn test_debug() {
 805    assert_eq!(
 806      format!("{:?}", Scalar::zero()),
 807      "0x0000000000000000000000000000000000000000000000000000000000000000"
 808    );
 809    assert_eq!(
 810      format!("{:?}", Scalar::one()),
 811      "0x0000000000000000000000000000000000000000000000000000000000000001"
 812    );
 813    assert_eq!(
 814      format!("{:?}", R2),
 815      "0x1824b159acc5056f998c4fefecbc4ff55884b7fa0003480200000001fffffffe"
 816    );
 817  }
 818
 819  #[test]
 820  fn test_equality() {
 821    assert_eq!(Scalar::zero(), Scalar::zero());
 822    assert_eq!(Scalar::one(), Scalar::one());
 823    assert_eq!(R2, R2);
 824
 825    assert!(Scalar::zero() != Scalar::one());
 826    assert!(Scalar::one() != R2);
 827  }
 828
 829  #[test]
 830  fn test_to_bytes() {
 831    assert_eq!(
 832      Scalar::zero().to_bytes(),
 833      [
 834        0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
 835        0, 0
 836      ]
 837    );
 838
 839    assert_eq!(
 840      Scalar::one().to_bytes(),
 841      [
 842        1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
 843        0, 0
 844      ]
 845    );
 846
 847    assert_eq!(
 848      R2.to_bytes(),
 849      [
 850        29, 149, 152, 141, 116, 49, 236, 214, 112, 207, 125, 115, 244, 91, 239, 198, 254, 255, 255,
 851        255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 15
 852      ]
 853    );
 854
 855    assert_eq!(
 856      (-&Scalar::one()).to_bytes(),
 857      [
 858        236, 211, 245, 92, 26, 99, 18, 88, 214, 156, 247, 162, 222, 249, 222, 20, 0, 0, 0, 0, 0, 0,
 859        0, 0, 0, 0, 0, 0, 0, 0, 0, 16
 860      ]
 861    );
 862  }
 863
 864  #[test]
 865  fn test_from_bytes() {
 866    assert_eq!(
 867      Scalar::from_bytes(&[
 868        0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
 869        0, 0
 870      ])
 871      .unwrap(),
 872      Scalar::zero()
 873    );
 874
 875    assert_eq!(
 876      Scalar::from_bytes(&[
 877        1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
 878        0, 0
 879      ])
 880      .unwrap(),
 881      Scalar::one()
 882    );
 883
 884    assert_eq!(
 885      Scalar::from_bytes(&[
 886        29, 149, 152, 141, 116, 49, 236, 214, 112, 207, 125, 115, 244, 91, 239, 198, 254, 255, 255,
 887        255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 15
 888      ])
 889      .unwrap(),
 890      R2
 891    );
 892
 893    // -1 should work
 894    assert!(
 895      Scalar::from_bytes(&[
 896        236, 211, 245, 92, 26, 99, 18, 88, 214, 156, 247, 162, 222, 249, 222, 20, 0, 0, 0, 0, 0, 0,
 897        0, 0, 0, 0, 0, 0, 0, 0, 0, 16
 898      ])
 899      .is_some()
 900      .unwrap_u8()
 901        == 1
 902    );
 903
 904    // modulus is invalid
 905    assert!(
 906      Scalar::from_bytes(&[
 907        1, 0, 0, 0, 255, 255, 255, 255, 254, 91, 254, 255, 2, 164, 189, 83, 5, 216, 161, 9, 8, 216,
 908        57, 51, 72, 125, 157, 41, 83, 167, 237, 115
 909      ])
 910      .is_none()
 911      .unwrap_u8()
 912        == 1
 913    );
 914
 915    // Anything larger than the modulus is invalid
 916    assert!(
 917      Scalar::from_bytes(&[
 918        2, 0, 0, 0, 255, 255, 255, 255, 254, 91, 254, 255, 2, 164, 189, 83, 5, 216, 161, 9, 8, 216,
 919        57, 51, 72, 125, 157, 41, 83, 167, 237, 115
 920      ])
 921      .is_none()
 922      .unwrap_u8()
 923        == 1
 924    );
 925    assert!(
 926      Scalar::from_bytes(&[
 927        1, 0, 0, 0, 255, 255, 255, 255, 254, 91, 254, 255, 2, 164, 189, 83, 5, 216, 161, 9, 8, 216,
 928        58, 51, 72, 125, 157, 41, 83, 167, 237, 115
 929      ])
 930      .is_none()
 931      .unwrap_u8()
 932        == 1
 933    );
 934    assert!(
 935      Scalar::from_bytes(&[
 936        1, 0, 0, 0, 255, 255, 255, 255, 254, 91, 254, 255, 2, 164, 189, 83, 5, 216, 161, 9, 8, 216,
 937        57, 51, 72, 125, 157, 41, 83, 167, 237, 116
 938      ])
 939      .is_none()
 940      .unwrap_u8()
 941        == 1
 942    );
 943  }
 944
 945  #[test]
 946  fn test_from_u512_zero() {
 947    assert_eq!(
 948      Scalar::zero(),
 949      Scalar::from_u512([
 950        MODULUS.0[0],
 951        MODULUS.0[1],
 952        MODULUS.0[2],
 953        MODULUS.0[3],
 954        0,
 955        0,
 956        0,
 957        0
 958      ])
 959    );
 960  }
 961
 962  #[test]
 963  fn test_from_u512_r() {
 964    assert_eq!(R, Scalar::from_u512([1, 0, 0, 0, 0, 0, 0, 0]));
 965  }
 966
 967  #[test]
 968  fn test_from_u512_r2() {
 969    assert_eq!(R2, Scalar::from_u512([0, 0, 0, 0, 1, 0, 0, 0]));
 970  }
 971
 972  #[test]
 973  fn test_from_u512_max() {
 974    let max_u64 = 0xffffffffffffffff;
 975    assert_eq!(
 976      R3 - R,
 977      Scalar::from_u512([max_u64, max_u64, max_u64, max_u64, max_u64, max_u64, max_u64, max_u64])
 978    );
 979  }
 980
 981  #[test]
 982  fn test_from_bytes_wide_r2() {
 983    assert_eq!(
 984      R2,
 985      Scalar::from_bytes_wide(&[
 986        29, 149, 152, 141, 116, 49, 236, 214, 112, 207, 125, 115, 244, 91, 239, 198, 254, 255, 255,
 987        255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 15, 0, 0, 0, 0, 0, 0, 0, 0, 0,
 988        0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
 989      ])
 990    );
 991  }
 992
 993  #[test]
 994  fn test_from_bytes_wide_negative_one() {
 995    assert_eq!(
 996      -&Scalar::one(),
 997      Scalar::from_bytes_wide(&[
 998        236, 211, 245, 92, 26, 99, 18, 88, 214, 156, 247, 162, 222, 249, 222, 20, 0, 0, 0, 0, 0, 0,
 999        0, 0, 0, 0, 0, 0, 0, 0, 0, 16, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
1000        0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
1001      ])
1002    );
1003  }
1004
1005  #[test]
1006  fn test_from_bytes_wide_maximum() {
1007    assert_eq!(
1008      Scalar::from_raw([
1009        0xa40611e3449c0f00,
1010        0xd00e1ba768859347,
1011        0xceec73d217f5be65,
1012        0x0399411b7c309a3d
1013      ]),
1014      Scalar::from_bytes_wide(&[0xff; 64])
1015    );
1016  }
1017
1018  #[test]
1019  fn test_zero() {
1020    assert_eq!(Scalar::zero(), -&Scalar::zero());
1021    assert_eq!(Scalar::zero(), Scalar::zero() + Scalar::zero());
1022    assert_eq!(Scalar::zero(), Scalar::zero() - Scalar::zero());
1023    assert_eq!(Scalar::zero(), Scalar::zero() * Scalar::zero());
1024  }
1025
1026  const LARGEST: Scalar = Scalar([
1027    0x5812631a5cf5d3ec,
1028    0x14def9dea2f79cd6,
1029    0x0000000000000000,
1030    0x1000000000000000,
1031  ]);
1032
1033  #[test]
1034  fn test_addition() {
1035    let mut tmp = LARGEST;
1036    tmp += &LARGEST;
1037
1038    assert_eq!(
1039      tmp,
1040      Scalar([
1041        0x5812631a5cf5d3eb,
1042        0x14def9dea2f79cd6,
1043        0x0000000000000000,
1044        0x1000000000000000,
1045      ])
1046    );
1047
1048    let mut tmp = LARGEST;
1049    tmp += &Scalar([1, 0, 0, 0]);
1050
1051    assert_eq!(tmp, Scalar::zero());
1052  }
1053
1054  #[test]
1055  fn test_negation() {
1056    let tmp = -&LARGEST;
1057
1058    assert_eq!(tmp, Scalar([1, 0, 0, 0]));
1059
1060    let tmp = -&Scalar::zero();
1061    assert_eq!(tmp, Scalar::zero());
1062    let tmp = -&Scalar([1, 0, 0, 0]);
1063    assert_eq!(tmp, LARGEST);
1064  }
1065
1066  #[test]
1067  fn test_subtraction() {
1068    let mut tmp = LARGEST;
1069    tmp -= &LARGEST;
1070
1071    assert_eq!(tmp, Scalar::zero());
1072
1073    let mut tmp = Scalar::zero();
1074    tmp -= &LARGEST;
1075
1076    let mut tmp2 = MODULUS;
1077    tmp2 -= &LARGEST;
1078
1079    assert_eq!(tmp, tmp2);
1080  }
1081
1082  #[test]
1083  fn test_multiplication() {
1084    let mut cur = LARGEST;
1085
1086    for _ in 0..100 {
1087      let mut tmp = cur;
1088      tmp *= &cur;
1089
1090      let mut tmp2 = Scalar::zero();
1091      for b in cur
1092        .to_bytes()
1093        .iter()
1094        .rev()
1095        .flat_map(|byte| (0..8).rev().map(move |i| ((byte >> i) & 1u8) == 1u8))
1096      {
1097        let tmp3 = tmp2;
1098        tmp2.add_assign(&tmp3);
1099
1100        if b {
1101          tmp2.add_assign(&cur);
1102        }
1103      }
1104
1105      assert_eq!(tmp, tmp2);
1106
1107      cur.add_assign(&LARGEST);
1108    }
1109  }
1110
1111  #[test]
1112  fn test_squaring() {
1113    let mut cur = LARGEST;
1114
1115    for _ in 0..100 {
1116      let mut tmp = cur;
1117      tmp = tmp.square();
1118
1119      let mut tmp2 = Scalar::zero();
1120      for b in cur
1121        .to_bytes()
1122        .iter()
1123        .rev()
1124        .flat_map(|byte| (0..8).rev().map(move |i| ((byte >> i) & 1u8) == 1u8))
1125      {
1126        let tmp3 = tmp2;
1127        tmp2.add_assign(&tmp3);
1128
1129        if b {
1130          tmp2.add_assign(&cur);
1131        }
1132      }
1133
1134      assert_eq!(tmp, tmp2);
1135
1136      cur.add_assign(&LARGEST);
1137    }
1138  }
1139
1140  #[test]
1141  fn test_inversion() {
1142    assert_eq!(Scalar::zero().invert().is_none().unwrap_u8(), 1);
1143    assert_eq!(Scalar::one().invert().unwrap(), Scalar::one());
1144    assert_eq!((-&Scalar::one()).invert().unwrap(), -&Scalar::one());
1145
1146    let mut tmp = R2;
1147
1148    for _ in 0..100 {
1149      let mut tmp2 = tmp.invert().unwrap();
1150      tmp2.mul_assign(&tmp);
1151
1152      assert_eq!(tmp2, Scalar::one());
1153
1154      tmp.add_assign(&R2);
1155    }
1156  }
1157
1158  #[test]
1159  fn test_invert_is_pow() {
1160    let q_minus_2 = [
1161      0x5812631a5cf5d3eb,
1162      0x14def9dea2f79cd6,
1163      0x0000000000000000,
1164      0x1000000000000000,
1165    ];
1166
1167    let mut r1 = R;
1168    let mut r2 = R;
1169    let mut r3 = R;
1170
1171    for _ in 0..100 {
1172      r1 = r1.invert().unwrap();
1173      r2 = r2.pow_vartime(&q_minus_2);
1174      r3 = r3.pow(&q_minus_2);
1175
1176      assert_eq!(r1, r2);
1177      assert_eq!(r2, r3);
1178      // Add R so we check something different next time around
1179      r1.add_assign(&R);
1180      r2 = r1;
1181      r3 = r1;
1182    }
1183  }
1184
1185  #[test]
1186  fn test_from_raw() {
1187    assert_eq!(
1188      Scalar::from_raw([
1189        0xd6ec31748d98951c,
1190        0xc6ef5bf4737dcf70,
1191        0xfffffffffffffffe,
1192        0x0fffffffffffffff
1193      ]),
1194      Scalar::from_raw([0xffffffffffffffff; 4])
1195    );
1196
1197    assert_eq!(Scalar::from_raw(MODULUS.0), Scalar::zero());
1198
1199    assert_eq!(Scalar::from_raw([1, 0, 0, 0]), R);
1200  }
1201
1202  #[test]
1203  fn test_double() {
1204    let a = Scalar::from_raw([
1205      0x1fff3231233ffffd,
1206      0x4884b7fa00034802,
1207      0x998c4fefecbc4ff3,
1208      0x1824b159acc50562,
1209    ]);
1210
1211    assert_eq!(a.double(), a + a);
1212  }
1213}