Skip to content

Instantly share code, notes, and snippets.

@ExpHP
Last active July 29, 2016 23:16
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save ExpHP/022f0cce07d443786a8c301494375b11 to your computer and use it in GitHub Desktop.
Save ExpHP/022f0cce07d443786a8c301494375b11 to your computer and use it in GitHub Desktop.
extern crate num_integer;
extern crate num_traits;
use num_integer::Integer;
use num_traits::{PrimInt,Signed,Zero,One};
#[derive(Copy,Clone,Debug,Eq,PartialEq)]
pub struct GcdData<X> {
// greatest common divisor
pub gcd: X,
// least common multiple
pub lcm: X,
// bezout coefficients
pub coeffs: (X, X),
// quotients of the inputs by the GCD
pub quotients: (X, X),
}
// NOTE:
// The Signed bound is unavoidable for Extended GCD because the Bezout
// coefficients can be negative. This is unfortunate for plain old gcd(),
// which technically shouldn't require the Signed bound.
// Since the bezout coefficients have no impact on each other or on the gcd,
// a sufficiently smart compiler can rip their computations out entirely.
// And as luck would have it, rustc is sufficiently smart!
#[inline(always)]
fn extended_gcd__inline<X>(a: X, b: X) -> GcdData<X> where
X: PrimInt + Integer + Signed,
{
let (a_sign, a) = (a.signum(), a.abs());
let (b_sign, b) = (b.signum(), b.abs());
// Silly trick because rust doesn't have true multiple assignment:
// Store the two targets in one variable! Order is (old, current).
let mut s = (X::one(), X::zero()); // a coefficient
let mut t = (X::zero(), X::one()); // b coefficient
let mut r = (a, b); // gcd
while r.1 != X::zero() {
let (div, rem) = (r.0/r.1, r.0%r.1);
r = (r.1, rem);
s = (s.1, s.0 - div * s.1);
t = (t.1, t.0 - div * t.1);
}
let quots = (a_sign * t.1.abs(), b_sign * s.1.abs());
GcdData {
gcd: r.0,
// FIXME think more about sign of LCM
// (current implementation preserves the property a*b == gcd*lcm
// which is nice, but I don't know if it is The Right Thing)
lcm: r.0*quots.0*quots.1,
coeffs: (a_sign*s.0, b_sign*t.0),
quotients: quots,
}
}
#[test]
mod tests {
// Gold standard for binary comparison.
#[inline(never)]
fn gcd__reference<X>(a: X, b: X) -> X where
X: PrimInt + Integer + Signed,
{
let mut a = a.abs();
let mut b = b.abs();
while b != X::zero() {
let tmp = b;
b = a % b;
a = tmp;
}
a
}
// Impressively, rustc grinds this down to a *byte-perfect match*
// against gcd__reference.
#[inline(never)]
fn gcd__optimized<X>(a: X, b: X) -> X where
X: PrimInt + Integer + Signed,
{
gcd(a, b)
}
}
Dump of assembler code for function _ZN9numtheory5tests14gcd__reference17hc1dcde6d8c1a542bE:
0x000000000000fd60 <+0>: push %rax
0x000000000000fd61 <+1>: mov %edi,%eax
0x000000000000fd63 <+3>: neg %eax
0x000000000000fd65 <+5>: cmovl %edi,%eax
0x000000000000fd68 <+8>: mov %esi,%edx
0x000000000000fd6a <+10>: neg %edx
0x000000000000fd6c <+12>: cmovl %esi,%edx
0x000000000000fd6f <+15>: jmp 0xfd85 <_ZN9numtheory5tests14gcd__reference17hc1dcde6d8c1a542bE+37>
0x000000000000fd71 <+17>: data32 data32 data32 data32 data32 nopw %cs:0x0(%rax,%rax,1)
0x000000000000fd80 <+32>: cltd
0x000000000000fd81 <+33>: idiv %ecx
0x000000000000fd83 <+35>: mov %ecx,%eax
0x000000000000fd85 <+37>: mov %edx,%ecx
0x000000000000fd87 <+39>: cmp $0xffffffff,%ecx
0x000000000000fd8a <+42>: je 0xfda0 <_ZN9numtheory5tests14gcd__reference17hc1dcde6d8c1a542bE+64>
0x000000000000fd8c <+44>: test %ecx,%ecx
0x000000000000fd8e <+46>: jne 0xfd80 <_ZN9numtheory5tests14gcd__reference17hc1dcde6d8c1a542bE+32>
0x000000000000fd90 <+48>: jmp 0xfdb3 <_ZN9numtheory5tests14gcd__reference17hc1dcde6d8c1a542bE+83>
0x000000000000fd92 <+50>: data32 data32 data32 data32 nopw %cs:0x0(%rax,%rax,1)
0x000000000000fda0 <+64>: cmp $0x80000000,%eax
0x000000000000fda5 <+69>: jne 0xfd80 <_ZN9numtheory5tests14gcd__reference17hc1dcde6d8c1a542bE+32>
0x000000000000fda7 <+71>: lea 0x2a456a(%rip),%rdi # 0x2b4318 <panic_loc9109>
0x000000000000fdae <+78>: callq 0x8e480 <_ZN4core9panicking5panic17hbfac80217e56ecbeE>
0x000000000000fdb3 <+83>: pop %rcx
0x000000000000fdb4 <+84>: retq
End of assembler dump.
Dump of assembler code for function _ZN9numtheory5tests14gcd__optimized17ha4ee88a7db73c270E:
0x000000000000fdc0 <+0>: push %rax
0x000000000000fdc1 <+1>: mov %edi,%eax
0x000000000000fdc3 <+3>: neg %eax
0x000000000000fdc5 <+5>: cmovl %edi,%eax
0x000000000000fdc8 <+8>: mov %esi,%edx
0x000000000000fdca <+10>: neg %edx
0x000000000000fdcc <+12>: cmovl %esi,%edx
0x000000000000fdcf <+15>: jmp 0xfde5 <_ZN9numtheory5tests14gcd__optimized17ha4ee88a7db73c270E+37>
0x000000000000fdd1 <+17>: data32 data32 data32 data32 data32 nopw %cs:0x0(%rax,%rax,1)
0x000000000000fde0 <+32>: cltd
0x000000000000fde1 <+33>: idiv %ecx
0x000000000000fde3 <+35>: mov %ecx,%eax
0x000000000000fde5 <+37>: mov %edx,%ecx
0x000000000000fde7 <+39>: cmp $0xffffffff,%ecx
0x000000000000fdea <+42>: je 0xfe00 <_ZN9numtheory5tests14gcd__optimized17ha4ee88a7db73c270E+64>
0x000000000000fdec <+44>: test %ecx,%ecx
0x000000000000fdee <+46>: jne 0xfde0 <_ZN9numtheory5tests14gcd__optimized17ha4ee88a7db73c270E+32>
0x000000000000fdf0 <+48>: jmp 0xfe13 <_ZN9numtheory5tests14gcd__optimized17ha4ee88a7db73c270E+83>
0x000000000000fdf2 <+50>: data32 data32 data32 data32 nopw %cs:0x0(%rax,%rax,1)
0x000000000000fe00 <+64>: cmp $0x80000000,%eax
0x000000000000fe05 <+69>: jne 0xfde0 <_ZN9numtheory5tests14gcd__optimized17ha4ee88a7db73c270E+32>
0x000000000000fe07 <+71>: lea 0x2a44e2(%rip),%rdi # 0x2b42f0 <panic_loc9105>
0x000000000000fe0e <+78>: callq 0x8e480 <_ZN4core9panicking5panic17hbfac80217e56ecbeE>
0x000000000000fe13 <+83>: pop %rcx
0x000000000000fe14 <+84>: retq
End of assembler dump.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment