Skip to main content

crypto_bigint/int/
gcd.rs

1//! This module implements (a constant variant of) the Optimized Extended Binary GCD algorithm,
2//! which is described by Pornin in "Optimized Binary GCD for Modular Inversion".
3//! Ref: <https://eprint.iacr.org/2020/972.pdf>
4
5use crate::primitives::u32_min;
6use crate::uint::gcd::{OddUintXgcdOutput, impl_gcd_unsigned_lhs, impl_gcd_unsigned_rhs};
7use crate::{Choice, Gcd, Int, NonZero, NonZeroInt, NonZeroUint, Odd, OddInt, OddUint, Uint, Xgcd};
8
9impl<const LIMBS: usize> Int<LIMBS> {
10    /// Compute the greatest common divisor of `self` and `rhs`.
11    #[must_use]
12    pub const fn gcd_unsigned(&self, rhs: &Uint<LIMBS>) -> Uint<LIMBS> {
13        self.abs().gcd(rhs)
14    }
15
16    /// Compute the greatest common divisor of `self` and `rhs`.
17    ///
18    /// Executes in variable time w.r.t. all input parameters.
19    #[must_use]
20    pub const fn gcd_unsigned_vartime(&self, rhs: &Uint<LIMBS>) -> Uint<LIMBS> {
21        self.abs().gcd_vartime(rhs)
22    }
23
24    /// Executes the Extended GCD algorithm.
25    ///
26    /// Given `(self, rhs)`, computes `(g, x, y)`, s.t. `self * x + rhs * y = g = gcd(self, rhs)`.
27    #[must_use]
28    pub const fn xgcd(&self, rhs: &Self) -> IntXgcdOutput<LIMBS> {
29        // Make sure `self` and `rhs` are nonzero.
30        let self_is_zero = self.is_nonzero().not();
31        let self_nz = Int::select(self, &Int::ONE, self_is_zero)
32            .to_nz()
33            .expect_copied("self is non zero by construction");
34        let rhs_is_zero = rhs.is_nonzero().not();
35        let rhs_nz = Int::select(rhs, &Int::ONE, rhs_is_zero)
36            .to_nz()
37            .expect_copied("rhs is non zero by construction");
38
39        let NonZeroIntXgcdOutput {
40            gcd,
41            mut x,
42            mut y,
43            mut lhs_on_gcd,
44            mut rhs_on_gcd,
45        } = self_nz.xgcd(&rhs_nz);
46
47        // Correct the gcd in case self and/or rhs was zero
48        let mut gcd = *gcd.as_ref();
49        gcd = Uint::select(&gcd, &rhs.abs(), self_is_zero);
50        gcd = Uint::select(&gcd, &self.abs(), rhs_is_zero);
51
52        // Correct the Bézout coefficients in case self and/or rhs was zero.
53        let signum_self =
54            Int::new_from_abs_sign(Uint::ONE, self.is_negative()).expect_copied("+/- 1");
55        let signum_rhs =
56            Int::new_from_abs_sign(Uint::ONE, rhs.is_negative()).expect_copied("+/- 1");
57        x = Int::select(&x, &Int::ZERO, self_is_zero);
58        y = Int::select(&y, &signum_rhs, self_is_zero);
59        x = Int::select(&x, &signum_self, rhs_is_zero);
60        y = Int::select(&y, &Int::ZERO, rhs_is_zero);
61
62        // Correct the quotients in case self and/or rhs was zero.
63        lhs_on_gcd = Int::select(&lhs_on_gcd, &signum_self, rhs_is_zero);
64        lhs_on_gcd = Int::select(&lhs_on_gcd, &Int::ZERO, self_is_zero);
65        rhs_on_gcd = Int::select(&rhs_on_gcd, &signum_rhs, self_is_zero);
66        rhs_on_gcd = Int::select(&rhs_on_gcd, &Int::ZERO, rhs_is_zero);
67
68        IntXgcdOutput {
69            gcd,
70            x,
71            y,
72            lhs_on_gcd,
73            rhs_on_gcd,
74        }
75    }
76}
77
78impl<const LIMBS: usize> NonZero<Int<LIMBS>> {
79    /// Compute the greatest common divisor of `self` and `rhs`.
80    #[must_use]
81    pub const fn gcd_unsigned(&self, rhs: &Uint<LIMBS>) -> NonZero<Uint<LIMBS>> {
82        self.abs().gcd_unsigned(rhs)
83    }
84
85    /// Compute the greatest common divisor of `self` and `rhs`.
86    ///
87    /// Executes in variable time w.r.t. all input parameters.
88    #[must_use]
89    pub const fn gcd_unsigned_vartime(&self, rhs: &Uint<LIMBS>) -> NonZeroUint<LIMBS> {
90        self.abs().gcd_unsigned_vartime(rhs)
91    }
92
93    /// Execute the Extended GCD algorithm.
94    ///
95    /// Given `(self, rhs)`, computes `(g, x, y)` s.t. `self * x + rhs * y = g = gcd(self, rhs)`.
96    #[must_use]
97    pub const fn xgcd(&self, rhs: &Self) -> NonZeroIntXgcdOutput<LIMBS> {
98        let (mut lhs, mut rhs) = (*self.as_ref(), *rhs.as_ref());
99
100        // Leverage the property that gcd(2^k * a, 2^k *b) = 2^k * gcd(a, b)
101        let i = lhs.0.trailing_zeros();
102        let j = rhs.0.trailing_zeros();
103        let k = u32_min(i, j);
104        lhs = lhs.shr(k);
105        rhs = rhs.shr(k);
106
107        // Note: at this point, either lhs or rhs is odd (or both).
108        // Swap to make sure lhs is odd.
109        let swap = Choice::from_u32_lt(j, i);
110        Int::conditional_swap(&mut lhs, &mut rhs, swap);
111        let lhs = lhs.to_odd().expect_copied("odd by construction");
112        let rhs = rhs.to_nz().expect_copied("non-zero by construction");
113
114        let OddIntXgcdOutput {
115            gcd,
116            mut x,
117            mut y,
118            mut lhs_on_gcd,
119            mut rhs_on_gcd,
120        } = lhs.xgcd(&rhs);
121
122        // Account for the parameter swap
123        Int::conditional_swap(&mut x, &mut y, swap);
124        Int::conditional_swap(&mut lhs_on_gcd, &mut rhs_on_gcd, swap);
125
126        // Reintroduce the factor 2^k to the gcd.
127        let gcd = gcd
128            .as_ref()
129            .shl(k)
130            .to_nz()
131            .expect_copied("is non-zero by construction");
132
133        NonZeroIntXgcdOutput {
134            gcd,
135            x,
136            y,
137            lhs_on_gcd,
138            rhs_on_gcd,
139        }
140    }
141}
142
143impl<const LIMBS: usize> Odd<Int<LIMBS>> {
144    /// Compute the greatest common divisor of `self` and `rhs`.
145    #[must_use]
146    pub const fn gcd_unsigned(&self, rhs: &Uint<LIMBS>) -> Odd<Uint<LIMBS>> {
147        self.abs().gcd_unsigned(rhs)
148    }
149
150    /// Compute the greatest common divisor of `self` and `rhs`.
151    ///
152    /// Executes in variable time w.r.t. all input parameters.
153    #[must_use]
154    pub const fn gcd_unsigned_vartime(&self, rhs: &Uint<LIMBS>) -> OddUint<LIMBS> {
155        self.abs().gcd_unsigned_vartime(rhs)
156    }
157
158    /// Execute the Extended GCD algorithm.
159    ///
160    /// Given `(self, rhs)`, computes `(g, x, y)` s.t. `self * x + rhs * y = g = gcd(self, rhs)`.
161    #[must_use]
162    pub const fn xgcd(&self, rhs: &NonZero<Int<LIMBS>>) -> OddIntXgcdOutput<LIMBS> {
163        let (abs_lhs, sgn_lhs) = self.abs_sign();
164        let (abs_rhs, sgn_rhs) = rhs.abs_sign();
165
166        let OddUintXgcdOutput {
167            gcd,
168            mut x,
169            mut y,
170            lhs_on_gcd: abs_lhs_on_gcd,
171            rhs_on_gcd: abs_rhs_on_gcd,
172        } = OddUintXgcdOutput::from_pattern_output(abs_lhs.binxgcd_nz(&abs_rhs));
173
174        x = x.wrapping_neg_if(sgn_lhs);
175        y = y.wrapping_neg_if(sgn_rhs);
176        let lhs_on_gcd =
177            Int::new_from_abs_sign(abs_lhs_on_gcd, sgn_lhs).expect_copied("no overflow");
178        let rhs_on_gcd =
179            Int::new_from_abs_sign(abs_rhs_on_gcd, sgn_rhs).expect_copied("no overflow");
180
181        OddIntXgcdOutput {
182            gcd,
183            x,
184            y,
185            lhs_on_gcd,
186            rhs_on_gcd,
187        }
188    }
189}
190
191/// Output of the Binary XGCD algorithm applied to two [`Int`]s.
192pub type IntXgcdOutput<const LIMBS: usize> = XgcdOutput<LIMBS, Uint<LIMBS>>;
193
194/// Output of the Binary XGCD algorithm applied to two [`NonZero<Int<LIMBS>>`]s.
195pub type NonZeroIntXgcdOutput<const LIMBS: usize> = XgcdOutput<LIMBS, NonZeroUint<LIMBS>>;
196
197/// Output of the Binary XGCD algorithm applied to two [`Odd<Int<LIMBS>>`]s.
198pub type OddIntXgcdOutput<const LIMBS: usize> = XgcdOutput<LIMBS, OddUint<LIMBS>>;
199
200#[derive(Debug, Copy, Clone)]
201pub struct XgcdOutput<const LIMBS: usize, T: Copy> {
202    pub gcd: T,
203    pub x: Int<LIMBS>,
204    pub y: Int<LIMBS>,
205    pub lhs_on_gcd: Int<LIMBS>,
206    pub rhs_on_gcd: Int<LIMBS>,
207}
208
209impl<const LIMBS: usize, T: Copy> XgcdOutput<LIMBS, T> {
210    /// Return the `gcd`.
211    pub fn gcd(&self) -> T {
212        self.gcd
213    }
214
215    /// Return the quotients `lhs.gcd` and `rhs/gcd`.
216    pub fn quotients(&self) -> (Int<LIMBS>, Int<LIMBS>) {
217        (self.lhs_on_gcd, self.rhs_on_gcd)
218    }
219
220    /// Return the Bézout coefficients `x` and `y` s.t. `lhs * x + rhs * y = gcd`.
221    pub fn bezout_coefficients(&self) -> (Int<LIMBS>, Int<LIMBS>) {
222        (self.x, self.y)
223    }
224}
225
226macro_rules! impl_int_gcd_abs_lhs {
227    ($slf:ty, [$($rhs:ty),+]) => {
228        $(
229            impl_int_gcd_abs_lhs!($slf, $rhs, $rhs);
230        )+
231    };
232    ($slf:ty, $rhs:ty, $out:ty) => {
233        impl<const LIMBS: usize> Gcd<$rhs> for $slf {
234            type Output = $out;
235
236            #[inline]
237            fn gcd(&self, rhs: &$rhs) -> Self::Output {
238                rhs.gcd_unsigned(&self.abs())
239            }
240
241            #[inline]
242            fn gcd_vartime(&self, rhs: &$rhs) -> Self::Output {
243                rhs.gcd_unsigned_vartime(&self.abs())
244            }
245        }
246    };
247}
248
249macro_rules! impl_int_gcd_abs_rhs {
250    ($slf:ty, [$($rhs:ty),+], $out:ty) => {
251        $(
252            impl_int_gcd_abs_rhs!($slf, $rhs, $out);
253        )+
254    };
255    ($slf:ty, $rhs:ty, $out:ty) => {
256        impl<const LIMBS: usize> Gcd<$rhs> for $slf {
257            type Output = $out;
258
259            #[inline]
260            fn gcd(&self, rhs: &$rhs) -> Self::Output {
261                self.gcd_unsigned(&rhs.abs())
262            }
263
264            #[inline]
265            fn gcd_vartime(&self, rhs: &$rhs) -> Self::Output {
266                self.gcd_unsigned_vartime(&rhs.abs())
267            }
268        }
269    }
270}
271
272// avoiding (NonZero|Odd)•(NonZero|Odd) combinations except for Self•Self to limit compilation time
273impl_gcd_unsigned_lhs!(Int<LIMBS>, Uint<LIMBS>, Uint<LIMBS>);
274impl_gcd_unsigned_lhs!(NonZeroInt<LIMBS>, Uint<LIMBS>, NonZeroUint<LIMBS>);
275impl_gcd_unsigned_lhs!(OddInt<LIMBS>, Uint<LIMBS>, OddUint<LIMBS>);
276impl_gcd_unsigned_rhs!(Uint<LIMBS>, Int<LIMBS>, Uint<LIMBS>);
277impl_gcd_unsigned_rhs!(Uint<LIMBS>, NonZeroInt<LIMBS>, NonZeroUint<LIMBS>);
278impl_gcd_unsigned_rhs!(Uint<LIMBS>, OddInt<LIMBS>, OddUint<LIMBS>);
279impl_int_gcd_abs_lhs!(Int<LIMBS>, Int<LIMBS>, Uint<LIMBS>);
280impl_int_gcd_abs_lhs!(Int<LIMBS>, [NonZeroUint<LIMBS>, OddUint<LIMBS>]);
281impl_int_gcd_abs_rhs!(NonZeroInt<LIMBS>, [Int<LIMBS>, NonZeroInt<LIMBS>], NonZeroUint<LIMBS>);
282impl_int_gcd_abs_rhs!(OddInt<LIMBS>, [Int<LIMBS>, OddInt<LIMBS>], OddUint<LIMBS>);
283
284impl<const LIMBS: usize> Xgcd for Int<LIMBS> {
285    type Output = IntXgcdOutput<LIMBS>;
286
287    fn xgcd(&self, rhs: &Int<LIMBS>) -> Self::Output {
288        self.xgcd(rhs)
289    }
290
291    fn xgcd_vartime(&self, rhs: &Int<LIMBS>) -> Self::Output {
292        // TODO(#853): implement vartime
293        self.xgcd(rhs)
294    }
295}
296
297impl<const LIMBS: usize> Xgcd for NonZeroInt<LIMBS> {
298    type Output = NonZeroIntXgcdOutput<LIMBS>;
299
300    fn xgcd(&self, rhs: &NonZeroInt<LIMBS>) -> Self::Output {
301        self.xgcd(rhs)
302    }
303
304    fn xgcd_vartime(&self, rhs: &NonZeroInt<LIMBS>) -> Self::Output {
305        // TODO(#853): implement vartime
306        self.xgcd(rhs)
307    }
308}
309
310impl<const LIMBS: usize> Xgcd for OddInt<LIMBS> {
311    type Output = OddIntXgcdOutput<LIMBS>;
312
313    fn xgcd(&self, rhs: &OddInt<LIMBS>) -> Self::Output {
314        self.xgcd(rhs.as_nz_ref())
315    }
316
317    fn xgcd_vartime(&self, rhs: &OddInt<LIMBS>) -> Self::Output {
318        // TODO(#853): implement vartime
319        self.xgcd(rhs.as_nz_ref())
320    }
321}
322
323#[cfg(all(test, not(miri)))]
324mod tests {
325    use crate::int::gcd::{IntXgcdOutput, NonZeroIntXgcdOutput, OddIntXgcdOutput};
326    use crate::{Concat, Gcd, Int, Uint};
327
328    impl<const LIMBS: usize> From<NonZeroIntXgcdOutput<LIMBS>> for IntXgcdOutput<LIMBS> {
329        fn from(value: NonZeroIntXgcdOutput<LIMBS>) -> Self {
330            let NonZeroIntXgcdOutput {
331                gcd,
332                x,
333                y,
334                lhs_on_gcd,
335                rhs_on_gcd,
336            } = value;
337            IntXgcdOutput {
338                gcd: *gcd.as_ref(),
339                x,
340                y,
341                lhs_on_gcd,
342                rhs_on_gcd,
343            }
344        }
345    }
346
347    impl<const LIMBS: usize> From<OddIntXgcdOutput<LIMBS>> for IntXgcdOutput<LIMBS> {
348        fn from(value: OddIntXgcdOutput<LIMBS>) -> Self {
349            let OddIntXgcdOutput {
350                gcd,
351                x,
352                y,
353                lhs_on_gcd,
354                rhs_on_gcd,
355            } = value;
356            IntXgcdOutput {
357                gcd: *gcd.as_ref(),
358                x,
359                y,
360                lhs_on_gcd,
361                rhs_on_gcd,
362            }
363        }
364    }
365
366    mod gcd {
367
368        use crate::{Gcd, I64, I128, I256, I512, I1024, I2048, I4096, Int, Uint};
369
370        fn test<const LIMBS: usize>(lhs: Int<LIMBS>, rhs: Int<LIMBS>, target: Uint<LIMBS>) {
371            assert_eq!(lhs.gcd(&rhs), target);
372            assert_eq!(lhs.gcd_vartime(&rhs), target);
373        }
374
375        fn run_tests<const LIMBS: usize>() {
376            let abs_min = *Int::MIN.as_uint();
377            let max = *Int::MAX.as_uint();
378            test(Int::<LIMBS>::MIN, Int::MIN, abs_min);
379            test(Int::<LIMBS>::MIN, Int::MINUS_ONE, Uint::ONE);
380            test(Int::<LIMBS>::MIN, Int::ZERO, abs_min);
381            test(Int::<LIMBS>::MIN, Int::ONE, Uint::ONE);
382            test(Int::<LIMBS>::MIN, Int::MAX, Uint::ONE);
383            test(Int::<LIMBS>::MINUS_ONE, Int::MIN, Uint::ONE);
384            test(Int::<LIMBS>::MINUS_ONE, Int::MINUS_ONE, Uint::ONE);
385            test(Int::<LIMBS>::MINUS_ONE, Int::ZERO, Uint::ONE);
386            test(Int::<LIMBS>::MINUS_ONE, Int::ONE, Uint::ONE);
387            test(Int::<LIMBS>::MINUS_ONE, Int::MAX, Uint::ONE);
388            test(Int::<LIMBS>::ZERO, Int::MIN, abs_min);
389            test(Int::<LIMBS>::ZERO, Int::MINUS_ONE, Uint::ONE);
390            test(Int::<LIMBS>::ZERO, Int::ZERO, Uint::ZERO);
391            test(Int::<LIMBS>::ZERO, Int::ONE, Uint::ONE);
392            test(Int::<LIMBS>::ZERO, Int::MAX, max);
393            test(Int::<LIMBS>::ONE, Int::MIN, Uint::ONE);
394            test(Int::<LIMBS>::ONE, Int::MINUS_ONE, Uint::ONE);
395            test(Int::<LIMBS>::ONE, Int::ZERO, Uint::ONE);
396            test(Int::<LIMBS>::ONE, Int::ONE, Uint::ONE);
397            test(Int::<LIMBS>::ONE, Int::MAX, Uint::ONE);
398            test(Int::<LIMBS>::MAX, Int::MIN, Uint::ONE);
399            test(Int::<LIMBS>::MAX, Int::MINUS_ONE, Uint::ONE);
400            test(Int::<LIMBS>::MAX, Int::ZERO, max);
401            test(Int::<LIMBS>::MAX, Int::ONE, Uint::ONE);
402            test(Int::<LIMBS>::MAX, Int::MAX, max);
403        }
404
405        #[test]
406        fn gcd_sizes() {
407            run_tests::<{ I64::LIMBS }>();
408            run_tests::<{ I128::LIMBS }>();
409            run_tests::<{ I256::LIMBS }>();
410            run_tests::<{ I512::LIMBS }>();
411            run_tests::<{ I1024::LIMBS }>();
412            run_tests::<{ I2048::LIMBS }>();
413            run_tests::<{ I4096::LIMBS }>();
414        }
415    }
416
417    fn xgcd_test<const LIMBS: usize, const DOUBLE: usize>(
418        lhs: Int<LIMBS>,
419        rhs: Int<LIMBS>,
420        output: IntXgcdOutput<LIMBS>,
421    ) where
422        Uint<LIMBS>: Concat<LIMBS, Output = Uint<DOUBLE>>,
423    {
424        let gcd = lhs.gcd(&rhs);
425        assert_eq!(gcd, output.gcd);
426
427        // Test quotients
428        let (lhs_on_gcd, rhs_on_gcd) = output.quotients();
429
430        if gcd.is_zero().to_bool() {
431            assert_eq!(lhs_on_gcd, Int::ZERO);
432            assert_eq!(rhs_on_gcd, Int::ZERO);
433        } else {
434            assert_eq!(lhs_on_gcd, lhs.div_unsigned(&gcd.to_nz().unwrap()));
435            assert_eq!(rhs_on_gcd, rhs.div_unsigned(&gcd.to_nz().unwrap()));
436        }
437
438        // Test the Bezout coefficients on minimality
439        let (x, y) = output.bezout_coefficients();
440        assert!(x.abs() <= rhs_on_gcd.abs() || rhs_on_gcd.is_zero().to_bool());
441        assert!(y.abs() <= lhs_on_gcd.abs() || lhs_on_gcd.is_zero().to_bool());
442        if lhs.abs() != rhs.abs() {
443            assert!(x.abs() <= rhs_on_gcd.abs().shr(1) || rhs_on_gcd.is_zero().to_bool());
444            assert!(y.abs() <= lhs_on_gcd.abs().shr(1) || lhs_on_gcd.is_zero().to_bool());
445        }
446
447        // Test the Bezout coefficients for correctness
448        assert_eq!(
449            x.concatenating_mul(&lhs)
450                .wrapping_add(&y.concatenating_mul(&rhs)),
451            *gcd.resize().as_int()
452        );
453    }
454
455    mod test_int_xgcd {
456        use crate::int::gcd::tests::xgcd_test;
457        use crate::{
458            Concat, Gcd, Int, U64, U128, U192, U256, U384, U512, U768, U1024, U2048, U4096, U8192,
459            Uint,
460        };
461
462        fn test<const LIMBS: usize, const DOUBLE: usize>(lhs: Int<LIMBS>, rhs: Int<LIMBS>)
463        where
464            Uint<LIMBS>: Concat<LIMBS, Output = Uint<DOUBLE>>,
465            Int<LIMBS>: Gcd<Output = Uint<LIMBS>>,
466        {
467            xgcd_test(lhs, rhs, lhs.xgcd(&rhs));
468        }
469
470        fn run_tests<const LIMBS: usize, const DOUBLE: usize>()
471        where
472            Uint<LIMBS>: Concat<LIMBS, Output = Uint<DOUBLE>>,
473            Int<LIMBS>: Gcd<Output = Uint<LIMBS>>,
474        {
475            test(Int::MIN, Int::MIN);
476            test(Int::MIN, Int::MINUS_ONE);
477            test(Int::MIN, Int::ZERO);
478            test(Int::MIN, Int::ONE);
479            test(Int::MIN, Int::MAX);
480            test(Int::ONE, Int::MIN);
481            test(Int::ONE, Int::MINUS_ONE);
482            test(Int::ONE, Int::ZERO);
483            test(Int::ONE, Int::ONE);
484            test(Int::ONE, Int::MAX);
485            test(Int::ZERO, Int::MIN);
486            test(Int::ZERO, Int::MINUS_ONE);
487            test(Int::ZERO, Int::ZERO);
488            test(Int::ZERO, Int::ONE);
489            test(Int::ZERO, Int::MAX);
490            test(Int::ONE, Int::MIN);
491            test(Int::ONE, Int::MINUS_ONE);
492            test(Int::ONE, Int::ZERO);
493            test(Int::ONE, Int::ONE);
494            test(Int::ONE, Int::MAX);
495            test(Int::MAX, Int::MIN);
496            test(Int::MAX, Int::MINUS_ONE);
497            test(Int::MAX, Int::ZERO);
498            test(Int::MAX, Int::ONE);
499            test(Int::MAX, Int::MAX);
500        }
501
502        #[test]
503        fn xgcd_sizes() {
504            run_tests::<{ U64::LIMBS }, { U128::LIMBS }>();
505            run_tests::<{ U128::LIMBS }, { U256::LIMBS }>();
506            run_tests::<{ U192::LIMBS }, { U384::LIMBS }>();
507            run_tests::<{ U256::LIMBS }, { U512::LIMBS }>();
508            run_tests::<{ U384::LIMBS }, { U768::LIMBS }>();
509            run_tests::<{ U512::LIMBS }, { U1024::LIMBS }>();
510            run_tests::<{ U1024::LIMBS }, { U2048::LIMBS }>();
511            run_tests::<{ U2048::LIMBS }, { U4096::LIMBS }>();
512            run_tests::<{ U4096::LIMBS }, { U8192::LIMBS }>();
513        }
514    }
515
516    mod test_nonzero_int_xgcd {
517        use crate::int::gcd::tests::xgcd_test;
518        use crate::{
519            Concat, Int, U64, U128, U192, U256, U384, U512, U768, U1024, U2048, U4096, U8192, Uint,
520        };
521
522        fn test<const LIMBS: usize, const DOUBLE: usize>(lhs: Int<LIMBS>, rhs: Int<LIMBS>)
523        where
524            Uint<LIMBS>: Concat<LIMBS, Output = Uint<DOUBLE>>,
525        {
526            let output = lhs.to_nz().unwrap().xgcd(&rhs.to_nz().unwrap());
527            xgcd_test(lhs, rhs, output.into());
528        }
529
530        fn run_tests<const LIMBS: usize, const DOUBLE: usize>()
531        where
532            Uint<LIMBS>: Concat<LIMBS, Output = Uint<DOUBLE>>,
533        {
534            test(Int::MIN, Int::MIN);
535            test(Int::MIN, Int::MINUS_ONE);
536            test(Int::MIN, Int::ONE);
537            test(Int::MIN, Int::MAX);
538            test(Int::MINUS_ONE, Int::MIN);
539            test(Int::MINUS_ONE, Int::MINUS_ONE);
540            test(Int::MINUS_ONE, Int::ONE);
541            test(Int::MINUS_ONE, Int::MAX);
542            test(Int::ONE, Int::MIN);
543            test(Int::ONE, Int::MINUS_ONE);
544            test(Int::ONE, Int::ONE);
545            test(Int::ONE, Int::MAX);
546            test(Int::MAX, Int::MIN);
547            test(Int::MAX, Int::MINUS_ONE);
548            test(Int::MAX, Int::ONE);
549            test(Int::MAX, Int::MAX);
550        }
551
552        #[test]
553        fn binxgcd() {
554            run_tests::<{ U64::LIMBS }, { U128::LIMBS }>();
555            run_tests::<{ U128::LIMBS }, { U256::LIMBS }>();
556            run_tests::<{ U192::LIMBS }, { U384::LIMBS }>();
557            run_tests::<{ U256::LIMBS }, { U512::LIMBS }>();
558            run_tests::<{ U384::LIMBS }, { U768::LIMBS }>();
559            run_tests::<{ U512::LIMBS }, { U1024::LIMBS }>();
560            run_tests::<{ U1024::LIMBS }, { U2048::LIMBS }>();
561            run_tests::<{ U2048::LIMBS }, { U4096::LIMBS }>();
562            run_tests::<{ U4096::LIMBS }, { U8192::LIMBS }>();
563        }
564    }
565
566    mod test_odd_int_xgcd {
567        use crate::int::gcd::tests::xgcd_test;
568        use crate::{
569            Concat, Int, U64, U128, U192, U256, U384, U512, U768, U1024, U2048, U4096, U8192, Uint,
570        };
571
572        fn test<const LIMBS: usize, const DOUBLE: usize>(lhs: Int<LIMBS>, rhs: Int<LIMBS>)
573        where
574            Uint<LIMBS>: Concat<LIMBS, Output = Uint<DOUBLE>>,
575        {
576            let output = lhs.to_odd().unwrap().xgcd(&rhs.to_nz().unwrap());
577            xgcd_test(lhs, rhs, output.into());
578        }
579
580        fn run_tests<const LIMBS: usize, const DOUBLE: usize>()
581        where
582            Uint<LIMBS>: Concat<LIMBS, Output = Uint<DOUBLE>>,
583        {
584            let neg_max = Int::MAX.wrapping_neg();
585            test(neg_max, neg_max);
586            test(neg_max, Int::MINUS_ONE);
587            test(neg_max, Int::ONE);
588            test(neg_max, Int::MAX);
589            test(Int::ONE, neg_max);
590            test(Int::ONE, Int::MINUS_ONE);
591            test(Int::ONE, Int::ONE);
592            test(Int::ONE, Int::MAX);
593            test(Int::MAX, neg_max);
594            test(Int::MAX, Int::MINUS_ONE);
595            test(Int::MAX, Int::ONE);
596            test(Int::MAX, Int::MAX);
597        }
598
599        #[test]
600        fn binxgcd() {
601            run_tests::<{ U64::LIMBS }, { U128::LIMBS }>();
602            run_tests::<{ U128::LIMBS }, { U256::LIMBS }>();
603            run_tests::<{ U192::LIMBS }, { U384::LIMBS }>();
604            run_tests::<{ U256::LIMBS }, { U512::LIMBS }>();
605            run_tests::<{ U384::LIMBS }, { U768::LIMBS }>();
606            run_tests::<{ U512::LIMBS }, { U1024::LIMBS }>();
607            run_tests::<{ U1024::LIMBS }, { U2048::LIMBS }>();
608            run_tests::<{ U2048::LIMBS }, { U4096::LIMBS }>();
609            run_tests::<{ U4096::LIMBS }, { U8192::LIMBS }>();
610        }
611    }
612
613    mod traits {
614        use crate::{Gcd, I256, U256};
615
616        #[test]
617        fn gcd_always_positive() {
618            // Two numbers with a shared factor of 61
619            let f = I256::from(59i32 * 61);
620            let g = I256::from(61i32 * 71);
621
622            assert_eq!(U256::from(61u32), f.gcd(&g));
623            assert_eq!(U256::from(61u32), f.wrapping_neg().gcd(&g));
624            assert_eq!(U256::from(61u32), f.gcd(&g.wrapping_neg()));
625            assert_eq!(U256::from(61u32), f.wrapping_neg().gcd(&g.wrapping_neg()));
626        }
627
628        #[test]
629        fn gcd_int_uint() {
630            // Two numbers with a shared factor of 61
631            let f = I256::from(59i32 * 61);
632            let g = U256::from(61u32 * 71);
633
634            assert_eq!(U256::from(61u32), f.gcd(&g));
635            assert_eq!(U256::from(61u32), f.wrapping_neg().gcd(&g));
636        }
637    }
638}