1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
//! Full-precision Euclidean step on the GCD state.
use super::GcdState;
impl GcdState {
/// Performs one full-precision Euclidean step.
///
/// This is the fallback path used when the quotient `r0 / r1` doesn't
/// fit in u64 (i.e. when the top-bit approximation `b == 0`). Computes
/// `q = r0 / r1` using full-width `div_rem`, updates the remainder pair,
/// and advances the Bézout coefficient pair using wrapping 256-bit
/// arithmetic. Toggles the swap parity.
///
/// Deliberately not `#[inline]`: this is a cold path (called 1–2 times
/// per inversion) and inlining pulls `div_rem` + `mul_u256`
/// into the caller's LTO context, hurting LLVM's codegen for unrelated
/// callers.
pub fn full_step(&mut self) {
let (q, r2) = self.r0.div_rem(&self.r1);
self.r0 = self.r1;
self.r1 = r2;
let qt1 = self.mul_u256(&q);
let new_t = self.t0.overflowing_sub(&qt1).0;
self.t0 = self.t1;
self.t1 = new_t;
self.even = !self.even;
}
}
#[cfg(test)]
mod ai_tests {
use crate::u256::U256;
use super::*;
/// Single full step on small values: gcd(100, 37).
#[test]
fn small_values() {
let mut state = GcdState::new(
U256::from_be_limbs([0, 0, 0, 100]),
U256::from_be_limbs([0, 0, 0, 37]),
);
// q = 100/37 = 2, r2 = 100 - 2*37 = 26
state.full_step();
assert_eq!(state.r0, U256::from_be_limbs([0, 0, 0, 37]));
assert_eq!(state.r1, U256::from_be_limbs([0, 0, 0, 26]));
assert!(!state.even);
}
/// Parity toggles on each full step.
#[test]
fn parity_toggle() {
let mut state = GcdState::new(
U256::from_be_limbs([0, 0, 0, 100]),
U256::from_be_limbs([0, 0, 0, 37]),
);
assert!(state.even);
state.full_step();
assert!(!state.even);
state.full_step();
assert!(state.even);
}
/// Coefficient tracking: t0 and t1 advance correctly.
#[test]
fn coefficient_update() {
let mut state = GcdState::new(
U256::from_be_limbs([0, 0, 0, 100]),
U256::from_be_limbs([0, 0, 0, 37]),
);
// Initially t0=0, t1=1. After step: q=2.
// new_t = t0 - q*t1 = 0 - 2*1 = -2 (wrapping)
// t0' = old t1 = 1, t1' = -2 (wrapping)
state.full_step();
assert_eq!(state.t0, U256::ONE);
// -2 mod 2^256
let neg2 = U256::ZERO
.overflowing_sub(&U256::from_be_limbs([0, 0, 0, 2]))
.0;
assert_eq!(state.t1, neg2);
}
/// Full step terminates when r1 becomes zero.
#[test]
fn terminates_at_exact_divisor() {
let mut state = GcdState::new(
U256::from_be_limbs([0, 0, 0, 6]),
U256::from_be_limbs([0, 0, 0, 3]),
);
// q = 6/3 = 2, r2 = 0
state.full_step();
assert_eq!(state.r0, U256::from_be_limbs([0, 0, 0, 3]));
assert_eq!(state.r1, U256::ZERO);
}
}