@@ -18,7 +18,8 @@ pub fn ExtendedGreatestCommonDivisor(S: anytype) type {
1818inline fn egcd_helper (other : anytype , odd : anytype , shift : anytype ) [3 ]@TypeOf (other , odd ) {
1919 const S = @TypeOf (other , odd );
2020 const toinv = @shrExact (other , @intCast (shift ));
21- const ctrl = @shrExact (odd , @intCast (shift ));
21+ const ctrl = @shrExact (odd , @intCast (shift )); // Invariant: |s|, |t|, |ctrl| < |MIN_OF(S)|
22+ const hctrl = 1 + @shrExact (ctrl - 1 , 1 );
2223
2324 var s : S = std .math .sign (toinv );
2425 var t : S = 0 ;
@@ -29,8 +30,11 @@ inline fn egcd_helper(other: anytype, odd: anytype, shift: anytype) [3]@TypeOf(o
2930 {
3031 const xz = @ctz (x );
3132 x = @shrExact (x , @intCast (xz ));
32- for (0.. xz ) | _ |
33- s = @shrExact (if (s & 1 == 0 ) s else s + ctrl , 1 );
33+ for (0.. xz ) | _ | {
34+ const s_odd = s & 1 != 0 ;
35+ s = @divFloor (s , 2 );
36+ if (s_odd ) s += hctrl ;
37+ }
3438 }
3539
3640 var y_minus_x = y -% x ;
@@ -50,8 +54,11 @@ inline fn egcd_helper(other: anytype, odd: anytype, shift: anytype) [3]@TypeOf(o
5054 t = copy_s ;
5155 }
5256 x = @shrExact (x , @intCast (xz ));
53- for (0.. xz ) | _ |
54- s = @shrExact (if (s & 1 == 0 ) s else s + ctrl , 1 );
57+ for (0.. xz ) | _ | {
58+ const s_odd = s & 1 != 0 ;
59+ s = @divFloor (s , 2 );
60+ if (s_odd ) s += hctrl ;
61+ }
5562 }
5663
5764 y = @shlExact (y , @intCast (shift ));
@@ -96,6 +103,15 @@ pub fn egcd(a: anytype, b: anytype) ExtendedGreatestCommonDivisor(@TypeOf(a, b))
96103}
97104
98105test {
106+ {
107+ const a : i8 = 127 ;
108+ const b : i8 = 126 ;
109+ const r = egcd (a , b );
110+ const g = r .gcd ;
111+ const s : i16 = r .bezout_coeff_1 ;
112+ const t : i16 = r .bezout_coeff_2 ;
113+ try std .testing .expect (s * a + t * b == g );
114+ }
99115 {
100116 const a : i4 = -8 ;
101117 const b : i4 = 1 ;
0 commit comments