daelk-cryptography curve25519-dalek源码解析——之Field表示

https://github.com/dalek-cryptography/curve25519-dalek

1. Scalar结构

针对p Scalar { self.unpack().invert().pack() } /// Unpack this `Scalar` to an `UnpackedScalar` for faster arithmetic. pub(crate) fn unpack(&self) -> UnpackedScalar { UnpackedScalar::from_bytes(&self.bytes) } 3.1 Scalar转换为UnpackedScalar

Scalar转换为UnpackedScalar的代码细节为：

	/// Unpack a 32 byte / 256 bit scalar into 5 52-bit limbs.
    pub fn from_bytes(bytes: &[u8; 32]) -> Scalar52 {
        let mut words = [0u64; 4];
        for i in 0..4 {
            for j in 0..8 {
                words[i] |= (bytes[(i * 8) + j] as u64)  40) | (words[2] > 28) | (words[3] > 16)                     & top_mask;

        s
    }

3.2 invert()操作

有限域内的乘法具有以下特征：

x(p-2) * x = x(p-1) = 1 (mod p)

由此可推测出，求有限域的x值的倒数可转换为求x(p-2)的值。

程序中，对Scalar值求倒数，是先通过unpack()函数将Scalar转换为UnpackedScalar，然后对UnpackedScalar求倒数，最后通过pack()函数将UnpackedScalar转换为Scalar值。

impl Scalar {
	/// let inv_X: Scalar = X.invert();
    /// assert!(XINV == inv_X);
    /// let should_be_one: Scalar = &inv_X * &X;
    /// assert!(should_be_one == Scalar::one());
    /// ```
    pub fn invert(&self) -> Scalar {
        self.unpack().invert().pack()
    }
	/// Unpack this `Scalar` to an `UnpackedScalar` for faster arithmetic.
    pub(crate) fn unpack(&self) -> UnpackedScalar {
        UnpackedScalar::from_bytes(&self.bytes)
    }
}

impl UnpackedScalar {
	/// Inverts an UnpackedScalar not in Montgomery form.
    pub fn invert(&self) -> UnpackedScalar {
        self.to_montgomery().montgomery_invert().from_montgomery()
    }
	/// Pack the limbs of this `UnpackedScalar` into a `Scalar`.
    fn pack(&self) -> Scalar {
        Scalar{ bytes: self.to_bytes() }
    }
}

对于u64_backend feature, 有 type UnpackedScalar = backend::serial::u64::scalar::Scalar52;，所以对于 to_montgomery()的具体实现如下：

impl Scalar52 {
	/// Puts a Scalar52 in to Montgomery form, i.e. computes `a*R (mod l)`
    #[inline(never)]
    pub fn to_montgomery(&self) -> Scalar52 {
        Scalar52::montgomery_mul(self, &constants::RR) //将数组中52*5=260，260bit所有位数都用上。pub struct Scalar52(pub [u64; 5]);
    }
    
	/// Compute `(a * b) / R` (mod l), where R is the Montgomery modulus 2^260
    #[inline(never)]
    pub fn montgomery_mul(a: &Scalar52, b: &Scalar52) -> Scalar52 {
        Scalar52::montgomery_reduce(&Scalar52::mul_internal(a, b))
    }

	/// Compute `a * b`
    #[inline(always)]
    pub (crate) fn mul_internal(a: &Scalar52, b: &Scalar52) -> [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
    }
    
    /// u64 * u64 = u128 multiply helper
	#[inline(always)]
	fn m(x: u64, y: u64) -> u128 {
 	   (x as u128) * (y as u128)
	}
	
 	/// Compute `limbs/R` (mod l), where R is the Montgomery modulus 2^260
    #[inline(always)]
    pub (crate) fn montgomery_reduce(limbs: &[u128; 9]) -> Scalar52 {

        #[inline(always)]
        fn part1(sum: u128) -> (u128, u64) {
            let p = (sum as u64).wrapping_mul(constants::LFACTOR) & ((1u64 > 52, p)
        }

        #[inline(always)]
        fn part2(sum: u128) -> (u128, u64) {
            let w = (sum as u64) & ((1u64 > 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
        Scalar52::sub(&Scalar52([r0,r1,r2,r3,r4]), l)
    }
}

4. constant.rs中常量值sage验证

/// constant.rs中有记录一些常量值。
/// `L` is the order of base point, i.e. 2^252 + 27742317777372353535851937790883648493
pub(crate) const L: Scalar52 = Scalar52([ 0x0002631a5cf5d3ed, 0x000dea2f79cd6581, 0x000000000014def9, 0x0000000000000000, 0x0000100000000000 ]);

/// 其实即为L[0]*LFACTOR = -1 (mod 2^52) = 2^52-1 (mod 2^52)
/// (L[i]