admin/mruby-mruby
1
0
Fork 0
mirror of https://github.com/mruby/mruby synced 2026-06-08 16:11:16 +00:00
Code Issues Packages Projects Releases Wiki Activity

mruby-bigint: rewrite mpz_gcd main loop as binary Stein algorithm

Browse Source 
Replace the classical Euclidean main loop (mpz_mod per iteration)
with Stein's binary GCD: subtract + factor out trailing 2s on
odd-maintained operands. Keep an mpz_mod fallback for heavily
unbalanced pairs (the smaller operand has at least two fewer limbs)
where one long division replaces many Stein subtracts.

The old loop allocated a temporary mpz_t every iteration to hold the
mod result; the Stein loop is allocation-free thanks to the in-place
paths in mpz_sub and mpz_div_2exp. This matches the "memory first"
priority and also happens to be faster on typical inputs because
mpz_mod's setup cost dominates when the quotient is small (the
classical Fibonacci-neighbor worst case).

Also add a small mpz_swap helper used by the new loop.

Benchmark (bin/mruby benchmark/bm_bigint_gcd.rb, median of 3):

  case                       before    after    ratio
  single-limb                 42 ms    50 ms    1.19  (fast-path
                                                      unchanged;
                                                      noise)
  fib(200) vs fib(201)        82 ms    23 ms    0.28
  balanced ~700-bit shared    46 ms    32 ms    0.70
  unbalanced big vs small     35 ms    28 ms    0.80
  power-of-2 path             38 ms    42 ms    1.11  (fast-path
                                                      unchanged;
                                                      noise)
  balanced ~2800-bit shared   19 ms    15 ms    0.79

Co-authored-by: Claude <noreply@anthropic.com>
This commit is contained in:
Yukihiro "Matz" Matsumoto
2026-04-23 22:52:05 +09:00
parent ec55e0d0f1
commit 8e91554c6d
1 changed files with 35 additions and 10 deletions
Download Patch File Download Diff File Expand all files Collapse all files
mrbgems/mruby-bigint/core/bigint.c
+35 -10
View File
@@ -348,6 +348,14 @@ mpz_move(mpz_ctx_t *ctx, mpz_t *y, mpz_t *x)
x->sz = 0;
}
static inline void
mpz_swap(mpz_t *a, mpz_t *b)
{
mpz_t tmp = *a;
*a = *b;
*b = tmp;
}
static size_t
digits(mpz_t *x)
{
@@ -4821,9 +4829,11 @@ mpz_power_of_2_p(mpz_t *x)
return (limb != 0) && ((limb & (limb - 1)) == 0);
}
/* Hybrid GCD: a binary (Stein) prelude factors out common powers of 2
and handles single-limb / power-of-2 fast paths; the multi-limb main
loop is classical Euclidean using mpz_mod. */
/* Binary GCD (Stein's algorithm): factor out common powers of 2,
then iterate on odd operands with subtract + trailing-zero shift.
For heavily unbalanced pairs (one operand has at least two more
limbs than the other) a single Euclidean step via mpz_mod replaces
many Stein subtracts. */
static void
mpz_gcd(mpz_ctx_t *ctx, mpz_t *gg, mpz_t *aa, mpz_t *bb)
{
@@ -4895,14 +4905,29 @@ mpz_gcd(mpz_ctx_t *ctx, mpz_t *gg, mpz_t *aa, mpz_t *bb)
mpz_div_2exp(ctx, &a, &a, a_zeros);
mpz_div_2exp(ctx, &b, &b, b_zeros);
/* Euclidean algorithm for multi-limb numbers */
/* Stein main loop. Invariant: a and b are positive and odd.
Euclidean fallback when b has >=2 more limbs than a. */
while (!zero_p(&b)) {
mpz_t temp;
mpz_init_temp(ctx, &temp, a.sz);
mpz_mod(ctx, &temp, &a, &b);
mpz_move(ctx, &a, &b);
mpz_move(ctx, &b, &temp);
mpz_clear(ctx, &temp);
if (mpz_cmp(ctx, &a, &b) > 0) {
mpz_swap(&a, &b);
}
if (b.sz >= a.sz + 2) {
mpz_t temp;
mpz_init_temp(ctx, &temp, a.sz);
mpz_mod(ctx, &temp, &b, &a);
mpz_move(ctx, &b, &temp);
mpz_clear(ctx, &temp);
if (zero_p(&b)) break;
size_t bz = mpz_trailing_zeros(&b);
if (bz > 0)
mpz_div_2exp(ctx, &b, &b, bz);
}
else {
mpz_sub(ctx, &b, &b, &a);
if (zero_p(&b)) break;
size_t bz = mpz_trailing_zeros(&b);
mpz_div_2exp(ctx, &b, &b, bz);
}
}
mpz_mul_2exp(ctx, gg, &a, shift);
mpz_clear(ctx, &a);