mirror/guile
mirror/guile
1
Fork 0
mirror of https://git.savannah.gnu.org/git/guile.git synced 2025-05-20 03:30:27 +02:00
Code Activity

Implement gcd with new integer lib

Browse source 
* libguile/integers.c (scm_integer_gcd_ii)
(scm_integer_gcd_zi, scm_integer_gcd_zz): New internal functions.
* libguile/integers.h: Declare internal functions.
* libguile/numbers.c (scm_gcd): Use the new functions.
This commit is contained in:
Andy Wingo 2021-12-19 08:57:06 +01:00
parent 025c7c8045
commit 1e0797db72
3 changed files with 86 additions and 81 deletions
Download patch file Download diff file Expand all files Collapse all files

75
libguile/integers.c
View file

@ -1742,3 +1742,78 @@ scm_integer_round_divide_zz (SCM x, SCM y, SCM *qp, SCM *rp)
{
integer_round_divide_zz (scm_bignum (x), scm_bignum (y), qp, rp);
}
SCM
scm_integer_gcd_ii (scm_t_inum x, scm_t_inum y)
{
scm_t_inum u = x < 0 ? -x : x;
scm_t_inum v = y < 0 ? -y : y;
scm_t_inum result;
if (x == 0)
result = v;
else if (y == 0)
result = u;
else
{
int k = 0;
/* Determine a common factor 2^k */
while (((u | v) & 1) == 0)
{
k++;
u >>= 1;
v >>= 1;
}
/* Now, any factor 2^n can be eliminated */
if ((u & 1) == 0)
while ((u & 1) == 0)
u >>= 1;
else
while ((v & 1) == 0)
v >>= 1;
/* Both u and v are now odd. Subtract the smaller one
from the larger one to produce an even number, remove
more factors of two, and repeat. */
while (u != v)
{
if (u > v)
{
u -= v;
while ((u & 1) == 0)
u >>= 1;
}
else
{
v -= u;
while ((v & 1) == 0)
v >>= 1;
}
}
result = u << k;
}
return ulong_to_scm (result);
}
SCM
scm_integer_gcd_zi (SCM x, scm_t_inum y)
{
scm_t_bits result;
if (y == 0)
return scm_abs (x);
if (y < 0)
y = -y;
result = mpz_gcd_ui (NULL, SCM_I_BIG_MPZ (x), y);
scm_remember_upto_here_1 (x);
return ulong_to_scm (result);
}
SCM
scm_integer_gcd_zz (SCM x, SCM y)
{
mpz_t result, zx, zy;
mpz_init (result);
alias_bignum_to_mpz (scm_bignum (x), zx);
alias_bignum_to_mpz (scm_bignum (y), zy);
mpz_gcd (result, zx, zy);
scm_remember_upto_here_2 (x, y);
return take_mpz (result);
}

4
libguile/integers.h
View file

@ -124,6 +124,10 @@ SCM_INTERNAL void scm_integer_round_divide_zi (SCM x, scm_t_inum y,
SCM_INTERNAL void scm_integer_round_divide_zz (SCM x, SCM y,
SCM *qp, SCM *rp);
SCM_INTERNAL SCM scm_integer_gcd_ii (scm_t_inum x, scm_t_inum y);
SCM_INTERNAL SCM scm_integer_gcd_zi (SCM x, scm_t_inum y);
SCM_INTERNAL SCM scm_integer_gcd_zz (SCM x, SCM y);
#endif /* SCM_INTEGERS_H */

88
libguile/numbers.c
View file

@ -2782,68 +2782,15 @@ SCM_PRIMITIVE_GENERIC (scm_i_gcd, "gcd", 0, 2, 1,
SCM
scm_gcd (SCM x, SCM y)
{
if (SCM_UNLIKELY (SCM_UNBNDP (y)))
if (SCM_UNBNDP (y))
return SCM_UNBNDP (x) ? SCM_INUM0 : scm_abs (x);
if (SCM_LIKELY (SCM_I_INUMP (x)))
if (SCM_I_INUMP (x))
{
if (SCM_LIKELY (SCM_I_INUMP (y)))
{
scm_t_inum xx = SCM_I_INUM (x);
scm_t_inum yy = SCM_I_INUM (y);
scm_t_inum u = xx < 0 ? -xx : xx;
scm_t_inum v = yy < 0 ? -yy : yy;
scm_t_inum result;
if (SCM_UNLIKELY (xx == 0))
result = v;
else if (SCM_UNLIKELY (yy == 0))
result = u;
else
{
int k = 0;
/* Determine a common factor 2^k */
while (((u | v) & 1) == 0)
{
k++;
u >>= 1;
v >>= 1;
}
/* Now, any factor 2^n can be eliminated */
if ((u & 1) == 0)
while ((u & 1) == 0)
u >>= 1;
else
while ((v & 1) == 0)
v >>= 1;
/* Both u and v are now odd. Subtract the smaller one
from the larger one to produce an even number, remove
more factors of two, and repeat. */
while (u != v)
{
if (u > v)
{
u -= v;
while ((u & 1) == 0)
u >>= 1;
}
else
{
v -= u;
while ((v & 1) == 0)
v >>= 1;
}
}
result = u << k;
}
return (SCM_POSFIXABLE (result)
? SCM_I_MAKINUM (result)
: scm_i_inum2big (result));
}
if (SCM_I_INUMP (y))
return scm_integer_gcd_ii (SCM_I_INUM (x), SCM_I_INUM (y));
else if (SCM_BIGP (y))
{
SCM_SWAP (x, y);
goto big_inum;
}
return scm_integer_gcd_zi (y, SCM_I_INUM (x));
else if (SCM_REALP (y) && scm_is_integer (y))
goto handle_inexacts;
else
@ -2852,30 +2799,9 @@ scm_gcd (SCM x, SCM y)
else if (SCM_BIGP (x))
{
if (SCM_I_INUMP (y))
{
scm_t_bits result;
scm_t_inum yy;
big_inum:
yy = SCM_I_INUM (y);
if (yy == 0)
return scm_abs (x);
if (yy < 0)
yy = -yy;
result = mpz_gcd_ui (NULL, SCM_I_BIG_MPZ (x), yy);
scm_remember_upto_here_1 (x);
return (SCM_POSFIXABLE (result)
? SCM_I_MAKINUM (result)
: scm_from_unsigned_integer (result));
}
return scm_integer_gcd_zi (x, SCM_I_INUM (y));
else if (SCM_BIGP (y))
{
SCM result = scm_i_mkbig ();
mpz_gcd (SCM_I_BIG_MPZ (result),
SCM_I_BIG_MPZ (x),
SCM_I_BIG_MPZ (y));
scm_remember_upto_here_2 (x, y);
return scm_i_normbig (result);
}
return scm_integer_gcd_zz (x, y);
else if (SCM_REALP (y) && scm_is_integer (y))
goto handle_inexacts;
else