diff --git a/mrbgems/mruby-bigint/core/bigint.c b/mrbgems/mruby-bigint/core/bigint.c index a36469f0c..49412c735 100644 --- a/mrbgems/mruby-bigint/core/bigint.c +++ b/mrbgems/mruby-bigint/core/bigint.c @@ -37,7 +37,6 @@ typedef struct mpz_pool { /* Forward declarations */ static int mpz_mul_sliding_window(mrb_state *mrb, mpz_t *result, mpz_t *first, mpz_t *second); static int udiv_pool(mrb_state *mrb, mpz_t *qq, mpz_t *rr, mpz_t *xx, mpz_t *yy); -static int mpz_sqrt_pool(mrb_state *mrb, mpz_t *z, mpz_t *x); /* Memory allocation tracking for benchmarking */ typedef struct allocation_stats { @@ -2708,12 +2707,6 @@ mpz_barrett_reduce(mrb_state *mrb, mpz_t *r, mpz_t *x, mpz_t *m, mpz_t *mu) static void mpz_sqrt(mrb_state *mrb, mpz_t *z, mpz_t *x) { - /* Try pool-based square root first for eligible operands */ - if (mpz_sqrt_pool(mrb, z, x)) { - return; /* Success with pool-based square root */ - } - - /* Fallback to traditional algorithm */ mrb_assert(x->sn >= 0); if (x->sz == 0) { @@ -2751,203 +2744,6 @@ mpz_sqrt(mrb_state *mrb, mpz_t *z, mpz_t *x) mpz_clear(mrb, &t); } -/* Pool-based square root using Newton-Raphson method with stack memory */ -static int -mpz_sqrt_pool(mrb_state *mrb, mpz_t *z, mpz_t *x) -{ - mrb_assert(x->sn >= 0); - - if (x->sz == 0) { - // sqrt(0) = 0 - z->sn = 0; - z->sz = 0; - return 1; /* Success - trivial case */ - } - - /* Only use pools for medium-sized operands that will benefit from stack allocation */ - size_t x_limbs = x->sz; - if (x_limbs < 4 || x_limbs > 64) { - return 0; /* Use traditional square root */ - } - - /* Estimate space needed: s (~x_limbs/2 + 1), t (~x_limbs/2 + 1), quotient (~x_limbs), remainder (~x_limbs) */ - size_t estimated_s_size = (x_limbs + 1) / 2 + 2; /* sqrt result size + margin */ - size_t temp_space = estimated_s_size * 3 + x_limbs * 2; /* s, t, division temps */ - if (temp_space > BIGINT_POOL_DEFAULT_SIZE / 2) { - return 0; /* Pool too small for all temporaries */ - } - - do { - mpz_pool_t pool_storage = {0}; - pool_storage.capacity = BIGINT_POOL_DEFAULT_SIZE; - pool_storage.active = 1; - mpz_pool_t *pool = &pool_storage; - - mpz_t s, t, quotient, remainder; - int pool_success = 0; - - /* Initialize temporary variables in pool */ - mpz_init_pool(mrb, &s, pool, estimated_s_size); - mpz_init_pool(mrb, &t, pool, estimated_s_size); - mpz_init_pool(mrb, "ient, pool, x_limbs); - mpz_init_pool(mrb, &remainder, pool, x_limbs); - - /* Verify all pool allocations succeeded */ - if (MPZ_POOL_VERIFY_4(s, t, quotient, remainder, pool)) { - - /* Estimate initial value: 1 << (bit_length(x) / 2) */ - size_t xbits = mpz_bits(x); - size_t sbit = (xbits + 1) / 2; - - /* Initialize s = 1 << sbit using pool memory */ - if (sbit == 0) { - s.p[0] = 1; - s.sz = 1; - s.sn = 1; - } - else { - size_t limb_shift = sbit / DIG_SIZE; - size_t bit_shift = sbit % DIG_SIZE; - - /* Clear s array */ - for (size_t i = 0; i < s.sz; i++) { - s.p[i] = 0; - } - - if (limb_shift < s.sz) { - s.p[limb_shift] = ((mp_limb)1) << bit_shift; - s.sz = limb_shift + 1; - s.sn = 1; - } - else { - /* Fallback if shift too large */ - s.p[0] = 1; - s.sz = 1; - s.sn = 1; - } - } - - /* Newton-Raphson iteration: s = (s + x / s) / 2 */ - int max_iterations = 100; /* Safety limit */ - for (int iter = 0; iter < max_iterations; iter++) { - /* t = x / s using pool-based division */ - if (udiv_pool(mrb, "ient, &remainder, x, &s)) { - /* Pool division succeeded - copy quotient to t */ - for (size_t i = 0; i < quotient.sz && i < t.sz; i++) { - t.p[i] = quotient.p[i]; - } - t.sz = (quotient.sz < t.sz) ? quotient.sz : t.sz; - t.sn = quotient.sn; - } - else { - /* Pool division failed - fallback to traditional */ - pool_success = 0; - break; - } - - /* t = t + s using pool-based addition */ - mpz_t temp_sum; - mpz_init_pool(mrb, &temp_sum, pool, estimated_s_size + 1); - if (is_pool_memory(&temp_sum, pool)) { - mpz_add(mrb, &temp_sum, &t, &s); - /* Copy sum back to t */ - for (size_t i = 0; i < temp_sum.sz && i < t.sz; i++) { - t.p[i] = temp_sum.p[i]; - } - t.sz = (temp_sum.sz < t.sz) ? temp_sum.sz : t.sz; - t.sn = temp_sum.sn; - MPZ_POOL_CLEANUP(mrb, temp_sum, pool); - } - else { - /* Pool addition failed - fallback to traditional */ - MPZ_POOL_CLEANUP(mrb, temp_sum, pool); - pool_success = 0; - break; - } - - /* t = t / 2 (right shift by 1 bit) */ - mp_limb carry = 0; - for (size_t i = t.sz; i > 0; i--) { - size_t idx = i - 1; - mp_limb current = t.p[idx]; - t.p[idx] = (current >> 1) | carry; - carry = (current & 1) << (DIG_SIZE - 1); - } - - /* Trim leading zeros */ - while (t.sz > 1 && t.p[t.sz - 1] == 0) { - t.sz--; - } - if (t.sz == 0) { - t.sz = 1; - t.p[0] = 0; - t.sn = 0; - } - - /* Check convergence: if t >= s, we're done */ - int cmp_result = 0; - if (t.sz > s.sz) { - cmp_result = 1; - } - else if (t.sz < s.sz) { - cmp_result = -1; - } - else { - for (size_t i = t.sz; i > 0; i--) { - size_t idx = i - 1; - if (t.p[idx] > s.p[idx]) { - cmp_result = 1; - break; - } - else if (t.p[idx] < s.p[idx]) { - cmp_result = -1; - break; - } - } - } - - if (cmp_result >= 0) { - /* Converged: t >= s */ - pool_success = 1; - break; - } - - /* s = t for next iteration */ - for (size_t i = 0; i < t.sz && i < s.sz; i++) { - s.p[i] = t.p[i]; - } - s.sz = (t.sz < s.sz) ? t.sz : s.sz; - s.sn = t.sn; - } - - if (pool_success) { - /* Copy final result from pool to heap-allocated output */ - mpz_realloc(mrb, z, s.sz); - for (size_t i = 0; i < s.sz; i++) { - z->p[i] = s.p[i]; - } - z->sz = s.sz; - z->sn = s.sn; - } - } - - /* Pool cleanup is automatic */ - MPZ_POOL_CLEANUP(mrb, s, pool); - MPZ_POOL_CLEANUP(mrb, t, pool); - MPZ_POOL_CLEANUP(mrb, quotient, pool); - MPZ_POOL_CLEANUP(mrb, remainder, pool); - pool_storage.active = 0; - - if (pool_success) { - return 1; /* Success - used pool memory for square root! */ - } - else { - return 0; /* Pool allocation failed, fallback */ - } - } while(0); - - return 0; /* Should not reach here */ -} /* Barrett reduction for efficient modular arithmetic with repeated operations */