ark_ff_optimized/
fp64.rs

1//! An implementation of a 64-bit STARK-friendly prime field with modulus `2^64 -
2//! 2^32 + 1`. The implementation follows <https://eprint.iacr.org/2022/274.pdf>
3//! and the code for the majority of functions was taken and adapted from
4//! <https://github.com/novifinancial/winterfell>
5//!
6//! This field and its implementation has many attractive properties:
7//! * Multiplication of two 32-bit values does not overflow field modulus.
8//! * Field arithmetic in this field can be implemented using a few 32-bit
9//!   addition, subtractions, and shifts.
10//! * 8 is the 64th root of unity which opens up potential for optimized FFT
11//!   implementations.
12
13use ark_ff::{fields::Fp64, BigInt, PrimeField, SqrtPrecomputation, Zero};
14use core::marker::PhantomData;
15
16/// Field modulus `p = 2^64 - 2^32 + 1`
17const MODULUS: u64 = 18_446_744_069_414_584_321;
18
19/// Square of auxiliary modulus R for Montgomery reduction `R2 ≡ (2^64)^2 mod p`
20const R2: u64 = 18_446_744_065_119_617_025;
21
22pub struct FpParams;
23impl ark_ff::FpConfig<1> for FpParams {
24    const MODULUS: ark_ff::BigInt<1> = BigInt([MODULUS]);
25    const GENERATOR: Fp64<Self> = into_mont(7);
26    const ZERO: Fp64<Self> = into_mont(0);
27    const ONE: Fp64<Self> = into_mont(1);
28    const TWO_ADICITY: u32 = 32;
29    const TWO_ADIC_ROOT_OF_UNITY: Fp64<Self> = into_mont(1_753_635_133_440_165_772);
30    const SQRT_PRECOMP: Option<ark_ff::SqrtPrecomputation<Fp64<Self>>> =
31        Some(SqrtPrecomputation::TonelliShanks {
32            two_adicity: Self::TWO_ADICITY,
33            quadratic_nonresidue_to_trace: Self::TWO_ADIC_ROOT_OF_UNITY,
34            trace_of_modulus_minus_one_div_two: &<Fp64<Self>>::TRACE_MINUS_ONE_DIV_TWO.0,
35        });
36
37    fn add_assign(a: &mut Fp64<Self>, b: &Fp64<Self>) {
38        // We compute a + b = a - (p - b).
39        let (x1, c1) = (a.0).0[0].overflowing_sub(MODULUS - (b.0).0[0]);
40        let adj = 0u32.wrapping_sub(u32::from(c1));
41        (a.0).0[0] = x1.wrapping_sub(u64::from(adj));
42    }
43
44    fn sub_assign(a: &mut Fp64<Self>, b: &Fp64<Self>) {
45        let (x1, c1) = (a.0).0[0].overflowing_sub((b.0).0[0]);
46        let adj = 0u32.wrapping_sub(u32::from(c1));
47        (a.0).0[0] = x1.wrapping_sub(u64::from(adj));
48    }
49
50    fn double_in_place(a: &mut Fp64<Self>) {
51        Self::add_assign(a, &a.clone());
52    }
53
54    fn mul_assign(a: &mut Fp64<Self>, b: &Fp64<Self>) {
55        (a.0).0[0] = mont_red(u128::from((a.0).0[0]) * u128::from((b.0).0[0]));
56    }
57
58    fn sum_of_products<const T: usize>(a: &[Fp64<Self>; T], b: &[Fp64<Self>; T]) -> Fp64<Self> {
59        a.iter().zip(b).map(|(&a, b)| a * b).sum()
60    }
61
62    fn square_in_place(a: &mut Fp64<Self>) {
63        let temp = *a;
64        Self::mul_assign(a, &temp);
65    }
66
67    fn inverse(a: &Fp64<Self>) -> Option<Fp64<Self>> {
68        if a.is_zero() {
69            None
70        } else {
71            let a = (a.0).0[0];
72            let t2 = exp_acc::<1>(a, a);
73            let t3 = exp_acc::<1>(t2, a);
74            let t6 = exp_acc::<3>(t3, t3);
75            let t12 = exp_acc::<6>(t6, t6);
76            let t24 = exp_acc::<12>(t12, t12);
77            let t30 = exp_acc::<6>(t24, t6);
78            let t31 = exp_acc::<1>(t30, a);
79            let t63 = exp_acc::<32>(t31, t31);
80            let inv = exp_acc::<1>(t63, a);
81            Some(ark_ff::Fp(BigInt([inv]), PhantomData))
82        }
83    }
84
85    fn from_bigint(other: ark_ff::BigInt<1>) -> Option<Fp64<Self>> {
86        let inner = other.0[0];
87        if inner.is_zero() {
88            Some(Self::ZERO)
89        } else if inner < MODULUS {
90            Some(into_mont(other.0[0]))
91        } else {
92            None
93        }
94    }
95
96    fn into_bigint(other: Fp64<Self>) -> ark_ff::BigInt<1> {
97        BigInt([mont_red(u128::from((other.0).0[0]))])
98    }
99
100    fn neg_in_place(a: &mut Fp64<Self>) {
101        let mut tmp = Self::ZERO;
102        Self::sub_assign(&mut tmp, a);
103        a.0 = tmp.0;
104    }
105}
106
107/// An optimized implementation of a 64-bit prime field with modulus `2^64 -
108/// 2^32 + 1`
109pub type Fp = Fp64<FpParams>;
110
111/// Converts a value into Montgomery representation
112#[inline]
113const fn into_mont(value: u64) -> Fp {
114    ark_ff::Fp(BigInt([mont_red(value as u128 * R2 as u128)]), PhantomData)
115}
116
117/// Performs Montgomery reduction
118#[inline]
119const fn mont_red(x: u128) -> u64 {
120    // See reference above for a description of the following implementation.
121    #[allow(clippy::cast_possible_truncation)]
122    let xl = x as u64;
123    let xh = (x >> 64) as u64;
124    let (a, overflow) = xl.overflowing_add(xl << 32);
125    let b = a.wrapping_sub(a >> 32).wrapping_sub(overflow as u64);
126    let (r, underflow) = xh.overflowing_sub(b);
127    r.wrapping_sub(0u32.wrapping_sub(underflow as u32) as u64)
128}
129
130/// Squares `base` N times and multiplies the result by the tail value.
131#[inline]
132const fn exp_acc<const N: usize>(base: u64, tail: u64) -> u64 {
133    let mut result = base;
134    let mut i = 0;
135    while i < N {
136        result = mont_red(result as u128 * result as u128);
137        i += 1;
138    }
139    mont_red(result as u128 * tail as u128)
140}
141
142#[cfg(test)]
143mod tests {
144    use super::Fp as TestField;
145    use ark_algebra_test_templates::test_field;
146
147    test_field!(generated; TestField; prime);
148}