Skip to content

Commit f14dabc

Browse files
authored
math.big: add a new greatest-common-divisor-algo for big.Integer, also add a benchmark for it (#12261)
1 parent f62b2dc commit f14dabc

File tree

5 files changed

+633
-0
lines changed

5 files changed

+633
-0
lines changed

.gitignore

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -20,6 +20,7 @@ vdbg.exe
2020
*.dll
2121
*.lib
2222
*.bak
23+
*.dylib
2324
a.out
2425
.noprefix.vrepl_temp
2526

vlib/math/big/integer.v

Lines changed: 68 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -761,3 +761,71 @@ pub fn (a Integer) isqrt() Integer {
761761
}
762762
return result
763763
}
764+
765+
[inline]
766+
fn bi_min(a Integer, b Integer) Integer {
767+
return if a < b { a } else { b }
768+
}
769+
770+
[inline]
771+
fn bi_max(a Integer, b Integer) Integer {
772+
return if a > b { a } else { b }
773+
}
774+
775+
[direct_array_access]
776+
fn (bi Integer) msb() u32 {
777+
for idx := 0; idx < bi.digits.len; idx += 1 {
778+
word := bi.digits[idx]
779+
if word > 0 {
780+
return u32((idx * 32) + bits.trailing_zeros_32(word))
781+
}
782+
}
783+
return u32(32)
784+
}
785+
786+
// Greatest-Common-Divisor https://en.wikipedia.org/wiki/Binary_GCD_algorithm
787+
// The code below follows the 2013-christmas-special by D. Lemire & R. Corderoy
788+
// https://en.algorithmica.org/hpc/analyzing-performance/gcd/
789+
//
790+
// discussion & further info https://lemire.me/blog/2013/12/26/fastest-way-to-compute-the-greatest-common-divisor/
791+
792+
pub fn (x Integer) gcd_binary(y Integer) Integer {
793+
// Since standard-euclid-gcd is much faster on smaller sizes 4-8-Byte.
794+
// In such a case, one could delegate back to big.Integer.gcd()
795+
// Uncomment below and a all long-long goes to euclid-gcd.
796+
//
797+
// if x.digits.len + y.digits.len <= 4 {
798+
// return x.gcd( y )
799+
// }
800+
801+
if x.signum == 0 {
802+
return y.abs()
803+
}
804+
if y.signum == 0 {
805+
return x.abs()
806+
}
807+
808+
if x.signum < 0 {
809+
return x.neg().gcd(y)
810+
}
811+
if y.signum < 0 {
812+
return x.gcd(y.neg())
813+
}
814+
815+
mut a := x
816+
mut b := y
817+
818+
mut az := a.msb()
819+
bz := b.msb()
820+
shift := util.umin(az, bz)
821+
b = b.rshift(bz)
822+
823+
for a.signum != 0 {
824+
a = a.rshift(az)
825+
diff := b - a
826+
az = diff.msb()
827+
b = bi_min(a, b)
828+
a = diff.abs()
829+
}
830+
return b.lshift(shift)
831+
}

0 commit comments

Comments
 (0)