From 3630b027ff059a60ea5fa924b6d84e970fcaa900 Mon Sep 17 00:00:00 2001 From: Alexander Wu Date: Sun, 18 Dec 2022 15:11:16 -0800 Subject: [PATCH 1/8] Add minor clarification in README about tests --- README.md | 2 ++ 1 file changed, 2 insertions(+) diff --git a/README.md b/README.md index e62ca3593..a81196d4e 100644 --- a/README.md +++ b/README.md @@ -49,6 +49,8 @@ This library comes with comprehensive unit and integration tests for each of the cargo test --all ``` +Note that this leaves out some tests that would be run in individual crates (e.g. `test-curves`) under `cargo test --all-features`. + ## Benchmarks To run the benchmarks, install the nightly Rust toolchain, via `rustup install nightly`, and then run the following command: From 76f8b8ebe49ffdfcae87596882a3a3bc105173d1 Mon Sep 17 00:00:00 2001 From: Alexander Wu Date: Sun, 18 Dec 2022 16:22:24 -0800 Subject: [PATCH 2/8] Move sqrt algorithms to individual functions and rename Case3Mod4 -> PowerCase3Mod4 --- ff/src/fields/models/fp/montgomery_backend.rs | 2 +- ff/src/fields/sqrt.rs | 129 ++++++++++-------- test-templates/src/fields.rs | 2 +- 3 files changed, 75 insertions(+), 58 deletions(-) diff --git a/ff/src/fields/models/fp/montgomery_backend.rs b/ff/src/fields/models/fp/montgomery_backend.rs index 76344227c..4b3a80219 100644 --- a/ff/src/fields/models/fp/montgomery_backend.rs +++ b/ff/src/fields/models/fp/montgomery_backend.rs @@ -545,7 +545,7 @@ pub const fn sqrt_precomputation>( ) -> Option, N>>> { match T::MODULUS.mod_4() { 3 => match T::MODULUS_PLUS_ONE_DIV_FOUR.as_ref() { - Some(BigInt(modulus_plus_one_div_four)) => Some(SqrtPrecomputation::Case3Mod4 { + Some(BigInt(modulus_plus_one_div_four)) => Some(SqrtPrecomputation::PowerCase3Mod4 { modulus_plus_one_div_four, }), None => None, diff --git a/ff/src/fields/sqrt.rs b/ff/src/fields/sqrt.rs index 6f3525988..f4f87fe89 100644 --- a/ff/src/fields/sqrt.rs +++ b/ff/src/fields/sqrt.rs @@ -72,7 +72,7 @@ pub enum SqrtPrecomputation { trace_of_modulus_minus_one_div_two: &'static [u64], }, /// To be used when the modulus is 3 mod 4. - Case3Mod4 { + PowerCase3Mod4 { modulus_plus_one_div_four: &'static [u64], }, } @@ -84,66 +84,83 @@ impl SqrtPrecomputation { two_adicity, quadratic_nonresidue_to_trace, trace_of_modulus_minus_one_div_two, - } => { - // https://eprint.iacr.org/2012/685.pdf (page 12, algorithm 5) - // Actually this is just normal Tonelli-Shanks; since `P::Generator` - // is a quadratic non-residue, `P::ROOT_OF_UNITY = P::GENERATOR ^ t` - // is also a quadratic non-residue (since `t` is odd). - if elem.is_zero() { - return Some(F::zero()); - } - // Try computing the square root (x at the end of the algorithm) - // Check at the end of the algorithm if x was a square root - // Begin Tonelli-Shanks - let mut z = *quadratic_nonresidue_to_trace; - let mut w = elem.pow(trace_of_modulus_minus_one_div_two); - let mut x = w * elem; - let mut b = x * &w; + } => tonelli_shanks( + elem, + two_adicity, + quadratic_nonresidue_to_trace, + trace_of_modulus_minus_one_div_two, + ), + Self::PowerCase3Mod4 { + modulus_plus_one_div_four, + } => power_case_three_mod_four(elem, modulus_plus_one_div_four), + } + } +} - let mut v = *two_adicity as usize; +/// https://eprint.iacr.org/2012/685.pdf (page 12, algorithm 5) +fn tonelli_shanks( + elem: &F, + two_adicity: &u32, + quadratic_nonresidue_to_trace: &F, + trace_of_modulus_minus_one_div_two: &'static [u64], +) -> Option { + // Actually this is just normal Tonelli-Shanks; since `P::Generator` + // is a quadratic non-residue, `P::ROOT_OF_UNITY = P::GENERATOR ^ t` + // is also a quadratic non-residue (since `t` is odd). + if elem.is_zero() { + return Some(F::zero()); + } + // Try computing the square root (x at the end of the algorithm) + // Check at the end of the algorithm if x was a square root + // Begin Tonelli-Shanks + let mut z = *quadratic_nonresidue_to_trace; + let mut w = elem.pow(trace_of_modulus_minus_one_div_two); + let mut x = w * elem; + let mut b = x * &w; - while !b.is_one() { - let mut k = 0usize; + let mut v = *two_adicity as usize; - let mut b2k = b; - while !b2k.is_one() { - // invariant: b2k = b^(2^k) after entering this loop - b2k.square_in_place(); - k += 1; - } + while !b.is_one() { + let mut k = 0usize; - if k == (*two_adicity as usize) { - // We are in the case where self^(T * 2^k) = x^(P::MODULUS - 1) = 1, - // which means that no square root exists. - return None; - } - let j = v - k; - w = z; - for _ in 1..j { - w.square_in_place(); - } + let mut b2k = b; + while !b2k.is_one() { + // invariant: b2k = b^(2^k) after entering this loop + b2k.square_in_place(); + k += 1; + } - z = w.square(); - b *= &z; - x *= &w; - v = k; - } - // Is x the square root? If so, return it. - if x.square() == *elem { - Some(x) - } else { - // Consistency check that if no square root is found, - // it is because none exists. - debug_assert!(!matches!(elem.legendre(), LegendreSymbol::QuadraticResidue)); - None - } - }, - Self::Case3Mod4 { - modulus_plus_one_div_four, - } => { - let result = elem.pow(modulus_plus_one_div_four.as_ref()); - (result.square() == *elem).then_some(result) - }, + if k == (*two_adicity as usize) { + // We are in the case where self^(T * 2^k) = x^(P::MODULUS - 1) = 1, + // which means that no square root exists. + return None; } + let j = v - k; + w = z; + for _ in 1..j { + w.square_in_place(); + } + + z = w.square(); + b *= &z; + x *= &w; + v = k; } + // Is x the square root? If so, return it. + if x.square() == *elem { + Some(x) + } else { + // Consistency check that if no square root is found, + // it is because none exists. + debug_assert!(!matches!(elem.legendre(), LegendreSymbol::QuadraticResidue)); + None + } +} + +fn power_case_three_mod_four( + elem: &F, + modulus_plus_one_div_four: &&'static [u64], +) -> Option { + let result = elem.pow(modulus_plus_one_div_four.as_ref()); + (result.square() == *elem).then_some(result) } diff --git a/test-templates/src/fields.rs b/test-templates/src/fields.rs index 309ea79ce..15bb18005 100644 --- a/test-templates/src/fields.rs +++ b/test-templates/src/fields.rs @@ -430,7 +430,7 @@ macro_rules! __test_field { let modulus_minus_one = &modulus - 1u8; assert_eq!(BigUint::from(<$field>::MODULUS_MINUS_ONE_DIV_TWO), &modulus_minus_one / 2u32); assert_eq!(<$field>::MODULUS_BIT_SIZE as u64, modulus.bits()); - if let Some(SqrtPrecomputation::Case3Mod4 { modulus_plus_one_div_four }) = <$field>::SQRT_PRECOMP { + if let Some(SqrtPrecomputation::PowerCase3Mod4 { modulus_plus_one_div_four }) = <$field>::SQRT_PRECOMP { // Handle the case where `(MODULUS + 1) / 4` // has fewer limbs than `MODULUS`. let check = ((&modulus + 1u8) / 4u8).to_u64_digits(); From 6ec9d448f92e57d154a153483b7f0ea6089b2976 Mon Sep 17 00:00:00 2001 From: Alexander Wu Date: Thu, 5 Jan 2023 12:36:44 -0500 Subject: [PATCH 3/8] Add functions that implement Shanks, Atkin, and Kong sqrt --- ff/src/fields/sqrt.rs | 145 +++++++++++++++++++++++++++++++++++++++++- 1 file changed, 142 insertions(+), 3 deletions(-) diff --git a/ff/src/fields/sqrt.rs b/ff/src/fields/sqrt.rs index f4f87fe89..083a04f07 100644 --- a/ff/src/fields/sqrt.rs +++ b/ff/src/fields/sqrt.rs @@ -65,13 +65,20 @@ impl LegendreSymbol { /// the corresponding condition. #[non_exhaustive] pub enum SqrtPrecomputation { - // Tonelli-Shanks algorithm works for all elements, no matter what the modulus is. + /// https://eprint.iacr.org/2012/685.pdf (page 12, algorithm 5). + /// Tonelli-Shanks algorithm works for all elements, no matter what the modulus is. + /// With _q_ as field order, _p_ as characteristic, and _m_ as extension degree: + /// * First factor _q - 1 = 2^s t_ where _t_ is odd. + /// * `two_adicity` - _s_. + /// * `quadratic_nonresidue_to_trace` - _c^t_, with random _c_ such that _c^2^(s - 1) = 1_. + /// * `trace_of_modulus_minus_one_div_two` - _(t - 1)/2 + 1_. TonelliShanks { two_adicity: u32, quadratic_nonresidue_to_trace: F, trace_of_modulus_minus_one_div_two: &'static [u64], }, - /// To be used when the modulus is 3 mod 4. + /// In the case of 3 mod 4, we can find the square root via an exponentiation, + /// sqrt(a) = a^(p+1)/4. This can be proved using Euler's criterion, a^(p-1)/2 = 1 mod p. PowerCase3Mod4 { modulus_plus_one_div_four: &'static [u64], }, @@ -97,7 +104,6 @@ impl SqrtPrecomputation { } } -/// https://eprint.iacr.org/2012/685.pdf (page 12, algorithm 5) fn tonelli_shanks( elem: &F, two_adicity: &u32, @@ -157,6 +163,139 @@ fn tonelli_shanks( } } +/// https://eprint.iacr.org/2012/685.pdf (page 9, algorithm 2). +/// With _q_ as field order, _p_ as characteristic, and _m_ as extension degree: +/// * `char_minus_three_div_four` - _(p - 3)/4_. +/// * `deg_minus_three_div_two_plus_one` - _(m - 3)/2 + 1_. +fn shanks( + elem: &F, + char_minus_three_div_four: &'static [u64], + deg_minus_three_div_two_plus_one: usize, +) -> Option { + // Computing a1 = Using decomposition of (q-3)/4 = a + p[pa + (3a+2)] * sum_i=1^(m-3)/2 p^2i + // where a = (p - 3) / 4. + // factor1 = elem^a + let factor1 = elem.pow(char_minus_three_div_four); + // elem_p = elem^p + let mut elem_p = elem.frobenius_map_pure(1); + // factor2_base = elem^(p^2)a * elem^3pa * elem^2p + let mut factor2_base = elem_p.frobenius_map_pure(1).pow(char_minus_three_div_four) + * elem_p.pow(&[3u64]).pow(char_minus_three_div_four) + * elem_p.square(); + // factor2 = prod_i=1^(m-3)/2 factor2_base^(p^2i) + let mut factor2 = F::one(); + for i in 1..deg_minus_three_div_two_plus_one { + factor2 *= factor2_base.frobenius_map_pure(i * 2 as usize); + } + let mut a1 = factor1 * factor2; + + let mut a0 = a1 * a1 * elem; + if a0 == -F::one() { + return None; + } + + let x = a1 * elem; + Some(x) +} + +/// https://eprint.iacr.org/2012/685.pdf (page 10, algorithm 3). +/// With _q_ as field order, _p_ as characteristic, and _m_ as extension degree: +/// * `trace` - _2^(q - 5)/8_. +/// * `char_minus_five_div_eight` - _(p - 5)/8_. +/// * `deg_minus_three_div_two_plus_one` - _(m - 3)/2 + 1_. +fn atkin( + elem: &F, + trace: &F, + char_minus_five_div_eight: &'static [u64], + deg_minus_three_div_two_plus_one: usize, +) -> Option { + // Computing a1 = elem^(q-5)/8 using decomposition of + // (q-5)/8 = a + p[pa + (5a+3)] * sum_i=1^(m-3)/2 p^2i + // where a = (p - 5) / 8. + // factor1 = elem^a + let factor1 = elem.pow(char_minus_five_div_eight); + // elem_p = elem^p + let mut elem_p = elem.frobenius_map_pure(1); + // factor2_base = elem^(p^2)a * elem^5pa * elem^3p + let factor2_base = elem_p.frobenius_map_pure(1).pow(char_minus_five_div_eight) + * elem_p.pow(&[5u64]).pow(char_minus_five_div_eight) + * elem_p.pow(&[3u64]); + // factor2 = prod_i=1^(m-3)/2 factor2_base^(p^2i) + let mut factor2 = F::one(); + for i in 1..deg_minus_three_div_two_plus_one { + factor2 *= factor2_base.frobenius_map_pure(2 * i); + } + let mut a1 = factor1 * factor2; + + let mut a0 = a1 * a1 * elem; + a0 *= a0; + if a0 == -F::one() { + return None; + } + + let b = a1 * trace; + let i = elem.double() * b * b; + let x = b * elem * (i - F::one()); + + Some(x) +} + +/// https://eprint.iacr.org/2012/685.pdf (page 11, algorithm 4). +/// With _q_ as field order, _p_ as characteristic, and _m_ as extension degree: +/// * `trace` - _2^(q - 9)/16_. +/// * `c` - nonzero value such that _chi_q(c) != 1_. +/// * `d` - _c^(q - 9)/8_. +/// * `c_squared` - _c^2_. +/// * `char_minus_nine_div_sixteen` - _(p - 9)/16_. +/// * `deg_minus_three_div_two_plus_one` - _(m - 3)/2 + 1_. +fn kong( + elem: &F, + trace: &F, + c: &F, + d: &F, + c_squared: &F, + char_minus_nine_div_sixteen: &'static [u64], + deg_minus_three_div_two_plus_one: usize, +) -> Option { + // Using decomposition of (q-9)/16 = a + p[pa + (9a+5)] * sum_i=1^(m-3)/2 p^2i + // a = (p - 9) / 16 + // factor1 = elem^a + let factor1 = elem.pow(char_minus_nine_div_sixteen); + // elem_p = elem^p + let elem_p = elem.frobenius_map_pure(1); + // factor2_base = elem^(p^2)a * elem^9pa * elem^5p + let factor2_base = elem_p + .frobenius_map_pure(1) + .pow(char_minus_nine_div_sixteen) + * elem_p.pow(&[9u64]).pow(char_minus_nine_div_sixteen) + * elem_p.pow(&[5u64]); + // factor2 = prod_i=1^(m-3)/2 factor2_base^(p^2i) + let mut factor2 = F::one(); + for i in 1..deg_minus_three_div_two_plus_one { + factor2 *= factor2_base.frobenius_map_pure(2 * i); + } + let a1 = factor1 * factor2; + + let mut a0 = a1 * a1 * elem; + a0 = a0.pow(&[4u64]); + if a0 == -F::one() { + return None; + } + + let b = a1 * trace; + let i = elem.double() * b * b; + let r = i * i; + if r == -F::one() { + let x = b * elem * (i - F::one()); + return Some(x); + } + + let u = b * d; + let i = u.double() * u * c_squared * elem; + let x = u * c * elem * (i - F::one()); + Some(x) +} + fn power_case_three_mod_four( elem: &F, modulus_plus_one_div_four: &&'static [u64], From 644794fec306813cc8d98e43abf5f3409e9f83e9 Mon Sep 17 00:00:00 2001 From: Alexander Wu Date: Thu, 5 Jan 2023 12:43:02 -0500 Subject: [PATCH 4/8] Rebase master and use updated frobenius_map --- ff/src/fields/sqrt.rs | 26 ++++++++++++-------------- 1 file changed, 12 insertions(+), 14 deletions(-) diff --git a/ff/src/fields/sqrt.rs b/ff/src/fields/sqrt.rs index 083a04f07..8899d545f 100644 --- a/ff/src/fields/sqrt.rs +++ b/ff/src/fields/sqrt.rs @@ -177,19 +177,19 @@ fn shanks( // factor1 = elem^a let factor1 = elem.pow(char_minus_three_div_four); // elem_p = elem^p - let mut elem_p = elem.frobenius_map_pure(1); + let elem_p = elem.frobenius_map(1); // factor2_base = elem^(p^2)a * elem^3pa * elem^2p - let mut factor2_base = elem_p.frobenius_map_pure(1).pow(char_minus_three_div_four) + let factor2_base = elem_p.frobenius_map(1).pow(char_minus_three_div_four) * elem_p.pow(&[3u64]).pow(char_minus_three_div_four) * elem_p.square(); // factor2 = prod_i=1^(m-3)/2 factor2_base^(p^2i) let mut factor2 = F::one(); for i in 1..deg_minus_three_div_two_plus_one { - factor2 *= factor2_base.frobenius_map_pure(i * 2 as usize); + factor2 *= factor2_base.frobenius_map(i * 2 as usize); } - let mut a1 = factor1 * factor2; + let a1 = factor1 * factor2; - let mut a0 = a1 * a1 * elem; + let a0 = a1 * a1 * elem; if a0 == -F::one() { return None; } @@ -215,17 +215,17 @@ fn atkin( // factor1 = elem^a let factor1 = elem.pow(char_minus_five_div_eight); // elem_p = elem^p - let mut elem_p = elem.frobenius_map_pure(1); + let elem_p = elem.frobenius_map(1); // factor2_base = elem^(p^2)a * elem^5pa * elem^3p - let factor2_base = elem_p.frobenius_map_pure(1).pow(char_minus_five_div_eight) + let factor2_base = elem_p.frobenius_map(1).pow(char_minus_five_div_eight) * elem_p.pow(&[5u64]).pow(char_minus_five_div_eight) * elem_p.pow(&[3u64]); // factor2 = prod_i=1^(m-3)/2 factor2_base^(p^2i) let mut factor2 = F::one(); for i in 1..deg_minus_three_div_two_plus_one { - factor2 *= factor2_base.frobenius_map_pure(2 * i); + factor2 *= factor2_base.frobenius_map(2 * i); } - let mut a1 = factor1 * factor2; + let a1 = factor1 * factor2; let mut a0 = a1 * a1 * elem; a0 *= a0; @@ -262,17 +262,15 @@ fn kong( // factor1 = elem^a let factor1 = elem.pow(char_minus_nine_div_sixteen); // elem_p = elem^p - let elem_p = elem.frobenius_map_pure(1); + let elem_p = elem.frobenius_map(1); // factor2_base = elem^(p^2)a * elem^9pa * elem^5p - let factor2_base = elem_p - .frobenius_map_pure(1) - .pow(char_minus_nine_div_sixteen) + let factor2_base = elem_p.frobenius_map(1).pow(char_minus_nine_div_sixteen) * elem_p.pow(&[9u64]).pow(char_minus_nine_div_sixteen) * elem_p.pow(&[5u64]); // factor2 = prod_i=1^(m-3)/2 factor2_base^(p^2i) let mut factor2 = F::one(); for i in 1..deg_minus_three_div_two_plus_one { - factor2 *= factor2_base.frobenius_map_pure(2 * i); + factor2 *= factor2_base.frobenius_map(2 * i); } let a1 = factor1 * factor2; From 427ce35c3b4bc1fda680b21760473a43bfa74912 Mon Sep 17 00:00:00 2001 From: Alexander Wu Date: Thu, 5 Jan 2023 15:12:12 -0500 Subject: [PATCH 5/8] Update algorithms to reuse previous results and .square(), and rename elem_p --- ff/src/fields/sqrt.rs | 56 +++++++++++++++++++++---------------------- 1 file changed, 28 insertions(+), 28 deletions(-) diff --git a/ff/src/fields/sqrt.rs b/ff/src/fields/sqrt.rs index 8899d545f..1dfb299e8 100644 --- a/ff/src/fields/sqrt.rs +++ b/ff/src/fields/sqrt.rs @@ -176,12 +176,12 @@ fn shanks( // where a = (p - 3) / 4. // factor1 = elem^a let factor1 = elem.pow(char_minus_three_div_four); - // elem_p = elem^p - let elem_p = elem.frobenius_map(1); + // elem_to_p = elem^p + let elem_to_p = elem.frobenius_map(1); // factor2_base = elem^(p^2)a * elem^3pa * elem^2p - let factor2_base = elem_p.frobenius_map(1).pow(char_minus_three_div_four) - * elem_p.pow(&[3u64]).pow(char_minus_three_div_four) - * elem_p.square(); + let factor2_base = elem_to_p.frobenius_map(1).pow(char_minus_three_div_four) + * elem_to_p.pow(&[3u64]).pow(char_minus_three_div_four) + * elem_to_p.square(); // factor2 = prod_i=1^(m-3)/2 factor2_base^(p^2i) let mut factor2 = F::one(); for i in 1..deg_minus_three_div_two_plus_one { @@ -189,13 +189,13 @@ fn shanks( } let a1 = factor1 * factor2; - let a0 = a1 * a1 * elem; + let a1_elem = a1 * elem; + let a0 = a1 * a1_elem; if a0 == -F::one() { return None; } - let x = a1 * elem; - Some(x) + Some(a1_elem) } /// https://eprint.iacr.org/2012/685.pdf (page 10, algorithm 3). @@ -214,12 +214,12 @@ fn atkin( // where a = (p - 5) / 8. // factor1 = elem^a let factor1 = elem.pow(char_minus_five_div_eight); - // elem_p = elem^p - let elem_p = elem.frobenius_map(1); + // elem_to_p = elem^p + let elem_to_p = elem.frobenius_map(1); // factor2_base = elem^(p^2)a * elem^5pa * elem^3p - let factor2_base = elem_p.frobenius_map(1).pow(char_minus_five_div_eight) - * elem_p.pow(&[5u64]).pow(char_minus_five_div_eight) - * elem_p.pow(&[3u64]); + let factor2_base = elem_to_p.frobenius_map(1).pow(char_minus_five_div_eight) + * elem_to_p.pow(&[5u64]).pow(char_minus_five_div_eight) + * elem_to_p.pow(&[3u64]); // factor2 = prod_i=1^(m-3)/2 factor2_base^(p^2i) let mut factor2 = F::one(); for i in 1..deg_minus_three_div_two_plus_one { @@ -227,15 +227,15 @@ fn atkin( } let a1 = factor1 * factor2; - let mut a0 = a1 * a1 * elem; - a0 *= a0; + let a0 = (a1.square() * elem).square(); if a0 == -F::one() { return None; } let b = a1 * trace; - let i = elem.double() * b * b; - let x = b * elem * (i - F::one()); + let elem_b = b * elem; + let i = elem_b.double() * b; + let x = elem_b * (i - F::one()); Some(x) } @@ -261,12 +261,12 @@ fn kong( // a = (p - 9) / 16 // factor1 = elem^a let factor1 = elem.pow(char_minus_nine_div_sixteen); - // elem_p = elem^p - let elem_p = elem.frobenius_map(1); + // elem_to_p = elem^p + let elem_to_p = elem.frobenius_map(1); // factor2_base = elem^(p^2)a * elem^9pa * elem^5p - let factor2_base = elem_p.frobenius_map(1).pow(char_minus_nine_div_sixteen) - * elem_p.pow(&[9u64]).pow(char_minus_nine_div_sixteen) - * elem_p.pow(&[5u64]); + let factor2_base = elem_to_p.frobenius_map(1).pow(char_minus_nine_div_sixteen) + * elem_to_p.pow(&[9u64]).pow(char_minus_nine_div_sixteen) + * elem_to_p.pow(&[5u64]); // factor2 = prod_i=1^(m-3)/2 factor2_base^(p^2i) let mut factor2 = F::one(); for i in 1..deg_minus_three_div_two_plus_one { @@ -274,22 +274,22 @@ fn kong( } let a1 = factor1 * factor2; - let mut a0 = a1 * a1 * elem; - a0 = a0.pow(&[4u64]); + let a0 = (a1.square() * elem).pow(&[4u64]); if a0 == -F::one() { return None; } let b = a1 * trace; - let i = elem.double() * b * b; - let r = i * i; + let elem_b = b * elem; + let mut i = elem_b.double() * b; + let r = i.square(); if r == -F::one() { - let x = b * elem * (i - F::one()); + let x = elem_b * (i - F::one()); return Some(x); } let u = b * d; - let i = u.double() * u * c_squared * elem; + i = u.square().double() * c_squared * elem; let x = u * c * elem * (i - F::one()); Some(x) } From b8b0176009776bf37097925f63d2cdea977e1ec5 Mon Sep 17 00:00:00 2001 From: Alexander Wu Date: Sun, 15 Jan 2023 16:44:17 -0800 Subject: [PATCH 6/8] Add enum cases for new sqrt algorithm functions --- ff/src/fields/sqrt.rs | 95 +++++++++++++++++++++++++++++++++---------- 1 file changed, 73 insertions(+), 22 deletions(-) diff --git a/ff/src/fields/sqrt.rs b/ff/src/fields/sqrt.rs index 1dfb299e8..778893246 100644 --- a/ff/src/fields/sqrt.rs +++ b/ff/src/fields/sqrt.rs @@ -77,6 +77,40 @@ pub enum SqrtPrecomputation { quadratic_nonresidue_to_trace: F, trace_of_modulus_minus_one_div_two: &'static [u64], }, + /// https://eprint.iacr.org/2012/685.pdf (page 9, algorithm 2). + /// With _q_ as field order, _p_ as characteristic, and _m_ as extension degree: + /// * `char_minus_three_div_four` - _(p - 3)/4_. + /// * `deg_minus_three_div_two_plus_one` - _(m - 3)/2 + 1_. + ShanksCase3Mod4 { + char_minus_three_div_four: &'static [u64], + deg_minus_three_div_two_plus_one: usize, + }, + /// https://eprint.iacr.org/2012/685.pdf (page 10, algorithm 3). + /// With _q_ as field order, _p_ as characteristic, and _m_ as extension degree: + /// * `trace` - _2^(q - 5)/8_. + /// * `char_minus_five_div_eight` - _(p - 5)/8_. + /// * `deg_minus_three_div_two_plus_one` - _(m - 3)/2 + 1_. + AtkinCase5Mod8 { + trace: F, + char_minus_five_div_eight: &'static [u64], + deg_minus_three_div_two_plus_one: usize, + }, + /// https://eprint.iacr.org/2012/685.pdf (page 11, algorithm 4). + /// With _q_ as field order, _p_ as characteristic, and _m_ as extension degree: + /// * `trace` - _2^(q - 9)/16_. + /// * `c` - nonzero value such that _chi_q(c) != 1_. + /// * `d` - _c^(q - 9)/8_. + /// * `c_squared` - _c^2_. + /// * `char_minus_nine_div_sixteen` - _(p - 9)/16_. + /// * `deg_minus_three_div_two_plus_one` - _(m - 3)/2 + 1_. + KongCase9Mod16 { + trace: F, + c: F, + d: F, + c_squared: F, + char_minus_nine_div_sixteen: &'static [u64], + deg_minus_three_div_two_plus_one: usize, + }, /// In the case of 3 mod 4, we can find the square root via an exponentiation, /// sqrt(a) = a^(p+1)/4. This can be proved using Euler's criterion, a^(p-1)/2 = 1 mod p. PowerCase3Mod4 { @@ -97,6 +131,40 @@ impl SqrtPrecomputation { quadratic_nonresidue_to_trace, trace_of_modulus_minus_one_div_two, ), + SqrtPrecomputation::ShanksCase3Mod4 { + char_minus_three_div_four, + deg_minus_three_div_two_plus_one, + } => shanks( + elem, + char_minus_three_div_four, + *deg_minus_three_div_two_plus_one, + ), + SqrtPrecomputation::AtkinCase5Mod8 { + trace, + char_minus_five_div_eight, + deg_minus_three_div_two_plus_one, + } => atkin( + elem, + trace, + char_minus_five_div_eight, + *deg_minus_three_div_two_plus_one, + ), + SqrtPrecomputation::KongCase9Mod16 { + trace, + c, + d, + c_squared, + char_minus_nine_div_sixteen, + deg_minus_three_div_two_plus_one, + } => kong( + elem, + trace, + c, + d, + c_squared, + char_minus_nine_div_sixteen, + *deg_minus_three_div_two_plus_one, + ), Self::PowerCase3Mod4 { modulus_plus_one_div_four, } => power_case_three_mod_four(elem, modulus_plus_one_div_four), @@ -108,7 +176,7 @@ fn tonelli_shanks( elem: &F, two_adicity: &u32, quadratic_nonresidue_to_trace: &F, - trace_of_modulus_minus_one_div_two: &'static [u64], + trace_of_modulus_minus_one_div_two: &[u64], ) -> Option { // Actually this is just normal Tonelli-Shanks; since `P::Generator` // is a quadratic non-residue, `P::ROOT_OF_UNITY = P::GENERATOR ^ t` @@ -163,13 +231,9 @@ fn tonelli_shanks( } } -/// https://eprint.iacr.org/2012/685.pdf (page 9, algorithm 2). -/// With _q_ as field order, _p_ as characteristic, and _m_ as extension degree: -/// * `char_minus_three_div_four` - _(p - 3)/4_. -/// * `deg_minus_three_div_two_plus_one` - _(m - 3)/2 + 1_. fn shanks( elem: &F, - char_minus_three_div_four: &'static [u64], + char_minus_three_div_four: &[u64], deg_minus_three_div_two_plus_one: usize, ) -> Option { // Computing a1 = Using decomposition of (q-3)/4 = a + p[pa + (3a+2)] * sum_i=1^(m-3)/2 p^2i @@ -198,15 +262,10 @@ fn shanks( Some(a1_elem) } -/// https://eprint.iacr.org/2012/685.pdf (page 10, algorithm 3). -/// With _q_ as field order, _p_ as characteristic, and _m_ as extension degree: -/// * `trace` - _2^(q - 5)/8_. -/// * `char_minus_five_div_eight` - _(p - 5)/8_. -/// * `deg_minus_three_div_two_plus_one` - _(m - 3)/2 + 1_. fn atkin( elem: &F, trace: &F, - char_minus_five_div_eight: &'static [u64], + char_minus_five_div_eight: &[u64], deg_minus_three_div_two_plus_one: usize, ) -> Option { // Computing a1 = elem^(q-5)/8 using decomposition of @@ -240,21 +299,13 @@ fn atkin( Some(x) } -/// https://eprint.iacr.org/2012/685.pdf (page 11, algorithm 4). -/// With _q_ as field order, _p_ as characteristic, and _m_ as extension degree: -/// * `trace` - _2^(q - 9)/16_. -/// * `c` - nonzero value such that _chi_q(c) != 1_. -/// * `d` - _c^(q - 9)/8_. -/// * `c_squared` - _c^2_. -/// * `char_minus_nine_div_sixteen` - _(p - 9)/16_. -/// * `deg_minus_three_div_two_plus_one` - _(m - 3)/2 + 1_. fn kong( elem: &F, trace: &F, c: &F, d: &F, c_squared: &F, - char_minus_nine_div_sixteen: &'static [u64], + char_minus_nine_div_sixteen: &[u64], deg_minus_three_div_two_plus_one: usize, ) -> Option { // Using decomposition of (q-9)/16 = a + p[pa + (9a+5)] * sum_i=1^(m-3)/2 p^2i @@ -296,7 +347,7 @@ fn kong( fn power_case_three_mod_four( elem: &F, - modulus_plus_one_div_four: &&'static [u64], + modulus_plus_one_div_four: &[u64], ) -> Option { let result = elem.pow(modulus_plus_one_div_four.as_ref()); (result.square() == *elem).then_some(result) From 00447872da433b546077cab7a9ae3ba590971920 Mon Sep 17 00:00:00 2001 From: Alexander Wu Date: Sun, 15 Jan 2023 17:54:09 -0800 Subject: [PATCH 7/8] Fix typo in docs Co-authored-by: Pratyush Mishra --- ff/src/fields/sqrt.rs | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/ff/src/fields/sqrt.rs b/ff/src/fields/sqrt.rs index 778893246..5ff43aad4 100644 --- a/ff/src/fields/sqrt.rs +++ b/ff/src/fields/sqrt.rs @@ -71,7 +71,7 @@ pub enum SqrtPrecomputation { /// * First factor _q - 1 = 2^s t_ where _t_ is odd. /// * `two_adicity` - _s_. /// * `quadratic_nonresidue_to_trace` - _c^t_, with random _c_ such that _c^2^(s - 1) = 1_. - /// * `trace_of_modulus_minus_one_div_two` - _(t - 1)/2 + 1_. + /// * `trace_of_modulus_minus_one_div_two` - _(t - 1)/2_. TonelliShanks { two_adicity: u32, quadratic_nonresidue_to_trace: F, From 06b51ffaa92bbede861c8779350717c971054688 Mon Sep 17 00:00:00 2001 From: Alexander Wu Date: Sun, 19 Feb 2023 22:02:33 -0800 Subject: [PATCH 8/8] Remove precompute parameters for deg_minus_three_div_two_plus_one --- ff/src/fields/sqrt.rs | 61 ++++++++++--------------------------------- 1 file changed, 14 insertions(+), 47 deletions(-) diff --git a/ff/src/fields/sqrt.rs b/ff/src/fields/sqrt.rs index 5ff43aad4..7912abb44 100644 --- a/ff/src/fields/sqrt.rs +++ b/ff/src/fields/sqrt.rs @@ -80,20 +80,16 @@ pub enum SqrtPrecomputation { /// https://eprint.iacr.org/2012/685.pdf (page 9, algorithm 2). /// With _q_ as field order, _p_ as characteristic, and _m_ as extension degree: /// * `char_minus_three_div_four` - _(p - 3)/4_. - /// * `deg_minus_three_div_two_plus_one` - _(m - 3)/2 + 1_. ShanksCase3Mod4 { char_minus_three_div_four: &'static [u64], - deg_minus_three_div_two_plus_one: usize, }, /// https://eprint.iacr.org/2012/685.pdf (page 10, algorithm 3). /// With _q_ as field order, _p_ as characteristic, and _m_ as extension degree: /// * `trace` - _2^(q - 5)/8_. /// * `char_minus_five_div_eight` - _(p - 5)/8_. - /// * `deg_minus_three_div_two_plus_one` - _(m - 3)/2 + 1_. AtkinCase5Mod8 { trace: F, char_minus_five_div_eight: &'static [u64], - deg_minus_three_div_two_plus_one: usize, }, /// https://eprint.iacr.org/2012/685.pdf (page 11, algorithm 4). /// With _q_ as field order, _p_ as characteristic, and _m_ as extension degree: @@ -102,14 +98,12 @@ pub enum SqrtPrecomputation { /// * `d` - _c^(q - 9)/8_. /// * `c_squared` - _c^2_. /// * `char_minus_nine_div_sixteen` - _(p - 9)/16_. - /// * `deg_minus_three_div_two_plus_one` - _(m - 3)/2 + 1_. KongCase9Mod16 { trace: F, c: F, d: F, c_squared: F, char_minus_nine_div_sixteen: &'static [u64], - deg_minus_three_div_two_plus_one: usize, }, /// In the case of 3 mod 4, we can find the square root via an exponentiation, /// sqrt(a) = a^(p+1)/4. This can be proved using Euler's criterion, a^(p-1)/2 = 1 mod p. @@ -133,38 +127,18 @@ impl SqrtPrecomputation { ), SqrtPrecomputation::ShanksCase3Mod4 { char_minus_three_div_four, - deg_minus_three_div_two_plus_one, - } => shanks( - elem, - char_minus_three_div_four, - *deg_minus_three_div_two_plus_one, - ), + } => shanks(elem, char_minus_three_div_four), SqrtPrecomputation::AtkinCase5Mod8 { trace, char_minus_five_div_eight, - deg_minus_three_div_two_plus_one, - } => atkin( - elem, - trace, - char_minus_five_div_eight, - *deg_minus_three_div_two_plus_one, - ), + } => atkin(elem, trace, char_minus_five_div_eight), SqrtPrecomputation::KongCase9Mod16 { trace, c, d, c_squared, char_minus_nine_div_sixteen, - deg_minus_three_div_two_plus_one, - } => kong( - elem, - trace, - c, - d, - c_squared, - char_minus_nine_div_sixteen, - *deg_minus_three_div_two_plus_one, - ), + } => kong(elem, trace, c, d, c_squared, char_minus_nine_div_sixteen), Self::PowerCase3Mod4 { modulus_plus_one_div_four, } => power_case_three_mod_four(elem, modulus_plus_one_div_four), @@ -231,11 +205,7 @@ fn tonelli_shanks( } } -fn shanks( - elem: &F, - char_minus_three_div_four: &[u64], - deg_minus_three_div_two_plus_one: usize, -) -> Option { +fn shanks(elem: &F, char_minus_three_div_four: &[u64]) -> Option { // Computing a1 = Using decomposition of (q-3)/4 = a + p[pa + (3a+2)] * sum_i=1^(m-3)/2 p^2i // where a = (p - 3) / 4. // factor1 = elem^a @@ -248,8 +218,9 @@ fn shanks( * elem_to_p.square(); // factor2 = prod_i=1^(m-3)/2 factor2_base^(p^2i) let mut factor2 = F::one(); - for i in 1..deg_minus_three_div_two_plus_one { - factor2 *= factor2_base.frobenius_map(i * 2 as usize); + let n = (F::extension_degree() as usize - 3) / 2; + for i in 1..(n + 1) { + factor2 *= factor2_base.frobenius_map(i * 2); } let a1 = factor1 * factor2; @@ -262,12 +233,7 @@ fn shanks( Some(a1_elem) } -fn atkin( - elem: &F, - trace: &F, - char_minus_five_div_eight: &[u64], - deg_minus_three_div_two_plus_one: usize, -) -> Option { +fn atkin(elem: &F, trace: &F, char_minus_five_div_eight: &[u64]) -> Option { // Computing a1 = elem^(q-5)/8 using decomposition of // (q-5)/8 = a + p[pa + (5a+3)] * sum_i=1^(m-3)/2 p^2i // where a = (p - 5) / 8. @@ -281,8 +247,9 @@ fn atkin( * elem_to_p.pow(&[3u64]); // factor2 = prod_i=1^(m-3)/2 factor2_base^(p^2i) let mut factor2 = F::one(); - for i in 1..deg_minus_three_div_two_plus_one { - factor2 *= factor2_base.frobenius_map(2 * i); + let n = (F::extension_degree() as usize - 3) / 2; + for i in 1..(n + 1) { + factor2 *= factor2_base.frobenius_map(i * 2); } let a1 = factor1 * factor2; @@ -306,7 +273,6 @@ fn kong( d: &F, c_squared: &F, char_minus_nine_div_sixteen: &[u64], - deg_minus_three_div_two_plus_one: usize, ) -> Option { // Using decomposition of (q-9)/16 = a + p[pa + (9a+5)] * sum_i=1^(m-3)/2 p^2i // a = (p - 9) / 16 @@ -320,8 +286,9 @@ fn kong( * elem_to_p.pow(&[5u64]); // factor2 = prod_i=1^(m-3)/2 factor2_base^(p^2i) let mut factor2 = F::one(); - for i in 1..deg_minus_three_div_two_plus_one { - factor2 *= factor2_base.frobenius_map(2 * i); + let n = (F::extension_degree() as usize - 3) / 2; + for i in 1..(n + 1) { + factor2 *= factor2_base.frobenius_map(i * 2); } let a1 = factor1 * factor2;