1 /* Copyright (C) 1995,1996,1997,1998,1999,2000,2001,2002,2003 Free Software Foundation, Inc.
3 * Portions Copyright 1990, 1991, 1992, 1993 by AT&T Bell Laboratories
4 * and Bellcore. See scm_divide.
7 * This library is free software; you can redistribute it and/or
8 * modify it under the terms of the GNU Lesser General Public
9 * License as published by the Free Software Foundation; either
10 * version 2.1 of the License, or (at your option) any later version.
12 * This library is distributed in the hope that it will be useful,
13 * but WITHOUT ANY WARRANTY; without even the implied warranty of
14 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
15 * Lesser General Public License for more details.
17 * You should have received a copy of the GNU Lesser General Public
18 * License along with this library; if not, write to the Free Software
19 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
23 /* General assumptions:
24 * All objects satisfying SCM_COMPLEXP() have a non-zero complex component.
25 * All objects satisfying SCM_BIGP() are too large to fit in a fixnum.
26 * If an object satisfies integer?, it's either an inum, a bignum, or a real.
27 * If floor (r) == r, r is an int, and mpz_set_d will DTRT.
28 * All objects satisfying SCM_FRACTIONP are never an integer.
33 - see if special casing bignums and reals in integer-exponent when
34 possible (to use mpz_pow and mpf_pow_ui) is faster.
36 - look in to better short-circuiting of common cases in
37 integer-expt and elsewhere.
39 - see if direct mpz operations can help in ash and elsewhere.
43 /* tell glibc (2.3) to give prototype for C99 trunc() */
55 #include "libguile/_scm.h"
56 #include "libguile/feature.h"
57 #include "libguile/ports.h"
58 #include "libguile/root.h"
59 #include "libguile/smob.h"
60 #include "libguile/strings.h"
62 #include "libguile/validate.h"
63 #include "libguile/numbers.h"
64 #include "libguile/deprecation.h"
66 #include "libguile/eq.h"
71 Wonder if this might be faster for some of our code? A switch on
72 the numtag would jump directly to the right case, and the
73 SCM_I_NUMTAG code might be faster than repeated SCM_FOOP tests...
75 #define SCM_I_NUMTAG_NOTNUM 0
76 #define SCM_I_NUMTAG_INUM 1
77 #define SCM_I_NUMTAG_BIG scm_tc16_big
78 #define SCM_I_NUMTAG_REAL scm_tc16_real
79 #define SCM_I_NUMTAG_COMPLEX scm_tc16_complex
80 #define SCM_I_NUMTAG(x) \
81 (SCM_INUMP(x) ? SCM_I_NUMTAG_INUM \
82 : (SCM_IMP(x) ? SCM_I_NUMTAG_NOTNUM \
83 : (((0xfcff & SCM_CELL_TYPE (x)) == scm_tc7_number) ? SCM_TYP16(x) \
84 : SCM_I_NUMTAG_NOTNUM)))
86 /* the macro above will not work as is with fractions */
89 #define SCM_SWAP(x, y) do { SCM __t = x; x = y; y = __t; } while (0)
91 /* FLOBUFLEN is the maximum number of characters neccessary for the
92 * printed or scm_string representation of an inexact number.
94 #define FLOBUFLEN (10+2*(sizeof(double)/sizeof(char)*SCM_CHAR_BIT*3+9)/10)
97 #if ! defined (HAVE_ISNAN)
102 return (IsNANorINF (x
) && NaN (x
) && ! IsINF (x
)) ? 1 : 0;
105 #if ! defined (HAVE_ISINF)
110 return (IsNANorINF (x
) && IsINF (x
)) ? 1 : 0;
117 /* mpz_cmp_d only recognises infinities in gmp 4.2 and up.
118 For prior versions use an explicit check here. */
119 #if __GNU_MP_VERSION < 4 \
120 || (__GNU_MP_VERSION == 4 && __GNU_MP_VERSION_MINOR < 2)
121 #define xmpz_cmp_d(z, d) \
122 (xisinf (d) ? (d < 0.0 ? 1 : -1) : mpz_cmp_d (z, d))
124 #define xmpz_cmp_d(z, d) mpz_cmp_d (z, d)
130 #if defined (HAVE_ISINF)
132 #elif defined (HAVE_FINITE) && defined (HAVE_ISNAN)
133 return (! (finite (x
) || isnan (x
)));
142 #if defined (HAVE_ISNAN)
151 static SCM abs_most_negative_fixnum
;
152 static mpz_t z_negative_one
;
156 static const char s_bignum
[] = "bignum";
158 SCM_C_INLINE_KEYWORD SCM
161 /* Return a newly created bignum. */
162 SCM z
= scm_double_cell (scm_tc16_big
, 0, 0, 0);
163 mpz_init (SCM_I_BIG_MPZ (z
));
167 SCM_C_INLINE_KEYWORD
static SCM
168 scm_i_clonebig (SCM src_big
, int same_sign_p
)
170 /* Copy src_big's value, negate it if same_sign_p is false, and return. */
171 SCM z
= scm_double_cell (scm_tc16_big
, 0, 0, 0);
172 mpz_init_set (SCM_I_BIG_MPZ (z
), SCM_I_BIG_MPZ (src_big
));
174 mpz_neg (SCM_I_BIG_MPZ (z
), SCM_I_BIG_MPZ (z
));
178 SCM_C_INLINE_KEYWORD
int
179 scm_i_bigcmp (SCM x
, SCM y
)
181 /* Return neg if x < y, pos if x > y, and 0 if x == y */
182 /* presume we already know x and y are bignums */
183 int result
= mpz_cmp (SCM_I_BIG_MPZ (x
), SCM_I_BIG_MPZ (y
));
184 scm_remember_upto_here_2 (x
, y
);
188 SCM_C_INLINE_KEYWORD SCM
189 scm_i_dbl2big (double d
)
191 /* results are only defined if d is an integer */
192 SCM z
= scm_double_cell (scm_tc16_big
, 0, 0, 0);
193 mpz_init_set_d (SCM_I_BIG_MPZ (z
), d
);
197 /* Convert a integer in double representation to a SCM number. */
199 SCM_C_INLINE_KEYWORD SCM
200 scm_i_dbl2num (double u
)
202 /* SCM_MOST_POSITIVE_FIXNUM+1 and SCM_MOST_NEGATIVE_FIXNUM are both
203 powers of 2, so there's no rounding when making "double" values
204 from them. If plain SCM_MOST_POSITIVE_FIXNUM was used it could
205 get rounded on a 64-bit machine, hence the "+1".
207 The use of floor() to force to an integer value ensures we get a
208 "numerically closest" value without depending on how a
209 double->long cast or how mpz_set_d will round. For reference,
210 double->long probably follows the hardware rounding mode,
211 mpz_set_d truncates towards zero. */
213 /* XXX - what happens when SCM_MOST_POSITIVE_FIXNUM etc is not
214 representable as a double? */
216 if (u
< (double) (SCM_MOST_POSITIVE_FIXNUM
+1)
217 && u
>= (double) SCM_MOST_NEGATIVE_FIXNUM
)
218 return SCM_MAKINUM ((long) u
);
220 return scm_i_dbl2big (u
);
223 /* scm_i_big2dbl() rounds to the closest representable double, in accordance
224 with R5RS exact->inexact.
226 The approach is to use mpz_get_d to pick out the high DBL_MANT_DIG bits
227 (ie. it truncates towards zero), then adjust to get the closest double by
228 examining the next lower bit and adding 1 if necessary.
230 Note that bignums exactly half way between representable doubles are
231 rounded to the next higher absolute value (ie. away from zero). This
232 seems like an adequate interpretation of R5RS "numerically closest", and
233 it's easier and faster than a full "nearest-even" style.
235 The bit test is done on the absolute value of the mpz_t, which means we
236 must use mpz_getlimbn. mpz_tstbit is not right, it treats negatives as
239 Prior to GMP 4.2, the rounding done by mpz_get_d was unspecified. It
240 happened to follow the hardware rounding mode, but on the absolute value
241 of its operand. This is not what we want, so we put the high
242 DBL_MANT_DIG bits into a temporary. This extra init/clear is a slowdown,
243 but doesn't matter too much since it's only for older GMP. */
246 scm_i_big2dbl (SCM b
)
251 bits
= mpz_sizeinbase (SCM_I_BIG_MPZ (b
), 2);
253 #if __GNU_MP_VERSION < 4 \
254 || (__GNU_MP_VERSION == 4 && __GNU_MP_VERSION_MINOR < 2)
256 /* GMP prior to 4.2, force truncate towards zero */
258 if (bits
> DBL_MANT_DIG
)
260 size_t shift
= bits
- DBL_MANT_DIG
;
261 mpz_init2 (tmp
, DBL_MANT_DIG
);
262 mpz_tdiv_q_2exp (tmp
, SCM_I_BIG_MPZ (b
), shift
);
263 result
= ldexp (mpz_get_d (tmp
), shift
);
268 result
= mpz_get_d (SCM_I_BIG_MPZ (b
));
273 result
= mpz_get_d (SCM_I_BIG_MPZ (b
));
276 if (bits
> DBL_MANT_DIG
)
278 unsigned long pos
= bits
- DBL_MANT_DIG
- 1;
279 /* test bit number "pos" in absolute value */
280 if (mpz_getlimbn (SCM_I_BIG_MPZ (b
), pos
/ GMP_NUMB_BITS
)
281 & ((mp_limb_t
) 1 << (pos
% GMP_NUMB_BITS
)))
283 result
+= ldexp ((double) mpz_sgn (SCM_I_BIG_MPZ (b
)), pos
+ 1);
287 scm_remember_upto_here_1 (b
);
291 SCM_C_INLINE_KEYWORD SCM
292 scm_i_normbig (SCM b
)
294 /* convert a big back to a fixnum if it'll fit */
295 /* presume b is a bignum */
296 if (mpz_fits_slong_p (SCM_I_BIG_MPZ (b
)))
298 long val
= mpz_get_si (SCM_I_BIG_MPZ (b
));
299 if (SCM_FIXABLE (val
))
300 b
= SCM_MAKINUM (val
);
305 static SCM_C_INLINE_KEYWORD SCM
306 scm_i_mpz2num (mpz_t b
)
308 /* convert a mpz number to a SCM number. */
309 if (mpz_fits_slong_p (b
))
311 long val
= mpz_get_si (b
);
312 if (SCM_FIXABLE (val
))
313 return SCM_MAKINUM (val
);
317 SCM z
= scm_double_cell (scm_tc16_big
, 0, 0, 0);
318 mpz_init_set (SCM_I_BIG_MPZ (z
), b
);
323 /* this is needed when we want scm_divide to make a float, not a ratio, even if passed two ints */
324 static SCM
scm_divide2real (SCM x
, SCM y
);
327 scm_make_ratio (SCM numerator
, SCM denominator
)
330 return scm_divide2real(numerator
, denominator
);
332 #define FUNC_NAME "make-ratio"
333 if (SCM_INUMP (denominator
))
335 if (SCM_EQ_P (denominator
, SCM_INUM0
))
336 scm_num_overflow ("make-ratio");
337 if (SCM_EQ_P (denominator
, SCM_MAKINUM(1)))
342 if (!(SCM_BIGP(denominator
)))
343 SCM_WRONG_TYPE_ARG (2, denominator
);
345 if (SCM_INUMP (numerator
))
347 if (SCM_EQ_P (numerator
, SCM_INUM0
))
349 if (SCM_INUMP (denominator
))
352 x
= SCM_INUM (numerator
);
353 y
= SCM_INUM (denominator
);
355 return SCM_MAKINUM(1);
357 return SCM_MAKINUM (x
/ y
);
359 return scm_double_cell (scm_tc16_fraction
, (scm_t_bits
)SCM_MAKINUM(-x
), (scm_t_bits
)SCM_MAKINUM(-y
), 0);
360 else return scm_double_cell (scm_tc16_fraction
, (scm_t_bits
)numerator
, (scm_t_bits
)denominator
, 0);
364 /* I assume bignums are actually big, so here there's no point in looking for a integer */
365 int sgn
= mpz_sgn (SCM_I_BIG_MPZ (denominator
));
366 if (sgn
< 0) /* if denominator negative, flip signs */
367 return scm_double_cell (scm_tc16_fraction
,
368 (scm_t_bits
)scm_difference (numerator
, SCM_UNDEFINED
),
369 (scm_t_bits
)scm_difference (denominator
, SCM_UNDEFINED
),
371 else return scm_double_cell (scm_tc16_fraction
, (scm_t_bits
)numerator
, (scm_t_bits
)denominator
, 0);
373 /* should this use SCM_UNPACK for the bignums? */
378 if (SCM_BIGP (numerator
))
380 /* can't use scm_divide to find integer here */
381 if (SCM_INUMP (denominator
))
383 long yy
= SCM_INUM (denominator
);
384 long abs_yy
= yy
< 0 ? -yy
: yy
;
385 int divisible_p
= mpz_divisible_ui_p (SCM_I_BIG_MPZ (numerator
), abs_yy
);
387 return scm_divide(numerator
, denominator
);
388 else return scm_double_cell (scm_tc16_fraction
, (scm_t_bits
)numerator
, (scm_t_bits
)denominator
, 0);
392 /* both are bignums */
393 if (SCM_EQ_P (numerator
, denominator
))
394 return SCM_MAKINUM(1);
395 int divisible_p
= mpz_divisible_p (SCM_I_BIG_MPZ (numerator
),
396 SCM_I_BIG_MPZ (denominator
));
398 return scm_divide(numerator
, denominator
);
399 else return scm_double_cell (scm_tc16_fraction
, (scm_t_bits
)numerator
, (scm_t_bits
)denominator
, 0);
402 else SCM_WRONG_TYPE_ARG (1, numerator
);
404 return SCM_BOOL_F
; /* won't happen */
409 static void scm_i_fraction_reduce (SCM z
)
411 if (!(SCM_FRACTION_REDUCED (z
)))
414 divisor
= scm_gcd (SCM_FRACTION_NUMERATOR (z
), SCM_FRACTION_DENOMINATOR (z
));
415 if (!(SCM_EQ_P (divisor
, SCM_MAKINUM(1))))
418 SCM_FRACTION_SET_NUMERATOR (z
, scm_divide (SCM_FRACTION_NUMERATOR (z
), divisor
));
419 SCM_FRACTION_SET_DENOMINATOR (z
, scm_divide (SCM_FRACTION_DENOMINATOR (z
), divisor
));
421 SCM_FRACTION_REDUCED_SET (z
);
426 scm_i_fraction2double (SCM z
)
428 return scm_num2dbl (scm_divide2real (SCM_FRACTION_NUMERATOR (z
),
429 SCM_FRACTION_DENOMINATOR (z
)),
433 SCM_DEFINE (scm_exact_p
, "exact?", 1, 0, 0,
435 "Return @code{#t} if @var{x} is an exact number, @code{#f}\n"
437 #define FUNC_NAME s_scm_exact_p
443 if (SCM_FRACTIONP (x
))
450 SCM_DEFINE (scm_odd_p
, "odd?", 1, 0, 0,
452 "Return @code{#t} if @var{n} is an odd number, @code{#f}\n"
454 #define FUNC_NAME s_scm_odd_p
458 long val
= SCM_INUM (n
);
459 return SCM_BOOL ((val
& 1L) != 0);
461 else if (SCM_BIGP (n
))
463 int odd_p
= mpz_odd_p (SCM_I_BIG_MPZ (n
));
464 scm_remember_upto_here_1 (n
);
465 return SCM_BOOL (odd_p
);
467 else if (!SCM_FALSEP (scm_inf_p (n
)))
469 else if (SCM_REALP (n
))
471 double rem
= fabs (fmod (SCM_REAL_VALUE(n
), 2.0));
477 SCM_WRONG_TYPE_ARG (1, n
);
480 SCM_WRONG_TYPE_ARG (1, n
);
485 SCM_DEFINE (scm_even_p
, "even?", 1, 0, 0,
487 "Return @code{#t} if @var{n} is an even number, @code{#f}\n"
489 #define FUNC_NAME s_scm_even_p
493 long val
= SCM_INUM (n
);
494 return SCM_BOOL ((val
& 1L) == 0);
496 else if (SCM_BIGP (n
))
498 int even_p
= mpz_even_p (SCM_I_BIG_MPZ (n
));
499 scm_remember_upto_here_1 (n
);
500 return SCM_BOOL (even_p
);
502 else if (!SCM_FALSEP (scm_inf_p (n
)))
504 else if (SCM_REALP (n
))
506 double rem
= fabs (fmod (SCM_REAL_VALUE(n
), 2.0));
512 SCM_WRONG_TYPE_ARG (1, n
);
515 SCM_WRONG_TYPE_ARG (1, n
);
519 SCM_DEFINE (scm_inf_p
, "inf?", 1, 0, 0,
521 "Return @code{#t} if @var{n} is infinite, @code{#f}\n"
523 #define FUNC_NAME s_scm_inf_p
526 return SCM_BOOL (xisinf (SCM_REAL_VALUE (n
)));
527 else if (SCM_COMPLEXP (n
))
528 return SCM_BOOL (xisinf (SCM_COMPLEX_REAL (n
))
529 || xisinf (SCM_COMPLEX_IMAG (n
)));
535 SCM_DEFINE (scm_nan_p
, "nan?", 1, 0, 0,
537 "Return @code{#t} if @var{n} is a NaN, @code{#f}\n"
539 #define FUNC_NAME s_scm_nan_p
542 return SCM_BOOL (xisnan (SCM_REAL_VALUE (n
)));
543 else if (SCM_COMPLEXP (n
))
544 return SCM_BOOL (xisnan (SCM_COMPLEX_REAL (n
))
545 || xisnan (SCM_COMPLEX_IMAG (n
)));
551 /* Guile's idea of infinity. */
552 static double guile_Inf
;
554 /* Guile's idea of not a number. */
555 static double guile_NaN
;
558 guile_ieee_init (void)
560 #if defined (HAVE_ISINF) || defined (HAVE_FINITE)
562 /* Some version of gcc on some old version of Linux used to crash when
563 trying to make Inf and NaN. */
567 guile_Inf
= 1.0 / (tmp
- tmp
);
568 #elif defined (__alpha__) && ! defined (linux)
569 extern unsigned int DINFINITY
[2];
570 guile_Inf
= (*(X_CAST(double *, DINFINITY
)));
577 if (guile_Inf
== tmp
)
585 #if defined (HAVE_ISNAN)
587 #if defined (__alpha__) && ! defined (linux)
588 extern unsigned int DQNAN
[2];
589 guile_NaN
= (*(X_CAST(double *, DQNAN
)));
591 guile_NaN
= guile_Inf
/ guile_Inf
;
597 SCM_DEFINE (scm_inf
, "inf", 0, 0, 0,
600 #define FUNC_NAME s_scm_inf
602 static int initialized
= 0;
608 return scm_make_real (guile_Inf
);
612 SCM_DEFINE (scm_nan
, "nan", 0, 0, 0,
615 #define FUNC_NAME s_scm_nan
617 static int initialized
= 0;
623 return scm_make_real (guile_NaN
);
628 SCM_PRIMITIVE_GENERIC (scm_abs
, "abs", 1, 0, 0,
630 "Return the absolute value of @var{x}.")
635 long int xx
= SCM_INUM (x
);
638 else if (SCM_POSFIXABLE (-xx
))
639 return SCM_MAKINUM (-xx
);
641 return scm_i_long2big (-xx
);
643 else if (SCM_BIGP (x
))
645 const int sgn
= mpz_sgn (SCM_I_BIG_MPZ (x
));
647 return scm_i_clonebig (x
, 0);
651 else if (SCM_REALP (x
))
652 return scm_make_real (fabs (SCM_REAL_VALUE (x
)));
653 else if (SCM_FRACTIONP (x
))
655 if (SCM_FALSEP (scm_negative_p (SCM_FRACTION_NUMERATOR (x
))))
657 return scm_make_ratio (scm_difference (SCM_FRACTION_NUMERATOR (x
), SCM_UNDEFINED
),
658 SCM_FRACTION_DENOMINATOR (x
));
661 SCM_WTA_DISPATCH_1 (g_scm_abs
, x
, 1, s_scm_abs
);
666 SCM_GPROC (s_quotient
, "quotient", 2, 0, 0, scm_quotient
, g_quotient
);
667 /* "Return the quotient of the numbers @var{x} and @var{y}."
670 scm_quotient (SCM x
, SCM y
)
674 long xx
= SCM_INUM (x
);
677 long yy
= SCM_INUM (y
);
679 scm_num_overflow (s_quotient
);
684 return SCM_MAKINUM (z
);
686 return scm_i_long2big (z
);
689 else if (SCM_BIGP (y
))
691 if ((SCM_INUM (x
) == SCM_MOST_NEGATIVE_FIXNUM
)
692 && (scm_i_bigcmp (abs_most_negative_fixnum
, y
) == 0))
693 /* Special case: x == fixnum-min && y == abs (fixnum-min) */
694 return SCM_MAKINUM (-1);
696 return SCM_MAKINUM (0);
699 SCM_WTA_DISPATCH_2 (g_quotient
, x
, y
, SCM_ARG2
, s_quotient
);
701 else if (SCM_BIGP (x
))
705 long yy
= SCM_INUM (y
);
707 scm_num_overflow (s_quotient
);
712 SCM result
= scm_i_mkbig ();
715 mpz_tdiv_q_ui (SCM_I_BIG_MPZ (result
),
718 mpz_neg (SCM_I_BIG_MPZ (result
), SCM_I_BIG_MPZ (result
));
721 mpz_tdiv_q_ui (SCM_I_BIG_MPZ (result
), SCM_I_BIG_MPZ (x
), yy
);
722 scm_remember_upto_here_1 (x
);
723 return scm_i_normbig (result
);
726 else if (SCM_BIGP (y
))
728 SCM result
= scm_i_mkbig ();
729 mpz_tdiv_q (SCM_I_BIG_MPZ (result
),
732 scm_remember_upto_here_2 (x
, y
);
733 return scm_i_normbig (result
);
736 SCM_WTA_DISPATCH_2 (g_quotient
, x
, y
, SCM_ARG2
, s_quotient
);
739 SCM_WTA_DISPATCH_2 (g_quotient
, x
, y
, SCM_ARG1
, s_quotient
);
742 SCM_GPROC (s_remainder
, "remainder", 2, 0, 0, scm_remainder
, g_remainder
);
743 /* "Return the remainder of the numbers @var{x} and @var{y}.\n"
745 * "(remainder 13 4) @result{} 1\n"
746 * "(remainder -13 4) @result{} -1\n"
750 scm_remainder (SCM x
, SCM y
)
756 long yy
= SCM_INUM (y
);
758 scm_num_overflow (s_remainder
);
761 long z
= SCM_INUM (x
) % yy
;
762 return SCM_MAKINUM (z
);
765 else if (SCM_BIGP (y
))
767 if ((SCM_INUM (x
) == SCM_MOST_NEGATIVE_FIXNUM
)
768 && (scm_i_bigcmp (abs_most_negative_fixnum
, y
) == 0))
769 /* Special case: x == fixnum-min && y == abs (fixnum-min) */
770 return SCM_MAKINUM (0);
775 SCM_WTA_DISPATCH_2 (g_remainder
, x
, y
, SCM_ARG2
, s_remainder
);
777 else if (SCM_BIGP (x
))
781 long yy
= SCM_INUM (y
);
783 scm_num_overflow (s_remainder
);
786 SCM result
= scm_i_mkbig ();
789 mpz_tdiv_r_ui (SCM_I_BIG_MPZ (result
), SCM_I_BIG_MPZ(x
), yy
);
790 scm_remember_upto_here_1 (x
);
791 return scm_i_normbig (result
);
794 else if (SCM_BIGP (y
))
796 SCM result
= scm_i_mkbig ();
797 mpz_tdiv_r (SCM_I_BIG_MPZ (result
),
800 scm_remember_upto_here_2 (x
, y
);
801 return scm_i_normbig (result
);
804 SCM_WTA_DISPATCH_2 (g_remainder
, x
, y
, SCM_ARG2
, s_remainder
);
807 SCM_WTA_DISPATCH_2 (g_remainder
, x
, y
, SCM_ARG1
, s_remainder
);
811 SCM_GPROC (s_modulo
, "modulo", 2, 0, 0, scm_modulo
, g_modulo
);
812 /* "Return the modulo of the numbers @var{x} and @var{y}.\n"
814 * "(modulo 13 4) @result{} 1\n"
815 * "(modulo -13 4) @result{} 3\n"
819 scm_modulo (SCM x
, SCM y
)
823 long xx
= SCM_INUM (x
);
826 long yy
= SCM_INUM (y
);
828 scm_num_overflow (s_modulo
);
831 /* FIXME: I think this may be a bug on some arches -- results
832 of % with negative second arg are undefined... */
850 return SCM_MAKINUM (result
);
853 else if (SCM_BIGP (y
))
855 int sgn_y
= mpz_sgn (SCM_I_BIG_MPZ (y
));
858 scm_num_overflow (s_modulo
);
866 SCM pos_y
= scm_i_clonebig (y
, 0);
867 /* do this after the last scm_op */
868 mpz_init_set_si (z_x
, xx
);
869 result
= pos_y
; /* re-use this bignum */
870 mpz_mod (SCM_I_BIG_MPZ (result
),
872 SCM_I_BIG_MPZ (pos_y
));
873 scm_remember_upto_here_1 (pos_y
);
877 result
= scm_i_mkbig ();
878 /* do this after the last scm_op */
879 mpz_init_set_si (z_x
, xx
);
880 mpz_mod (SCM_I_BIG_MPZ (result
),
883 scm_remember_upto_here_1 (y
);
886 if ((sgn_y
< 0) && mpz_sgn (SCM_I_BIG_MPZ (result
)) != 0)
887 mpz_add (SCM_I_BIG_MPZ (result
),
889 SCM_I_BIG_MPZ (result
));
890 scm_remember_upto_here_1 (y
);
891 /* and do this before the next one */
893 return scm_i_normbig (result
);
897 SCM_WTA_DISPATCH_2 (g_modulo
, x
, y
, SCM_ARG2
, s_modulo
);
899 else if (SCM_BIGP (x
))
903 long yy
= SCM_INUM (y
);
905 scm_num_overflow (s_modulo
);
908 SCM result
= scm_i_mkbig ();
909 mpz_mod_ui (SCM_I_BIG_MPZ (result
),
911 (yy
< 0) ? - yy
: yy
);
912 scm_remember_upto_here_1 (x
);
913 if ((yy
< 0) && (mpz_sgn (SCM_I_BIG_MPZ (result
)) != 0))
914 mpz_sub_ui (SCM_I_BIG_MPZ (result
),
915 SCM_I_BIG_MPZ (result
),
917 return scm_i_normbig (result
);
920 else if (SCM_BIGP (y
))
922 int sgn_y
= mpz_sgn (SCM_I_BIG_MPZ (y
));
924 scm_num_overflow (s_modulo
);
927 SCM result
= scm_i_mkbig ();
928 int y_sgn
= mpz_sgn (SCM_I_BIG_MPZ (y
));
929 SCM pos_y
= scm_i_clonebig (y
, y_sgn
>= 0);
930 mpz_mod (SCM_I_BIG_MPZ (result
),
932 SCM_I_BIG_MPZ (pos_y
));
934 scm_remember_upto_here_1 (x
);
935 if ((y_sgn
< 0) && (mpz_sgn (SCM_I_BIG_MPZ (result
)) != 0))
936 mpz_add (SCM_I_BIG_MPZ (result
),
938 SCM_I_BIG_MPZ (result
));
939 scm_remember_upto_here_2 (y
, pos_y
);
940 return scm_i_normbig (result
);
944 SCM_WTA_DISPATCH_2 (g_modulo
, x
, y
, SCM_ARG2
, s_modulo
);
947 SCM_WTA_DISPATCH_2 (g_modulo
, x
, y
, SCM_ARG1
, s_modulo
);
950 SCM_GPROC1 (s_gcd
, "gcd", scm_tc7_asubr
, scm_gcd
, g_gcd
);
951 /* "Return the greatest common divisor of all arguments.\n"
952 * "If called without arguments, 0 is returned."
955 scm_gcd (SCM x
, SCM y
)
958 return SCM_UNBNDP (x
) ? SCM_INUM0
: x
;
964 long xx
= SCM_INUM (x
);
965 long yy
= SCM_INUM (y
);
966 long u
= xx
< 0 ? -xx
: xx
;
967 long v
= yy
< 0 ? -yy
: yy
;
977 /* Determine a common factor 2^k */
978 while (!(1 & (u
| v
)))
984 /* Now, any factor 2^n can be eliminated */
1004 return (SCM_POSFIXABLE (result
)
1005 ? SCM_MAKINUM (result
)
1006 : scm_i_long2big (result
));
1008 else if (SCM_BIGP (y
))
1010 SCM result
= scm_i_mkbig ();
1011 SCM mx
= scm_i_mkbig ();
1012 mpz_set_si (SCM_I_BIG_MPZ (mx
), SCM_INUM (x
));
1013 scm_remember_upto_here_1 (x
);
1014 mpz_gcd (SCM_I_BIG_MPZ (result
),
1017 scm_remember_upto_here_2 (mx
, y
);
1018 return scm_i_normbig (result
);
1021 SCM_WTA_DISPATCH_2 (g_gcd
, x
, y
, SCM_ARG2
, s_gcd
);
1023 else if (SCM_BIGP (x
))
1027 unsigned long result
;
1028 long yy
= SCM_INUM (y
);
1033 result
= mpz_gcd_ui (NULL
, SCM_I_BIG_MPZ (x
), yy
);
1034 scm_remember_upto_here_1 (x
);
1035 return (SCM_POSFIXABLE (result
)
1036 ? SCM_MAKINUM (result
)
1037 : scm_ulong2num (result
));
1039 else if (SCM_BIGP (y
))
1041 SCM result
= scm_i_mkbig ();
1042 mpz_gcd (SCM_I_BIG_MPZ (result
),
1045 scm_remember_upto_here_2 (x
, y
);
1046 return scm_i_normbig (result
);
1049 SCM_WTA_DISPATCH_2 (g_gcd
, x
, y
, SCM_ARG2
, s_gcd
);
1052 SCM_WTA_DISPATCH_2 (g_gcd
, x
, y
, SCM_ARG1
, s_gcd
);
1055 SCM_GPROC1 (s_lcm
, "lcm", scm_tc7_asubr
, scm_lcm
, g_lcm
);
1056 /* "Return the least common multiple of the arguments.\n"
1057 * "If called without arguments, 1 is returned."
1060 scm_lcm (SCM n1
, SCM n2
)
1062 if (SCM_UNBNDP (n2
))
1064 if (SCM_UNBNDP (n1
))
1065 return SCM_MAKINUM (1L);
1066 n2
= SCM_MAKINUM (1L);
1069 SCM_GASSERT2 (SCM_INUMP (n1
) || SCM_BIGP (n1
),
1070 g_lcm
, n1
, n2
, SCM_ARG1
, s_lcm
);
1071 SCM_GASSERT2 (SCM_INUMP (n2
) || SCM_BIGP (n2
),
1072 g_lcm
, n1
, n2
, SCM_ARGn
, s_lcm
);
1078 SCM d
= scm_gcd (n1
, n2
);
1079 if (SCM_EQ_P (d
, SCM_INUM0
))
1082 return scm_abs (scm_product (n1
, scm_quotient (n2
, d
)));
1086 /* inum n1, big n2 */
1089 SCM result
= scm_i_mkbig ();
1090 long nn1
= SCM_INUM (n1
);
1091 if (nn1
== 0) return SCM_INUM0
;
1092 if (nn1
< 0) nn1
= - nn1
;
1093 mpz_lcm_ui (SCM_I_BIG_MPZ (result
), SCM_I_BIG_MPZ (n2
), nn1
);
1094 scm_remember_upto_here_1 (n2
);
1109 SCM result
= scm_i_mkbig ();
1110 mpz_lcm(SCM_I_BIG_MPZ (result
),
1112 SCM_I_BIG_MPZ (n2
));
1113 scm_remember_upto_here_2(n1
, n2
);
1114 /* shouldn't need to normalize b/c lcm of 2 bigs should be big */
1120 #ifndef scm_long2num
1121 #define SCM_LOGOP_RETURN(x) scm_ulong2num(x)
1123 #define SCM_LOGOP_RETURN(x) SCM_MAKINUM(x)
1126 /* Emulating 2's complement bignums with sign magnitude arithmetic:
1131 + + + x (map digit:logand X Y)
1132 + - + x (map digit:logand X (lognot (+ -1 Y)))
1133 - + + y (map digit:logand (lognot (+ -1 X)) Y)
1134 - - - (+ 1 (map digit:logior (+ -1 X) (+ -1 Y)))
1139 + + + (map digit:logior X Y)
1140 + - - y (+ 1 (map digit:logand (lognot X) (+ -1 Y)))
1141 - + - x (+ 1 (map digit:logand (+ -1 X) (lognot Y)))
1142 - - - x (+ 1 (map digit:logand (+ -1 X) (+ -1 Y)))
1147 + + + (map digit:logxor X Y)
1148 + - - (+ 1 (map digit:logxor X (+ -1 Y)))
1149 - + - (+ 1 (map digit:logxor (+ -1 X) Y))
1150 - - + (map digit:logxor (+ -1 X) (+ -1 Y))
1155 + + (any digit:logand X Y)
1156 + - (any digit:logand X (lognot (+ -1 Y)))
1157 - + (any digit:logand (lognot (+ -1 X)) Y)
1162 SCM_DEFINE1 (scm_logand
, "logand", scm_tc7_asubr
,
1164 "Return the bitwise AND of the integer arguments.\n\n"
1166 "(logand) @result{} -1\n"
1167 "(logand 7) @result{} 7\n"
1168 "(logand #b111 #b011 #\b001) @result{} 1\n"
1170 #define FUNC_NAME s_scm_logand
1174 if (SCM_UNBNDP (n2
))
1176 if (SCM_UNBNDP (n1
))
1177 return SCM_MAKINUM (-1);
1178 else if (!SCM_NUMBERP (n1
))
1179 SCM_WRONG_TYPE_ARG (SCM_ARG1
, n1
);
1180 else if (SCM_NUMBERP (n1
))
1183 SCM_WRONG_TYPE_ARG (SCM_ARG1
, n1
);
1188 nn1
= SCM_INUM (n1
);
1191 long nn2
= SCM_INUM (n2
);
1192 return SCM_MAKINUM (nn1
& nn2
);
1194 else if SCM_BIGP (n2
)
1200 SCM result_z
= scm_i_mkbig ();
1202 mpz_init_set_si (nn1_z
, nn1
);
1203 mpz_and (SCM_I_BIG_MPZ (result_z
), nn1_z
, SCM_I_BIG_MPZ (n2
));
1204 scm_remember_upto_here_1 (n2
);
1206 return scm_i_normbig (result_z
);
1210 SCM_WRONG_TYPE_ARG (SCM_ARG2
, n2
);
1212 else if (SCM_BIGP (n1
))
1217 nn1
= SCM_INUM (n1
);
1220 else if (SCM_BIGP (n2
))
1222 SCM result_z
= scm_i_mkbig ();
1223 mpz_and (SCM_I_BIG_MPZ (result_z
),
1225 SCM_I_BIG_MPZ (n2
));
1226 scm_remember_upto_here_2 (n1
, n2
);
1227 return scm_i_normbig (result_z
);
1230 SCM_WRONG_TYPE_ARG (SCM_ARG2
, n2
);
1233 SCM_WRONG_TYPE_ARG (SCM_ARG1
, n1
);
1238 SCM_DEFINE1 (scm_logior
, "logior", scm_tc7_asubr
,
1240 "Return the bitwise OR of the integer arguments.\n\n"
1242 "(logior) @result{} 0\n"
1243 "(logior 7) @result{} 7\n"
1244 "(logior #b000 #b001 #b011) @result{} 3\n"
1246 #define FUNC_NAME s_scm_logior
1250 if (SCM_UNBNDP (n2
))
1252 if (SCM_UNBNDP (n1
))
1254 else if (SCM_NUMBERP (n1
))
1257 SCM_WRONG_TYPE_ARG (SCM_ARG1
, n1
);
1262 nn1
= SCM_INUM (n1
);
1265 long nn2
= SCM_INUM (n2
);
1266 return SCM_MAKINUM (nn1
| nn2
);
1268 else if (SCM_BIGP (n2
))
1274 SCM result_z
= scm_i_mkbig ();
1276 mpz_init_set_si (nn1_z
, nn1
);
1277 mpz_ior (SCM_I_BIG_MPZ (result_z
), nn1_z
, SCM_I_BIG_MPZ (n2
));
1278 scm_remember_upto_here_1 (n2
);
1284 SCM_WRONG_TYPE_ARG (SCM_ARG2
, n2
);
1286 else if (SCM_BIGP (n1
))
1291 nn1
= SCM_INUM (n1
);
1294 else if (SCM_BIGP (n2
))
1296 SCM result_z
= scm_i_mkbig ();
1297 mpz_ior (SCM_I_BIG_MPZ (result_z
),
1299 SCM_I_BIG_MPZ (n2
));
1300 scm_remember_upto_here_2 (n1
, n2
);
1304 SCM_WRONG_TYPE_ARG (SCM_ARG2
, n2
);
1307 SCM_WRONG_TYPE_ARG (SCM_ARG1
, n1
);
1312 SCM_DEFINE1 (scm_logxor
, "logxor", scm_tc7_asubr
,
1314 "Return the bitwise XOR of the integer arguments. A bit is\n"
1315 "set in the result if it is set in an odd number of arguments.\n"
1317 "(logxor) @result{} 0\n"
1318 "(logxor 7) @result{} 7\n"
1319 "(logxor #b000 #b001 #b011) @result{} 2\n"
1320 "(logxor #b000 #b001 #b011 #b011) @result{} 1\n"
1322 #define FUNC_NAME s_scm_logxor
1326 if (SCM_UNBNDP (n2
))
1328 if (SCM_UNBNDP (n1
))
1330 else if (SCM_NUMBERP (n1
))
1333 SCM_WRONG_TYPE_ARG (SCM_ARG1
, n1
);
1338 nn1
= SCM_INUM (n1
);
1341 long nn2
= SCM_INUM (n2
);
1342 return SCM_MAKINUM (nn1
^ nn2
);
1344 else if (SCM_BIGP (n2
))
1348 SCM result_z
= scm_i_mkbig ();
1350 mpz_init_set_si (nn1_z
, nn1
);
1351 mpz_xor (SCM_I_BIG_MPZ (result_z
), nn1_z
, SCM_I_BIG_MPZ (n2
));
1352 scm_remember_upto_here_1 (n2
);
1354 return scm_i_normbig (result_z
);
1358 SCM_WRONG_TYPE_ARG (SCM_ARG2
, n2
);
1360 else if (SCM_BIGP (n1
))
1365 nn1
= SCM_INUM (n1
);
1368 else if (SCM_BIGP (n2
))
1370 SCM result_z
= scm_i_mkbig ();
1371 mpz_xor (SCM_I_BIG_MPZ (result_z
),
1373 SCM_I_BIG_MPZ (n2
));
1374 scm_remember_upto_here_2 (n1
, n2
);
1375 return scm_i_normbig (result_z
);
1378 SCM_WRONG_TYPE_ARG (SCM_ARG2
, n2
);
1381 SCM_WRONG_TYPE_ARG (SCM_ARG1
, n1
);
1386 SCM_DEFINE (scm_logtest
, "logtest", 2, 0, 0,
1389 "(logtest j k) @equiv{} (not (zero? (logand j k)))\n\n"
1390 "(logtest #b0100 #b1011) @result{} #f\n"
1391 "(logtest #b0100 #b0111) @result{} #t\n"
1393 #define FUNC_NAME s_scm_logtest
1402 long nk
= SCM_INUM (k
);
1403 return SCM_BOOL (nj
& nk
);
1405 else if (SCM_BIGP (k
))
1413 mpz_init_set_si (nj_z
, nj
);
1414 mpz_and (nj_z
, nj_z
, SCM_I_BIG_MPZ (k
));
1415 scm_remember_upto_here_1 (k
);
1416 result
= SCM_BOOL (mpz_sgn (nj_z
) != 0);
1422 SCM_WRONG_TYPE_ARG (SCM_ARG2
, k
);
1424 else if (SCM_BIGP (j
))
1432 else if (SCM_BIGP (k
))
1436 mpz_init (result_z
);
1440 scm_remember_upto_here_2 (j
, k
);
1441 result
= SCM_BOOL (mpz_sgn (result_z
) != 0);
1442 mpz_clear (result_z
);
1446 SCM_WRONG_TYPE_ARG (SCM_ARG2
, k
);
1449 SCM_WRONG_TYPE_ARG (SCM_ARG1
, j
);
1454 SCM_DEFINE (scm_logbit_p
, "logbit?", 2, 0, 0,
1457 "(logbit? index j) @equiv{} (logtest (integer-expt 2 index) j)\n\n"
1458 "(logbit? 0 #b1101) @result{} #t\n"
1459 "(logbit? 1 #b1101) @result{} #f\n"
1460 "(logbit? 2 #b1101) @result{} #t\n"
1461 "(logbit? 3 #b1101) @result{} #t\n"
1462 "(logbit? 4 #b1101) @result{} #f\n"
1464 #define FUNC_NAME s_scm_logbit_p
1466 unsigned long int iindex
;
1468 SCM_VALIDATE_INUM_MIN (SCM_ARG1
, index
, 0);
1469 iindex
= (unsigned long int) SCM_INUM (index
);
1472 return SCM_BOOL ((1L << iindex
) & SCM_INUM (j
));
1473 else if (SCM_BIGP (j
))
1475 int val
= mpz_tstbit (SCM_I_BIG_MPZ (j
), iindex
);
1476 scm_remember_upto_here_1 (j
);
1477 return SCM_BOOL (val
);
1480 SCM_WRONG_TYPE_ARG (SCM_ARG2
, j
);
1485 SCM_DEFINE (scm_lognot
, "lognot", 1, 0, 0,
1487 "Return the integer which is the ones-complement of the integer\n"
1491 "(number->string (lognot #b10000000) 2)\n"
1492 " @result{} \"-10000001\"\n"
1493 "(number->string (lognot #b0) 2)\n"
1494 " @result{} \"-1\"\n"
1496 #define FUNC_NAME s_scm_lognot
1498 if (SCM_INUMP (n
)) {
1499 /* No overflow here, just need to toggle all the bits making up the inum.
1500 Enhancement: No need to strip the tag and add it back, could just xor
1501 a block of 1 bits, if that worked with the various debug versions of
1503 return SCM_MAKINUM (~ SCM_INUM (n
));
1505 } else if (SCM_BIGP (n
)) {
1506 SCM result
= scm_i_mkbig ();
1507 mpz_com (SCM_I_BIG_MPZ (result
), SCM_I_BIG_MPZ (n
));
1508 scm_remember_upto_here_1 (n
);
1512 SCM_WRONG_TYPE_ARG (SCM_ARG1
, n
);
1517 SCM_DEFINE (scm_integer_expt
, "integer-expt", 2, 0, 0,
1519 "Return @var{n} raised to the non-negative integer exponent\n"
1523 "(integer-expt 2 5)\n"
1525 "(integer-expt -3 3)\n"
1528 #define FUNC_NAME s_scm_integer_expt
1531 SCM z_i2
= SCM_BOOL_F
;
1533 SCM acc
= SCM_MAKINUM (1L);
1535 /* 0^0 == 1 according to R5RS */
1536 if (SCM_EQ_P (n
, SCM_INUM0
) || SCM_EQ_P (n
, acc
))
1537 return SCM_FALSEP (scm_zero_p(k
)) ? n
: acc
;
1538 else if (SCM_EQ_P (n
, SCM_MAKINUM (-1L)))
1539 return SCM_FALSEP (scm_even_p (k
)) ? n
: acc
;
1543 else if (SCM_BIGP (k
))
1545 z_i2
= scm_i_clonebig (k
, 1);
1546 mpz_init_set (SCM_I_BIG_MPZ (z_i2
), SCM_I_BIG_MPZ (k
));
1547 scm_remember_upto_here_1 (k
);
1550 else if (SCM_REALP (k
))
1552 double r
= SCM_REAL_VALUE (k
);
1554 SCM_WRONG_TYPE_ARG (2, k
);
1555 if ((r
> SCM_MOST_POSITIVE_FIXNUM
) || (r
< SCM_MOST_NEGATIVE_FIXNUM
))
1557 z_i2
= scm_i_mkbig ();
1558 mpz_init_set_d (SCM_I_BIG_MPZ (z_i2
), r
);
1567 SCM_WRONG_TYPE_ARG (2, k
);
1571 if (mpz_sgn(SCM_I_BIG_MPZ (z_i2
)) == -1)
1573 mpz_neg (SCM_I_BIG_MPZ (z_i2
), SCM_I_BIG_MPZ (z_i2
));
1574 n
= scm_divide (n
, SCM_UNDEFINED
);
1578 if (mpz_sgn(SCM_I_BIG_MPZ (z_i2
)) == 0)
1580 mpz_clear (SCM_I_BIG_MPZ (z_i2
));
1583 if (mpz_cmp_ui(SCM_I_BIG_MPZ (z_i2
), 1) == 0)
1585 mpz_clear (SCM_I_BIG_MPZ (z_i2
));
1586 return scm_product (acc
, n
);
1588 if (mpz_tstbit(SCM_I_BIG_MPZ (z_i2
), 0))
1589 acc
= scm_product (acc
, n
);
1590 n
= scm_product (n
, n
);
1591 mpz_fdiv_q_2exp (SCM_I_BIG_MPZ (z_i2
), SCM_I_BIG_MPZ (z_i2
), 1);
1599 n
= scm_divide (n
, SCM_UNDEFINED
);
1606 return scm_product (acc
, n
);
1608 acc
= scm_product (acc
, n
);
1609 n
= scm_product (n
, n
);
1616 SCM_DEFINE (scm_ash
, "ash", 2, 0, 0,
1618 "Return @var{n} shifted left by @var{cnt} bits, or shifted right\n"
1619 "if @var{cnt} is negative. This is an ``arithmetic'' shift.\n"
1621 "This is effectively a multiplication by 2^@var{cnt}}, and when\n"
1622 "@var{cnt} is negative it's a division, rounded towards negative\n"
1623 "infinity. (Note that this is not the same rounding as\n"
1624 "@code{quotient} does.)\n"
1626 "With @var{n} viewed as an infinite precision twos complement,\n"
1627 "@code{ash} means a left shift introducing zero bits, or a right\n"
1628 "shift dropping bits.\n"
1631 "(number->string (ash #b1 3) 2) @result{} \"1000\"\n"
1632 "(number->string (ash #b1010 -1) 2) @result{} \"101\"\n"
1634 ";; -23 is bits ...11101001, -6 is bits ...111010\n"
1635 "(ash -23 -2) @result{} -6\n"
1637 #define FUNC_NAME s_scm_ash
1641 SCM_VALIDATE_INUM (2, cnt
);
1643 bits_to_shift
= SCM_INUM (cnt
);
1645 if (bits_to_shift
< 0)
1647 /* Shift right by abs(cnt) bits. This is realized as a division
1648 by div:=2^abs(cnt). However, to guarantee the floor
1649 rounding, negative values require some special treatment.
1651 SCM div
= scm_integer_expt (SCM_MAKINUM (2),
1652 SCM_MAKINUM (-bits_to_shift
));
1654 /* scm_quotient assumes its arguments are integers, but it's legal to (ash 1/2 -1) */
1655 if (SCM_FALSEP (scm_negative_p (n
)))
1656 return scm_quotient (n
, div
);
1658 return scm_sum (SCM_MAKINUM (-1L),
1659 scm_quotient (scm_sum (SCM_MAKINUM (1L), n
), div
));
1662 /* Shift left is done by multiplication with 2^CNT */
1663 return scm_product (n
, scm_integer_expt (SCM_MAKINUM (2), cnt
));
1668 SCM_DEFINE (scm_bit_extract
, "bit-extract", 3, 0, 0,
1669 (SCM n
, SCM start
, SCM end
),
1670 "Return the integer composed of the @var{start} (inclusive)\n"
1671 "through @var{end} (exclusive) bits of @var{n}. The\n"
1672 "@var{start}th bit becomes the 0-th bit in the result.\n"
1675 "(number->string (bit-extract #b1101101010 0 4) 2)\n"
1676 " @result{} \"1010\"\n"
1677 "(number->string (bit-extract #b1101101010 4 9) 2)\n"
1678 " @result{} \"10110\"\n"
1680 #define FUNC_NAME s_scm_bit_extract
1682 unsigned long int istart
, iend
;
1683 SCM_VALIDATE_INUM_MIN_COPY (2, start
,0, istart
);
1684 SCM_VALIDATE_INUM_MIN_COPY (3, end
, 0, iend
);
1685 SCM_ASSERT_RANGE (3, end
, (iend
>= istart
));
1689 long int in
= SCM_INUM (n
);
1690 unsigned long int bits
= iend
- istart
;
1692 if (in
< 0 && bits
>= SCM_I_FIXNUM_BIT
)
1694 /* Since we emulate two's complement encoded numbers, this
1695 * special case requires us to produce a result that has
1696 * more bits than can be stored in a fixnum. Thus, we fall
1697 * back to the more general algorithm that is used for
1703 if (istart
< SCM_I_FIXNUM_BIT
)
1706 if (bits
< SCM_I_FIXNUM_BIT
)
1707 return SCM_MAKINUM (in
& ((1L << bits
) - 1));
1708 else /* we know: in >= 0 */
1709 return SCM_MAKINUM (in
);
1712 return SCM_MAKINUM (-1L & ((1L << bits
) - 1));
1714 return SCM_MAKINUM (0);
1716 else if (SCM_BIGP (n
))
1720 SCM num1
= SCM_MAKINUM (1L);
1721 SCM num2
= SCM_MAKINUM (2L);
1722 SCM bits
= SCM_MAKINUM (iend
- istart
);
1723 SCM mask
= scm_difference (scm_integer_expt (num2
, bits
), num1
);
1724 return scm_logand (mask
, scm_ash (n
, SCM_MAKINUM (-istart
)));
1728 SCM_WRONG_TYPE_ARG (SCM_ARG1
, n
);
1732 static const char scm_logtab
[] = {
1733 0, 1, 1, 2, 1, 2, 2, 3, 1, 2, 2, 3, 2, 3, 3, 4
1736 SCM_DEFINE (scm_logcount
, "logcount", 1, 0, 0,
1738 "Return the number of bits in integer @var{n}. If integer is\n"
1739 "positive, the 1-bits in its binary representation are counted.\n"
1740 "If negative, the 0-bits in its two's-complement binary\n"
1741 "representation are counted. If 0, 0 is returned.\n"
1744 "(logcount #b10101010)\n"
1751 #define FUNC_NAME s_scm_logcount
1755 unsigned long int c
= 0;
1756 long int nn
= SCM_INUM (n
);
1761 c
+= scm_logtab
[15 & nn
];
1764 return SCM_MAKINUM (c
);
1766 else if (SCM_BIGP (n
))
1768 unsigned long count
;
1769 if (mpz_sgn (SCM_I_BIG_MPZ (n
)) >= 0)
1770 count
= mpz_popcount (SCM_I_BIG_MPZ (n
));
1772 count
= mpz_hamdist (SCM_I_BIG_MPZ (n
), z_negative_one
);
1773 scm_remember_upto_here_1 (n
);
1774 return SCM_MAKINUM (count
);
1777 SCM_WRONG_TYPE_ARG (SCM_ARG1
, n
);
1782 static const char scm_ilentab
[] = {
1783 0, 1, 2, 2, 3, 3, 3, 3, 4, 4, 4, 4, 4, 4, 4, 4
1787 SCM_DEFINE (scm_integer_length
, "integer-length", 1, 0, 0,
1789 "Return the number of bits necessary to represent @var{n}.\n"
1792 "(integer-length #b10101010)\n"
1794 "(integer-length 0)\n"
1796 "(integer-length #b1111)\n"
1799 #define FUNC_NAME s_scm_integer_length
1803 unsigned long int c
= 0;
1805 long int nn
= SCM_INUM (n
);
1811 l
= scm_ilentab
[15 & nn
];
1814 return SCM_MAKINUM (c
- 4 + l
);
1816 else if (SCM_BIGP (n
))
1818 /* mpz_sizeinbase looks at the absolute value of negatives, whereas we
1819 want a ones-complement. If n is ...111100..00 then mpz_sizeinbase is
1820 1 too big, so check for that and adjust. */
1821 size_t size
= mpz_sizeinbase (SCM_I_BIG_MPZ (n
), 2);
1822 if (mpz_sgn (SCM_I_BIG_MPZ (n
)) < 0
1823 && mpz_scan0 (SCM_I_BIG_MPZ (n
), /* no 0 bits above the lowest 1 */
1824 mpz_scan1 (SCM_I_BIG_MPZ (n
), 0)) == ULONG_MAX
)
1826 scm_remember_upto_here_1 (n
);
1827 return SCM_MAKINUM (size
);
1830 SCM_WRONG_TYPE_ARG (SCM_ARG1
, n
);
1834 /*** NUMBERS -> STRINGS ***/
1836 static const double fx
[] =
1837 { 0.0, 5e-1, 5e-2, 5e-3, 5e-4, 5e-5,
1838 5e-6, 5e-7, 5e-8, 5e-9, 5e-10,
1839 5e-11, 5e-12, 5e-13, 5e-14, 5e-15,
1840 5e-16, 5e-17, 5e-18, 5e-19, 5e-20};
1843 idbl2str (double f
, char *a
)
1845 int efmt
, dpt
, d
, i
, wp
= scm_dblprec
;
1851 #ifdef HAVE_COPYSIGN
1852 double sgn
= copysign (1.0, f
);
1858 goto zero
; /*{a[0]='0'; a[1]='.'; a[2]='0'; return 3;} */
1864 strcpy (a
, "-inf.0");
1866 strcpy (a
, "+inf.0");
1869 else if (xisnan (f
))
1871 strcpy (a
, "+nan.0");
1881 #ifdef DBL_MIN_10_EXP /* Prevent unnormalized values, as from
1882 make-uniform-vector, from causing infinite loops. */
1886 if (exp
-- < DBL_MIN_10_EXP
)
1897 if (exp
++ > DBL_MAX_10_EXP
)
1917 if (f
+ fx
[wp
] >= 10.0)
1924 dpt
= (exp
+ 9999) % 3;
1928 efmt
= (exp
< -3) || (exp
> wp
+ 2);
1953 if (f
+ fx
[wp
] >= 1.0)
1967 if ((dpt
> 4) && (exp
> 6))
1969 d
= (a
[0] == '-' ? 2 : 1);
1970 for (i
= ch
++; i
> d
; i
--)
1983 if (a
[ch
- 1] == '.')
1984 a
[ch
++] = '0'; /* trailing zero */
1993 for (i
= 10; i
<= exp
; i
*= 10);
1994 for (i
/= 10; i
; i
/= 10)
1996 a
[ch
++] = exp
/ i
+ '0';
2005 iflo2str (SCM flt
, char *str
)
2008 if (SCM_REALP (flt
))
2009 i
= idbl2str (SCM_REAL_VALUE (flt
), str
);
2012 i
= idbl2str (SCM_COMPLEX_REAL (flt
), str
);
2013 if (SCM_COMPLEX_IMAG (flt
) != 0.0)
2015 double imag
= SCM_COMPLEX_IMAG (flt
);
2016 /* Don't output a '+' for negative numbers or for Inf and
2017 NaN. They will provide their own sign. */
2018 if (0 <= imag
&& !xisinf (imag
) && !xisnan (imag
))
2020 i
+= idbl2str (imag
, &str
[i
]);
2027 /* convert a long to a string (unterminated). returns the number of
2028 characters in the result.
2030 p is destination: worst case (base 2) is SCM_INTBUFLEN */
2032 scm_iint2str (long num
, int rad
, char *p
)
2036 unsigned long n
= (num
< 0) ? -num
: num
;
2038 for (n
/= rad
; n
> 0; n
/= rad
)
2055 p
[i
] = d
+ ((d
< 10) ? '0' : 'a' - 10);
2060 SCM_DEFINE (scm_number_to_string
, "number->string", 1, 1, 0,
2062 "Return a string holding the external representation of the\n"
2063 "number @var{n} in the given @var{radix}. If @var{n} is\n"
2064 "inexact, a radix of 10 will be used.")
2065 #define FUNC_NAME s_scm_number_to_string
2069 if (SCM_UNBNDP (radix
))
2073 SCM_VALIDATE_INUM (2, radix
);
2074 base
= SCM_INUM (radix
);
2075 /* FIXME: ask if range limit was OK, and if so, document */
2076 SCM_ASSERT_RANGE (2, radix
, (base
>= 2) && (base
<= 36));
2081 char num_buf
[SCM_INTBUFLEN
];
2082 size_t length
= scm_iint2str (SCM_INUM (n
), base
, num_buf
);
2083 return scm_mem2string (num_buf
, length
);
2085 else if (SCM_BIGP (n
))
2087 char *str
= mpz_get_str (NULL
, base
, SCM_I_BIG_MPZ (n
));
2088 scm_remember_upto_here_1 (n
);
2089 return scm_take0str (str
);
2091 else if (SCM_FRACTIONP (n
))
2093 scm_i_fraction_reduce (n
);
2094 return scm_string_append (scm_list_3 (scm_number_to_string (SCM_FRACTION_NUMERATOR (n
), radix
),
2095 scm_mem2string ("/", 1),
2096 scm_number_to_string (SCM_FRACTION_DENOMINATOR (n
), radix
)));
2098 else if (SCM_INEXACTP (n
))
2100 char num_buf
[FLOBUFLEN
];
2101 return scm_mem2string (num_buf
, iflo2str (n
, num_buf
));
2104 SCM_WRONG_TYPE_ARG (1, n
);
2109 /* These print routines used to be stubbed here so that scm_repl.c
2110 wouldn't need SCM_BIGDIG conditionals (pre GMP) */
2113 scm_print_real (SCM sexp
, SCM port
, scm_print_state
*pstate SCM_UNUSED
)
2115 char num_buf
[FLOBUFLEN
];
2116 scm_lfwrite (num_buf
, iflo2str (sexp
, num_buf
), port
);
2121 scm_print_complex (SCM sexp
, SCM port
, scm_print_state
*pstate SCM_UNUSED
)
2124 char num_buf
[FLOBUFLEN
];
2125 scm_lfwrite (num_buf
, iflo2str (sexp
, num_buf
), port
);
2130 scm_i_print_fraction (SCM sexp
, SCM port
, scm_print_state
*pstate SCM_UNUSED
)
2133 scm_i_fraction_reduce (sexp
);
2134 str
= scm_number_to_string (sexp
, SCM_UNDEFINED
);
2135 scm_lfwrite (SCM_STRING_CHARS (str
), SCM_STRING_LENGTH (str
), port
);
2136 scm_remember_upto_here_1 (str
);
2141 scm_bigprint (SCM exp
, SCM port
, scm_print_state
*pstate SCM_UNUSED
)
2143 char *str
= mpz_get_str (NULL
, 10, SCM_I_BIG_MPZ (exp
));
2144 scm_remember_upto_here_1 (exp
);
2145 scm_lfwrite (str
, (size_t) strlen (str
), port
);
2149 /*** END nums->strs ***/
2152 /*** STRINGS -> NUMBERS ***/
2154 /* The following functions implement the conversion from strings to numbers.
2155 * The implementation somehow follows the grammar for numbers as it is given
2156 * in R5RS. Thus, the functions resemble syntactic units (<ureal R>,
2157 * <uinteger R>, ...) that are used to build up numbers in the grammar. Some
2158 * points should be noted about the implementation:
2159 * * Each function keeps a local index variable 'idx' that points at the
2160 * current position within the parsed string. The global index is only
2161 * updated if the function could parse the corresponding syntactic unit
2163 * * Similarly, the functions keep track of indicators of inexactness ('#',
2164 * '.' or exponents) using local variables ('hash_seen', 'x'). Again, the
2165 * global exactness information is only updated after each part has been
2166 * successfully parsed.
2167 * * Sequences of digits are parsed into temporary variables holding fixnums.
2168 * Only if these fixnums would overflow, the result variables are updated
2169 * using the standard functions scm_add, scm_product, scm_divide etc. Then,
2170 * the temporary variables holding the fixnums are cleared, and the process
2171 * starts over again. If for example fixnums were able to store five decimal
2172 * digits, a number 1234567890 would be parsed in two parts 12345 and 67890,
2173 * and the result was computed as 12345 * 100000 + 67890. In other words,
2174 * only every five digits two bignum operations were performed.
2177 enum t_exactness
{NO_EXACTNESS
, INEXACT
, EXACT
};
2179 /* R5RS, section 7.1.1, lexical structure of numbers: <uinteger R>. */
2181 /* In non ASCII-style encodings the following macro might not work. */
2182 #define XDIGIT2UINT(d) (isdigit (d) ? (d) - '0' : tolower (d) - 'a' + 10)
2185 mem2uinteger (const char* mem
, size_t len
, unsigned int *p_idx
,
2186 unsigned int radix
, enum t_exactness
*p_exactness
)
2188 unsigned int idx
= *p_idx
;
2189 unsigned int hash_seen
= 0;
2190 scm_t_bits shift
= 1;
2192 unsigned int digit_value
;
2202 digit_value
= XDIGIT2UINT (c
);
2203 if (digit_value
>= radix
)
2207 result
= SCM_MAKINUM (digit_value
);
2215 digit_value
= XDIGIT2UINT (c
);
2216 if (digit_value
>= radix
)
2228 if (SCM_MOST_POSITIVE_FIXNUM
/ radix
< shift
)
2230 result
= scm_product (result
, SCM_MAKINUM (shift
));
2232 result
= scm_sum (result
, SCM_MAKINUM (add
));
2239 shift
= shift
* radix
;
2240 add
= add
* radix
+ digit_value
;
2245 result
= scm_product (result
, SCM_MAKINUM (shift
));
2247 result
= scm_sum (result
, SCM_MAKINUM (add
));
2251 *p_exactness
= INEXACT
;
2257 /* R5RS, section 7.1.1, lexical structure of numbers: <decimal 10>. Only
2258 * covers the parts of the rules that start at a potential point. The value
2259 * of the digits up to the point have been parsed by the caller and are given
2260 * in variable result. The content of *p_exactness indicates, whether a hash
2261 * has already been seen in the digits before the point.
2264 /* In non ASCII-style encodings the following macro might not work. */
2265 #define DIGIT2UINT(d) ((d) - '0')
2268 mem2decimal_from_point (SCM result
, const char* mem
, size_t len
,
2269 unsigned int *p_idx
, enum t_exactness
*p_exactness
)
2271 unsigned int idx
= *p_idx
;
2272 enum t_exactness x
= *p_exactness
;
2277 if (mem
[idx
] == '.')
2279 scm_t_bits shift
= 1;
2281 unsigned int digit_value
;
2282 SCM big_shift
= SCM_MAKINUM (1);
2293 digit_value
= DIGIT2UINT (c
);
2304 if (SCM_MOST_POSITIVE_FIXNUM
/ 10 < shift
)
2306 big_shift
= scm_product (big_shift
, SCM_MAKINUM (shift
));
2307 result
= scm_product (result
, SCM_MAKINUM (shift
));
2309 result
= scm_sum (result
, SCM_MAKINUM (add
));
2317 add
= add
* 10 + digit_value
;
2323 big_shift
= scm_product (big_shift
, SCM_MAKINUM (shift
));
2324 result
= scm_product (result
, SCM_MAKINUM (shift
));
2325 result
= scm_sum (result
, SCM_MAKINUM (add
));
2328 result
= scm_divide2real (result
, big_shift
);
2330 /* We've seen a decimal point, thus the value is implicitly inexact. */
2342 /* R5RS, section 7.1.1, lexical structure of numbers: <suffix> */
2373 exponent
= DIGIT2UINT (c
);
2380 if (exponent
<= SCM_MAXEXP
)
2381 exponent
= exponent
* 10 + DIGIT2UINT (c
);
2387 if (exponent
> SCM_MAXEXP
)
2389 size_t exp_len
= idx
- start
;
2390 SCM exp_string
= scm_mem2string (&mem
[start
], exp_len
);
2391 SCM exp_num
= scm_string_to_number (exp_string
, SCM_UNDEFINED
);
2392 scm_out_of_range ("string->number", exp_num
);
2395 e
= scm_integer_expt (SCM_MAKINUM (10), SCM_MAKINUM (exponent
));
2397 result
= scm_product (result
, e
);
2399 result
= scm_divide2real (result
, e
);
2401 /* We've seen an exponent, thus the value is implicitly inexact. */
2419 /* R5RS, section 7.1.1, lexical structure of numbers: <ureal R> */
2422 mem2ureal (const char* mem
, size_t len
, unsigned int *p_idx
,
2423 unsigned int radix
, enum t_exactness
*p_exactness
)
2425 unsigned int idx
= *p_idx
;
2431 if (idx
+5 <= len
&& !strncmp (mem
+idx
, "inf.0", 5))
2437 if (idx
+4 < len
&& !strncmp (mem
+idx
, "nan.", 4))
2439 enum t_exactness x
= EXACT
;
2441 /* Cobble up the fraction. We might want to set the NaN's
2442 mantissa from it. */
2444 mem2uinteger (mem
, len
, &idx
, 10, &x
);
2449 if (mem
[idx
] == '.')
2453 else if (idx
+ 1 == len
)
2455 else if (!isdigit (mem
[idx
+ 1]))
2458 result
= mem2decimal_from_point (SCM_MAKINUM (0), mem
, len
,
2459 p_idx
, p_exactness
);
2463 enum t_exactness x
= EXACT
;
2466 uinteger
= mem2uinteger (mem
, len
, &idx
, radix
, &x
);
2467 if (SCM_FALSEP (uinteger
))
2472 else if (mem
[idx
] == '/')
2478 divisor
= mem2uinteger (mem
, len
, &idx
, radix
, &x
);
2479 if (SCM_FALSEP (divisor
))
2482 /* both are int/big here, I assume */
2483 result
= scm_make_ratio (uinteger
, divisor
);
2485 else if (radix
== 10)
2487 result
= mem2decimal_from_point (uinteger
, mem
, len
, &idx
, &x
);
2488 if (SCM_FALSEP (result
))
2499 /* When returning an inexact zero, make sure it is represented as a
2500 floating point value so that we can change its sign.
2502 if (SCM_EQ_P (result
, SCM_MAKINUM(0)) && *p_exactness
== INEXACT
)
2503 result
= scm_make_real (0.0);
2509 /* R5RS, section 7.1.1, lexical structure of numbers: <complex R> */
2512 mem2complex (const char* mem
, size_t len
, unsigned int idx
,
2513 unsigned int radix
, enum t_exactness
*p_exactness
)
2537 ureal
= mem2ureal (mem
, len
, &idx
, radix
, p_exactness
);
2538 if (SCM_FALSEP (ureal
))
2540 /* input must be either +i or -i */
2545 if (mem
[idx
] == 'i' || mem
[idx
] == 'I')
2551 return scm_make_rectangular (SCM_MAKINUM (0), SCM_MAKINUM (sign
));
2558 if (sign
== -1 && SCM_FALSEP (scm_nan_p (ureal
)))
2559 ureal
= scm_difference (ureal
, SCM_UNDEFINED
);
2568 /* either +<ureal>i or -<ureal>i */
2575 return scm_make_rectangular (SCM_MAKINUM (0), ureal
);
2578 /* polar input: <real>@<real>. */
2603 angle
= mem2ureal (mem
, len
, &idx
, radix
, p_exactness
);
2604 if (SCM_FALSEP (angle
))
2609 if (sign
== -1 && SCM_FALSEP (scm_nan_p (ureal
)))
2610 angle
= scm_difference (angle
, SCM_UNDEFINED
);
2612 result
= scm_make_polar (ureal
, angle
);
2617 /* expecting input matching <real>[+-]<ureal>?i */
2624 int sign
= (c
== '+') ? 1 : -1;
2625 SCM imag
= mem2ureal (mem
, len
, &idx
, radix
, p_exactness
);
2627 if (SCM_FALSEP (imag
))
2628 imag
= SCM_MAKINUM (sign
);
2629 else if (sign
== -1 && SCM_FALSEP (scm_nan_p (ureal
)))
2630 imag
= scm_difference (imag
, SCM_UNDEFINED
);
2634 if (mem
[idx
] != 'i' && mem
[idx
] != 'I')
2641 return scm_make_rectangular (ureal
, imag
);
2650 /* R5RS, section 7.1.1, lexical structure of numbers: <number> */
2652 enum t_radix
{NO_RADIX
=0, DUAL
=2, OCT
=8, DEC
=10, HEX
=16};
2655 scm_i_mem2number (const char* mem
, size_t len
, unsigned int default_radix
)
2657 unsigned int idx
= 0;
2658 unsigned int radix
= NO_RADIX
;
2659 enum t_exactness forced_x
= NO_EXACTNESS
;
2660 enum t_exactness implicit_x
= EXACT
;
2663 /* R5RS, section 7.1.1, lexical structure of numbers: <prefix R> */
2664 while (idx
+ 2 < len
&& mem
[idx
] == '#')
2666 switch (mem
[idx
+ 1])
2669 if (radix
!= NO_RADIX
)
2674 if (radix
!= NO_RADIX
)
2679 if (forced_x
!= NO_EXACTNESS
)
2684 if (forced_x
!= NO_EXACTNESS
)
2689 if (radix
!= NO_RADIX
)
2694 if (radix
!= NO_RADIX
)
2704 /* R5RS, section 7.1.1, lexical structure of numbers: <complex R> */
2705 if (radix
== NO_RADIX
)
2706 result
= mem2complex (mem
, len
, idx
, default_radix
, &implicit_x
);
2708 result
= mem2complex (mem
, len
, idx
, (unsigned int) radix
, &implicit_x
);
2710 if (SCM_FALSEP (result
))
2716 if (SCM_INEXACTP (result
))
2717 /* FIXME: This may change the value. */
2718 return scm_inexact_to_exact (result
);
2722 if (SCM_INEXACTP (result
))
2725 return scm_exact_to_inexact (result
);
2728 if (implicit_x
== INEXACT
)
2730 if (SCM_INEXACTP (result
))
2733 return scm_exact_to_inexact (result
);
2741 SCM_DEFINE (scm_string_to_number
, "string->number", 1, 1, 0,
2742 (SCM string
, SCM radix
),
2743 "Return a number of the maximally precise representation\n"
2744 "expressed by the given @var{string}. @var{radix} must be an\n"
2745 "exact integer, either 2, 8, 10, or 16. If supplied, @var{radix}\n"
2746 "is a default radix that may be overridden by an explicit radix\n"
2747 "prefix in @var{string} (e.g. \"#o177\"). If @var{radix} is not\n"
2748 "supplied, then the default radix is 10. If string is not a\n"
2749 "syntactically valid notation for a number, then\n"
2750 "@code{string->number} returns @code{#f}.")
2751 #define FUNC_NAME s_scm_string_to_number
2755 SCM_VALIDATE_STRING (1, string
);
2756 SCM_VALIDATE_INUM_MIN_DEF_COPY (2, radix
,2,10, base
);
2757 answer
= scm_i_mem2number (SCM_STRING_CHARS (string
),
2758 SCM_STRING_LENGTH (string
),
2760 return scm_return_first (answer
, string
);
2765 /*** END strs->nums ***/
2769 scm_make_real (double x
)
2771 SCM z
= scm_double_cell (scm_tc16_real
, 0, 0, 0);
2773 SCM_REAL_VALUE (z
) = x
;
2779 scm_make_complex (double x
, double y
)
2782 return scm_make_real (x
);
2786 SCM_NEWSMOB (z
, scm_tc16_complex
, scm_gc_malloc (sizeof (scm_t_complex
),
2788 SCM_COMPLEX_REAL (z
) = x
;
2789 SCM_COMPLEX_IMAG (z
) = y
;
2796 scm_bigequal (SCM x
, SCM y
)
2798 int result
= mpz_cmp (SCM_I_BIG_MPZ (x
), SCM_I_BIG_MPZ (y
));
2799 scm_remember_upto_here_2 (x
, y
);
2800 return SCM_BOOL (0 == result
);
2804 scm_real_equalp (SCM x
, SCM y
)
2806 return SCM_BOOL (SCM_REAL_VALUE (x
) == SCM_REAL_VALUE (y
));
2810 scm_complex_equalp (SCM x
, SCM y
)
2812 return SCM_BOOL (SCM_COMPLEX_REAL (x
) == SCM_COMPLEX_REAL (y
)
2813 && SCM_COMPLEX_IMAG (x
) == SCM_COMPLEX_IMAG (y
));
2817 scm_i_fraction_equalp (SCM x
, SCM y
)
2819 scm_i_fraction_reduce (x
);
2820 scm_i_fraction_reduce (y
);
2821 return SCM_BOOL (scm_equal_p (SCM_FRACTION_NUMERATOR (x
), SCM_FRACTION_NUMERATOR (y
))
2822 && scm_equal_p (SCM_FRACTION_DENOMINATOR (x
), SCM_FRACTION_DENOMINATOR (y
)));
2826 SCM_REGISTER_PROC (s_number_p
, "number?", 1, 0, 0, scm_number_p
);
2827 /* "Return @code{#t} if @var{x} is a number, @code{#f}\n"
2828 * "else. Note that the sets of complex, real, rational and\n"
2829 * "integer values form subsets of the set of numbers, i. e. the\n"
2830 * "predicate will be fulfilled for any number."
2832 SCM_DEFINE (scm_number_p
, "complex?", 1, 0, 0,
2834 "Return @code{#t} if @var{x} is a complex number, @code{#f}\n"
2835 "otherwise. Note that the sets of real, rational and integer\n"
2836 "values form subsets of the set of complex numbers, i. e. the\n"
2837 "predicate will also be fulfilled if @var{x} is a real,\n"
2838 "rational or integer number.")
2839 #define FUNC_NAME s_scm_number_p
2841 return SCM_BOOL (SCM_NUMBERP (x
));
2846 SCM_DEFINE (scm_real_p
, "real?", 1, 0, 0,
2848 "Return @code{#t} if @var{x} is a real number, @code{#f}\n"
2849 "otherwise. Note that the set of integer values forms a subset of\n"
2850 "the set of real numbers, i. e. the predicate will also be\n"
2851 "fulfilled if @var{x} is an integer number.")
2852 #define FUNC_NAME s_scm_real_p
2854 /* we can't represent irrational numbers. */
2855 return scm_rational_p (x
);
2859 SCM_DEFINE (scm_rational_p
, "rational?", 1, 0, 0,
2861 "Return @code{#t} if @var{x} is a rational number, @code{#f}\n"
2862 "otherwise. Note that the set of integer values forms a subset of\n"
2863 "the set of rational numbers, i. e. the predicate will also be\n"
2864 "fulfilled if @var{x} is an integer number.")
2865 #define FUNC_NAME s_scm_rational_p
2869 else if (SCM_IMP (x
))
2871 else if (SCM_BIGP (x
))
2873 else if (SCM_FRACTIONP (x
))
2875 else if (SCM_REALP (x
))
2876 /* due to their limited precision, all floating point numbers are
2877 rational as well. */
2885 SCM_DEFINE (scm_integer_p
, "integer?", 1, 0, 0,
2887 "Return @code{#t} if @var{x} is an integer number, @code{#f}\n"
2889 #define FUNC_NAME s_scm_integer_p
2898 if (!SCM_INEXACTP (x
))
2900 if (SCM_COMPLEXP (x
))
2902 r
= SCM_REAL_VALUE (x
);
2910 SCM_DEFINE (scm_inexact_p
, "inexact?", 1, 0, 0,
2912 "Return @code{#t} if @var{x} is an inexact number, @code{#f}\n"
2914 #define FUNC_NAME s_scm_inexact_p
2916 return SCM_BOOL (SCM_INEXACTP (x
));
2921 SCM_GPROC1 (s_eq_p
, "=", scm_tc7_rpsubr
, scm_num_eq_p
, g_eq_p
);
2922 /* "Return @code{#t} if all parameters are numerically equal." */
2924 scm_num_eq_p (SCM x
, SCM y
)
2928 long xx
= SCM_INUM (x
);
2931 long yy
= SCM_INUM (y
);
2932 return SCM_BOOL (xx
== yy
);
2934 else if (SCM_BIGP (y
))
2936 else if (SCM_REALP (y
))
2937 return SCM_BOOL ((double) xx
== SCM_REAL_VALUE (y
));
2938 else if (SCM_COMPLEXP (y
))
2939 return SCM_BOOL (((double) xx
== SCM_COMPLEX_REAL (y
))
2940 && (0.0 == SCM_COMPLEX_IMAG (y
)));
2941 else if (SCM_FRACTIONP (y
))
2944 SCM_WTA_DISPATCH_2 (g_eq_p
, x
, y
, SCM_ARGn
, s_eq_p
);
2946 else if (SCM_BIGP (x
))
2950 else if (SCM_BIGP (y
))
2952 int cmp
= mpz_cmp (SCM_I_BIG_MPZ (x
), SCM_I_BIG_MPZ (y
));
2953 scm_remember_upto_here_2 (x
, y
);
2954 return SCM_BOOL (0 == cmp
);
2956 else if (SCM_REALP (y
))
2959 if (xisnan (SCM_REAL_VALUE (y
)))
2961 cmp
= xmpz_cmp_d (SCM_I_BIG_MPZ (x
), SCM_REAL_VALUE (y
));
2962 scm_remember_upto_here_1 (x
);
2963 return SCM_BOOL (0 == cmp
);
2965 else if (SCM_COMPLEXP (y
))
2968 if (0.0 != SCM_COMPLEX_IMAG (y
))
2970 if (xisnan (SCM_COMPLEX_REAL (y
)))
2972 cmp
= xmpz_cmp_d (SCM_I_BIG_MPZ (x
), SCM_COMPLEX_REAL (y
));
2973 scm_remember_upto_here_1 (x
);
2974 return SCM_BOOL (0 == cmp
);
2976 else if (SCM_FRACTIONP (y
))
2979 SCM_WTA_DISPATCH_2 (g_eq_p
, x
, y
, SCM_ARGn
, s_eq_p
);
2981 else if (SCM_REALP (x
))
2984 return SCM_BOOL (SCM_REAL_VALUE (x
) == (double) SCM_INUM (y
));
2985 else if (SCM_BIGP (y
))
2988 if (xisnan (SCM_REAL_VALUE (x
)))
2990 cmp
= xmpz_cmp_d (SCM_I_BIG_MPZ (y
), SCM_REAL_VALUE (x
));
2991 scm_remember_upto_here_1 (y
);
2992 return SCM_BOOL (0 == cmp
);
2994 else if (SCM_REALP (y
))
2995 return SCM_BOOL (SCM_REAL_VALUE (x
) == SCM_REAL_VALUE (y
));
2996 else if (SCM_COMPLEXP (y
))
2997 return SCM_BOOL ((SCM_REAL_VALUE (x
) == SCM_COMPLEX_REAL (y
))
2998 && (0.0 == SCM_COMPLEX_IMAG (y
)));
2999 else if (SCM_FRACTIONP (y
))
3000 return SCM_BOOL (SCM_REAL_VALUE (x
) == scm_i_fraction2double (y
));
3002 SCM_WTA_DISPATCH_2 (g_eq_p
, x
, y
, SCM_ARGn
, s_eq_p
);
3004 else if (SCM_COMPLEXP (x
))
3007 return SCM_BOOL ((SCM_COMPLEX_REAL (x
) == (double) SCM_INUM (y
))
3008 && (SCM_COMPLEX_IMAG (x
) == 0.0));
3009 else if (SCM_BIGP (y
))
3012 if (0.0 != SCM_COMPLEX_IMAG (x
))
3014 if (xisnan (SCM_COMPLEX_REAL (x
)))
3016 cmp
= xmpz_cmp_d (SCM_I_BIG_MPZ (y
), SCM_COMPLEX_REAL (x
));
3017 scm_remember_upto_here_1 (y
);
3018 return SCM_BOOL (0 == cmp
);
3020 else if (SCM_REALP (y
))
3021 return SCM_BOOL ((SCM_COMPLEX_REAL (x
) == SCM_REAL_VALUE (y
))
3022 && (SCM_COMPLEX_IMAG (x
) == 0.0));
3023 else if (SCM_COMPLEXP (y
))
3024 return SCM_BOOL ((SCM_COMPLEX_REAL (x
) == SCM_COMPLEX_REAL (y
))
3025 && (SCM_COMPLEX_IMAG (x
) == SCM_COMPLEX_IMAG (y
)));
3026 else if (SCM_FRACTIONP (y
))
3027 return SCM_BOOL ((SCM_COMPLEX_REAL (x
) == scm_i_fraction2double (y
))
3028 && (SCM_COMPLEX_IMAG (x
) == 0.0));
3030 SCM_WTA_DISPATCH_2 (g_eq_p
, x
, y
, SCM_ARGn
, s_eq_p
);
3032 else if (SCM_FRACTIONP (x
))
3036 else if (SCM_BIGP (y
))
3038 else if (SCM_REALP (y
))
3039 return SCM_BOOL (scm_i_fraction2double (x
) == SCM_REAL_VALUE (y
));
3040 else if (SCM_COMPLEXP (y
))
3041 return SCM_BOOL ((scm_i_fraction2double (x
) == SCM_COMPLEX_REAL (y
))
3042 && (0.0 == SCM_COMPLEX_IMAG (y
)));
3043 else if (SCM_FRACTIONP (y
))
3044 return scm_i_fraction_equalp (x
, y
);
3046 SCM_WTA_DISPATCH_2 (g_eq_p
, x
, y
, SCM_ARGn
, s_eq_p
);
3049 SCM_WTA_DISPATCH_2 (g_eq_p
, x
, y
, SCM_ARG1
, s_eq_p
);
3053 SCM_GPROC1 (s_less_p
, "<", scm_tc7_rpsubr
, scm_less_p
, g_less_p
);
3054 /* "Return @code{#t} if the list of parameters is monotonically\n"
3058 scm_less_p (SCM x
, SCM y
)
3062 long xx
= SCM_INUM (x
);
3065 long yy
= SCM_INUM (y
);
3066 return SCM_BOOL (xx
< yy
);
3068 else if (SCM_BIGP (y
))
3070 int sgn
= mpz_sgn (SCM_I_BIG_MPZ (y
));
3071 scm_remember_upto_here_1 (y
);
3072 return SCM_BOOL (sgn
> 0);
3074 else if (SCM_REALP (y
))
3075 return SCM_BOOL ((double) xx
< SCM_REAL_VALUE (y
));
3076 else if (SCM_FRACTIONP (y
))
3077 return SCM_BOOL ((double) xx
< scm_i_fraction2double (y
));
3079 SCM_WTA_DISPATCH_2 (g_less_p
, x
, y
, SCM_ARGn
, s_less_p
);
3081 else if (SCM_BIGP (x
))
3085 int sgn
= mpz_sgn (SCM_I_BIG_MPZ (x
));
3086 scm_remember_upto_here_1 (x
);
3087 return SCM_BOOL (sgn
< 0);
3089 else if (SCM_BIGP (y
))
3091 int cmp
= mpz_cmp (SCM_I_BIG_MPZ (x
), SCM_I_BIG_MPZ (y
));
3092 scm_remember_upto_here_2 (x
, y
);
3093 return SCM_BOOL (cmp
< 0);
3095 else if (SCM_REALP (y
))
3098 if (xisnan (SCM_REAL_VALUE (y
)))
3100 cmp
= xmpz_cmp_d (SCM_I_BIG_MPZ (x
), SCM_REAL_VALUE (y
));
3101 scm_remember_upto_here_1 (x
);
3102 return SCM_BOOL (cmp
< 0);
3104 else if (SCM_FRACTIONP (y
))
3107 cmp
= xmpz_cmp_d (SCM_I_BIG_MPZ (x
), scm_i_fraction2double (y
));
3108 scm_remember_upto_here_1 (x
);
3109 return SCM_BOOL (cmp
< 0);
3112 SCM_WTA_DISPATCH_2 (g_less_p
, x
, y
, SCM_ARGn
, s_less_p
);
3114 else if (SCM_REALP (x
))
3117 return SCM_BOOL (SCM_REAL_VALUE (x
) < (double) SCM_INUM (y
));
3118 else if (SCM_BIGP (y
))
3121 if (xisnan (SCM_REAL_VALUE (x
)))
3123 cmp
= xmpz_cmp_d (SCM_I_BIG_MPZ (y
), SCM_REAL_VALUE (x
));
3124 scm_remember_upto_here_1 (y
);
3125 return SCM_BOOL (cmp
> 0);
3127 else if (SCM_REALP (y
))
3128 return SCM_BOOL (SCM_REAL_VALUE (x
) < SCM_REAL_VALUE (y
));
3129 else if (SCM_FRACTIONP (y
))
3130 return SCM_BOOL (SCM_REAL_VALUE (x
) < scm_i_fraction2double (y
));
3132 SCM_WTA_DISPATCH_2 (g_less_p
, x
, y
, SCM_ARGn
, s_less_p
);
3134 else if (SCM_FRACTIONP (x
))
3137 return SCM_BOOL (scm_i_fraction2double (x
) < (double) SCM_INUM (y
));
3138 else if (SCM_BIGP (y
))
3141 if (xisnan (SCM_REAL_VALUE (x
)))
3143 cmp
= xmpz_cmp_d (SCM_I_BIG_MPZ (y
), scm_i_fraction2double (x
));
3144 scm_remember_upto_here_1 (y
);
3145 return SCM_BOOL (cmp
> 0);
3147 else if (SCM_REALP (y
))
3148 return SCM_BOOL (scm_i_fraction2double (x
) < SCM_REAL_VALUE (y
));
3149 else if (SCM_FRACTIONP (y
))
3150 return SCM_BOOL (scm_i_fraction2double (x
) < scm_i_fraction2double (y
));
3152 SCM_WTA_DISPATCH_2 (g_less_p
, x
, y
, SCM_ARGn
, s_less_p
);
3155 SCM_WTA_DISPATCH_2 (g_less_p
, x
, y
, SCM_ARG1
, s_less_p
);
3159 SCM_GPROC1 (s_scm_gr_p
, ">", scm_tc7_rpsubr
, scm_gr_p
, g_gr_p
);
3160 /* "Return @code{#t} if the list of parameters is monotonically\n"
3163 #define FUNC_NAME s_scm_gr_p
3165 scm_gr_p (SCM x
, SCM y
)
3167 if (!SCM_NUMBERP (x
))
3168 SCM_WTA_DISPATCH_2 (g_gr_p
, x
, y
, SCM_ARG1
, FUNC_NAME
);
3169 else if (!SCM_NUMBERP (y
))
3170 SCM_WTA_DISPATCH_2 (g_gr_p
, x
, y
, SCM_ARG2
, FUNC_NAME
);
3172 return scm_less_p (y
, x
);
3177 SCM_GPROC1 (s_scm_leq_p
, "<=", scm_tc7_rpsubr
, scm_leq_p
, g_leq_p
);
3178 /* "Return @code{#t} if the list of parameters is monotonically\n"
3181 #define FUNC_NAME s_scm_leq_p
3183 scm_leq_p (SCM x
, SCM y
)
3185 if (!SCM_NUMBERP (x
))
3186 SCM_WTA_DISPATCH_2 (g_leq_p
, x
, y
, SCM_ARG1
, FUNC_NAME
);
3187 else if (!SCM_NUMBERP (y
))
3188 SCM_WTA_DISPATCH_2 (g_leq_p
, x
, y
, SCM_ARG2
, FUNC_NAME
);
3189 else if (SCM_NFALSEP (scm_nan_p (x
)) || SCM_NFALSEP (scm_nan_p (y
)))
3192 return SCM_BOOL_NOT (scm_less_p (y
, x
));
3197 SCM_GPROC1 (s_scm_geq_p
, ">=", scm_tc7_rpsubr
, scm_geq_p
, g_geq_p
);
3198 /* "Return @code{#t} if the list of parameters is monotonically\n"
3201 #define FUNC_NAME s_scm_geq_p
3203 scm_geq_p (SCM x
, SCM y
)
3205 if (!SCM_NUMBERP (x
))
3206 SCM_WTA_DISPATCH_2 (g_geq_p
, x
, y
, SCM_ARG1
, FUNC_NAME
);
3207 else if (!SCM_NUMBERP (y
))
3208 SCM_WTA_DISPATCH_2 (g_geq_p
, x
, y
, SCM_ARG2
, FUNC_NAME
);
3209 else if (SCM_NFALSEP (scm_nan_p (x
)) || SCM_NFALSEP (scm_nan_p (y
)))
3212 return SCM_BOOL_NOT (scm_less_p (x
, y
));
3217 SCM_GPROC (s_zero_p
, "zero?", 1, 0, 0, scm_zero_p
, g_zero_p
);
3218 /* "Return @code{#t} if @var{z} is an exact or inexact number equal to\n"
3225 return SCM_BOOL (SCM_EQ_P (z
, SCM_INUM0
));
3226 else if (SCM_BIGP (z
))
3228 else if (SCM_REALP (z
))
3229 return SCM_BOOL (SCM_REAL_VALUE (z
) == 0.0);
3230 else if (SCM_COMPLEXP (z
))
3231 return SCM_BOOL (SCM_COMPLEX_REAL (z
) == 0.0
3232 && SCM_COMPLEX_IMAG (z
) == 0.0);
3233 else if (SCM_FRACTIONP (z
))
3236 SCM_WTA_DISPATCH_1 (g_zero_p
, z
, SCM_ARG1
, s_zero_p
);
3240 SCM_GPROC (s_positive_p
, "positive?", 1, 0, 0, scm_positive_p
, g_positive_p
);
3241 /* "Return @code{#t} if @var{x} is an exact or inexact number greater than\n"
3245 scm_positive_p (SCM x
)
3248 return SCM_BOOL (SCM_INUM (x
) > 0);
3249 else if (SCM_BIGP (x
))
3251 int sgn
= mpz_sgn (SCM_I_BIG_MPZ (x
));
3252 scm_remember_upto_here_1 (x
);
3253 return SCM_BOOL (sgn
> 0);
3255 else if (SCM_REALP (x
))
3256 return SCM_BOOL(SCM_REAL_VALUE (x
) > 0.0);
3257 else if (SCM_FRACTIONP (x
))
3258 return scm_positive_p (SCM_FRACTION_NUMERATOR (x
));
3260 SCM_WTA_DISPATCH_1 (g_positive_p
, x
, SCM_ARG1
, s_positive_p
);
3264 SCM_GPROC (s_negative_p
, "negative?", 1, 0, 0, scm_negative_p
, g_negative_p
);
3265 /* "Return @code{#t} if @var{x} is an exact or inexact number less than\n"
3269 scm_negative_p (SCM x
)
3272 return SCM_BOOL (SCM_INUM (x
) < 0);
3273 else if (SCM_BIGP (x
))
3275 int sgn
= mpz_sgn (SCM_I_BIG_MPZ (x
));
3276 scm_remember_upto_here_1 (x
);
3277 return SCM_BOOL (sgn
< 0);
3279 else if (SCM_REALP (x
))
3280 return SCM_BOOL(SCM_REAL_VALUE (x
) < 0.0);
3281 else if (SCM_FRACTIONP (x
))
3282 return scm_negative_p (SCM_FRACTION_NUMERATOR (x
));
3284 SCM_WTA_DISPATCH_1 (g_negative_p
, x
, SCM_ARG1
, s_negative_p
);
3288 SCM_GPROC1 (s_max
, "max", scm_tc7_asubr
, scm_max
, g_max
);
3289 /* "Return the maximum of all parameter values."
3292 scm_max (SCM x
, SCM y
)
3297 SCM_WTA_DISPATCH_0 (g_max
, s_max
);
3298 else if (SCM_NUMBERP (x
))
3301 SCM_WTA_DISPATCH_1 (g_max
, x
, SCM_ARG1
, s_max
);
3306 long xx
= SCM_INUM (x
);
3309 long yy
= SCM_INUM (y
);
3310 return (xx
< yy
) ? y
: x
;
3312 else if (SCM_BIGP (y
))
3314 int sgn
= mpz_sgn (SCM_I_BIG_MPZ (y
));
3315 scm_remember_upto_here_1 (y
);
3316 return (sgn
< 0) ? x
: y
;
3318 else if (SCM_REALP (y
))
3321 /* if y==NaN then ">" is false and we return NaN */
3322 return (z
> SCM_REAL_VALUE (y
)) ? scm_make_real (z
) : y
;
3324 else if (SCM_FRACTIONP (y
))
3327 return (z
> scm_i_fraction2double (y
)) ? x
: y
;
3330 SCM_WTA_DISPATCH_2 (g_max
, x
, y
, SCM_ARGn
, s_max
);
3332 else if (SCM_BIGP (x
))
3336 int sgn
= mpz_sgn (SCM_I_BIG_MPZ (x
));
3337 scm_remember_upto_here_1 (x
);
3338 return (sgn
< 0) ? y
: x
;
3340 else if (SCM_BIGP (y
))
3342 int cmp
= mpz_cmp (SCM_I_BIG_MPZ (x
), SCM_I_BIG_MPZ (y
));
3343 scm_remember_upto_here_2 (x
, y
);
3344 return (cmp
> 0) ? x
: y
;
3346 else if (SCM_REALP (y
))
3348 double yy
= SCM_REAL_VALUE (y
);
3352 cmp
= xmpz_cmp_d (SCM_I_BIG_MPZ (x
), yy
);
3353 scm_remember_upto_here_1 (x
);
3354 return (cmp
> 0) ? x
: y
;
3356 else if (SCM_FRACTIONP (y
))
3358 double yy
= scm_i_fraction2double (y
);
3360 cmp
= xmpz_cmp_d (SCM_I_BIG_MPZ (x
), yy
);
3361 scm_remember_upto_here_1 (x
);
3362 return (cmp
> 0) ? x
: y
;
3365 SCM_WTA_DISPATCH_2 (g_max
, x
, y
, SCM_ARGn
, s_max
);
3367 else if (SCM_REALP (x
))
3371 double z
= SCM_INUM (y
);
3372 /* if x==NaN then "<" is false and we return NaN */
3373 return (SCM_REAL_VALUE (x
) < z
) ? scm_make_real (z
) : x
;
3375 else if (SCM_BIGP (y
))
3377 double xx
= SCM_REAL_VALUE (x
);
3381 cmp
= xmpz_cmp_d (SCM_I_BIG_MPZ (y
), xx
);
3382 scm_remember_upto_here_1 (y
);
3383 return (cmp
< 0) ? x
: y
;
3385 else if (SCM_REALP (y
))
3387 /* if x==NaN then our explicit check means we return NaN
3388 if y==NaN then ">" is false and we return NaN
3389 calling isnan is unavoidable, since it's the only way to know
3390 which of x or y causes any compares to be false */
3391 double xx
= SCM_REAL_VALUE (x
);
3392 return (xisnan (xx
) || xx
> SCM_REAL_VALUE (y
)) ? x
: y
;
3394 else if (SCM_FRACTIONP (y
))
3396 double yy
= scm_i_fraction2double (y
);
3397 double xx
= SCM_REAL_VALUE (x
);
3398 return (xx
< yy
) ? scm_make_real (yy
) : x
;
3401 SCM_WTA_DISPATCH_2 (g_max
, x
, y
, SCM_ARGn
, s_max
);
3403 else if (SCM_FRACTIONP (x
))
3407 double z
= SCM_INUM (y
);
3408 return (scm_i_fraction2double (x
) < z
) ? y
: x
;
3410 else if (SCM_BIGP (y
))
3412 double xx
= scm_i_fraction2double (x
);
3414 cmp
= xmpz_cmp_d (SCM_I_BIG_MPZ (y
), xx
);
3415 scm_remember_upto_here_1 (y
);
3416 return (cmp
< 0) ? x
: y
;
3418 else if (SCM_REALP (y
))
3420 double xx
= scm_i_fraction2double (x
);
3421 return (xx
< SCM_REAL_VALUE (y
)) ? y
: scm_make_real (xx
);
3423 else if (SCM_FRACTIONP (y
))
3425 double yy
= scm_i_fraction2double (y
);
3426 double xx
= scm_i_fraction2double (x
);
3427 return (xx
< yy
) ? y
: x
;
3430 SCM_WTA_DISPATCH_2 (g_max
, x
, y
, SCM_ARGn
, s_max
);
3433 SCM_WTA_DISPATCH_2 (g_max
, x
, y
, SCM_ARG1
, s_max
);
3437 SCM_GPROC1 (s_min
, "min", scm_tc7_asubr
, scm_min
, g_min
);
3438 /* "Return the minium of all parameter values."
3441 scm_min (SCM x
, SCM y
)
3446 SCM_WTA_DISPATCH_0 (g_min
, s_min
);
3447 else if (SCM_NUMBERP (x
))
3450 SCM_WTA_DISPATCH_1 (g_min
, x
, SCM_ARG1
, s_min
);
3455 long xx
= SCM_INUM (x
);
3458 long yy
= SCM_INUM (y
);
3459 return (xx
< yy
) ? x
: y
;
3461 else if (SCM_BIGP (y
))
3463 int sgn
= mpz_sgn (SCM_I_BIG_MPZ (y
));
3464 scm_remember_upto_here_1 (y
);
3465 return (sgn
< 0) ? y
: x
;
3467 else if (SCM_REALP (y
))
3470 /* if y==NaN then "<" is false and we return NaN */
3471 return (z
< SCM_REAL_VALUE (y
)) ? scm_make_real (z
) : y
;
3473 else if (SCM_FRACTIONP (y
))
3476 return (z
< scm_i_fraction2double (y
)) ? x
: y
;
3479 SCM_WTA_DISPATCH_2 (g_min
, x
, y
, SCM_ARGn
, s_min
);
3481 else if (SCM_BIGP (x
))
3485 int sgn
= mpz_sgn (SCM_I_BIG_MPZ (x
));
3486 scm_remember_upto_here_1 (x
);
3487 return (sgn
< 0) ? x
: y
;
3489 else if (SCM_BIGP (y
))
3491 int cmp
= mpz_cmp (SCM_I_BIG_MPZ (x
), SCM_I_BIG_MPZ (y
));
3492 scm_remember_upto_here_2 (x
, y
);
3493 return (cmp
> 0) ? y
: x
;
3495 else if (SCM_REALP (y
))
3497 double yy
= SCM_REAL_VALUE (y
);
3501 cmp
= xmpz_cmp_d (SCM_I_BIG_MPZ (x
), yy
);
3502 scm_remember_upto_here_1 (x
);
3503 return (cmp
> 0) ? y
: x
;
3505 else if (SCM_FRACTIONP (y
))
3507 double yy
= scm_i_fraction2double (y
);
3509 cmp
= xmpz_cmp_d (SCM_I_BIG_MPZ (x
), yy
);
3510 scm_remember_upto_here_1 (x
);
3511 return (cmp
> 0) ? y
: x
;
3514 SCM_WTA_DISPATCH_2 (g_min
, x
, y
, SCM_ARGn
, s_min
);
3516 else if (SCM_REALP (x
))
3520 double z
= SCM_INUM (y
);
3521 /* if x==NaN then "<" is false and we return NaN */
3522 return (z
< SCM_REAL_VALUE (x
)) ? scm_make_real (z
) : x
;
3524 else if (SCM_BIGP (y
))
3526 double xx
= SCM_REAL_VALUE (x
);
3530 cmp
= xmpz_cmp_d (SCM_I_BIG_MPZ (y
), xx
);
3531 scm_remember_upto_here_1 (y
);
3532 return (cmp
< 0) ? y
: x
;
3534 else if (SCM_REALP (y
))
3536 /* if x==NaN then our explicit check means we return NaN
3537 if y==NaN then "<" is false and we return NaN
3538 calling isnan is unavoidable, since it's the only way to know
3539 which of x or y causes any compares to be false */
3540 double xx
= SCM_REAL_VALUE (x
);
3541 return (xisnan (xx
) || xx
< SCM_REAL_VALUE (y
)) ? x
: y
;
3543 else if (SCM_FRACTIONP (y
))
3545 double yy
= scm_i_fraction2double (y
);
3546 double xx
= SCM_REAL_VALUE (x
);
3547 return (yy
< xx
) ? scm_make_real (yy
) : x
;
3550 SCM_WTA_DISPATCH_2 (g_min
, x
, y
, SCM_ARGn
, s_min
);
3552 else if (SCM_FRACTIONP (x
))
3556 double z
= SCM_INUM (y
);
3557 return (scm_i_fraction2double (x
) < z
) ? x
: y
;
3559 else if (SCM_BIGP (y
))
3561 double xx
= scm_i_fraction2double (x
);
3563 cmp
= xmpz_cmp_d (SCM_I_BIG_MPZ (y
), xx
);
3564 scm_remember_upto_here_1 (y
);
3565 return (cmp
< 0) ? y
: x
;
3567 else if (SCM_REALP (y
))
3569 double xx
= scm_i_fraction2double (x
);
3570 return (SCM_REAL_VALUE (y
) < xx
) ? y
: scm_make_real (xx
);
3572 else if (SCM_FRACTIONP (y
))
3574 double yy
= scm_i_fraction2double (y
);
3575 double xx
= scm_i_fraction2double (x
);
3576 return (xx
< yy
) ? x
: y
;
3579 SCM_WTA_DISPATCH_2 (g_max
, x
, y
, SCM_ARGn
, s_max
);
3582 SCM_WTA_DISPATCH_2 (g_min
, x
, y
, SCM_ARG1
, s_min
);
3586 SCM_GPROC1 (s_sum
, "+", scm_tc7_asubr
, scm_sum
, g_sum
);
3587 /* "Return the sum of all parameter values. Return 0 if called without\n"
3591 scm_sum (SCM x
, SCM y
)
3595 if (SCM_NUMBERP (x
)) return x
;
3596 if (SCM_UNBNDP (x
)) return SCM_INUM0
;
3597 SCM_WTA_DISPATCH_1 (g_sum
, x
, SCM_ARG1
, s_sum
);
3604 long xx
= SCM_INUM (x
);
3605 long yy
= SCM_INUM (y
);
3606 long int z
= xx
+ yy
;
3607 return SCM_FIXABLE (z
) ? SCM_MAKINUM (z
) : scm_i_long2big (z
);
3609 else if (SCM_BIGP (y
))
3614 else if (SCM_REALP (y
))
3616 long int xx
= SCM_INUM (x
);
3617 return scm_make_real (xx
+ SCM_REAL_VALUE (y
));
3619 else if (SCM_COMPLEXP (y
))
3621 long int xx
= SCM_INUM (x
);
3622 return scm_make_complex (xx
+ SCM_COMPLEX_REAL (y
),
3623 SCM_COMPLEX_IMAG (y
));
3625 else if (SCM_FRACTIONP (y
))
3626 return scm_make_ratio (scm_sum (SCM_FRACTION_NUMERATOR (y
),
3627 scm_product (x
, SCM_FRACTION_DENOMINATOR (y
))),
3628 SCM_FRACTION_DENOMINATOR (y
));
3630 SCM_WTA_DISPATCH_2 (g_sum
, x
, y
, SCM_ARGn
, s_sum
);
3631 } else if (SCM_BIGP (x
))
3638 inum
= SCM_INUM (y
);
3641 bigsgn
= mpz_sgn (SCM_I_BIG_MPZ (x
));
3644 SCM result
= scm_i_mkbig ();
3645 mpz_sub_ui (SCM_I_BIG_MPZ (result
), SCM_I_BIG_MPZ (x
), - inum
);
3646 scm_remember_upto_here_1 (x
);
3647 /* we know the result will have to be a bignum */
3650 return scm_i_normbig (result
);
3654 SCM result
= scm_i_mkbig ();
3655 mpz_add_ui (SCM_I_BIG_MPZ (result
), SCM_I_BIG_MPZ (x
), inum
);
3656 scm_remember_upto_here_1 (x
);
3657 /* we know the result will have to be a bignum */
3660 return scm_i_normbig (result
);
3663 else if (SCM_BIGP (y
))
3665 SCM result
= scm_i_mkbig ();
3666 int sgn_x
= mpz_sgn (SCM_I_BIG_MPZ (x
));
3667 int sgn_y
= mpz_sgn (SCM_I_BIG_MPZ (y
));
3668 mpz_add (SCM_I_BIG_MPZ (result
),
3671 scm_remember_upto_here_2 (x
, y
);
3672 /* we know the result will have to be a bignum */
3675 return scm_i_normbig (result
);
3677 else if (SCM_REALP (y
))
3679 double result
= mpz_get_d (SCM_I_BIG_MPZ (x
)) + SCM_REAL_VALUE (y
);
3680 scm_remember_upto_here_1 (x
);
3681 return scm_make_real (result
);
3683 else if (SCM_COMPLEXP (y
))
3685 double real_part
= (mpz_get_d (SCM_I_BIG_MPZ (x
))
3686 + SCM_COMPLEX_REAL (y
));
3687 scm_remember_upto_here_1 (x
);
3688 return scm_make_complex (real_part
, SCM_COMPLEX_IMAG (y
));
3690 else if (SCM_FRACTIONP (y
))
3691 return scm_make_ratio (scm_sum (SCM_FRACTION_NUMERATOR (y
),
3692 scm_product (x
, SCM_FRACTION_DENOMINATOR (y
))),
3693 SCM_FRACTION_DENOMINATOR (y
));
3695 SCM_WTA_DISPATCH_2 (g_sum
, x
, y
, SCM_ARGn
, s_sum
);
3697 else if (SCM_REALP (x
))
3700 return scm_make_real (SCM_REAL_VALUE (x
) + SCM_INUM (y
));
3701 else if (SCM_BIGP (y
))
3703 double result
= mpz_get_d (SCM_I_BIG_MPZ (y
)) + SCM_REAL_VALUE (x
);
3704 scm_remember_upto_here_1 (y
);
3705 return scm_make_real (result
);
3707 else if (SCM_REALP (y
))
3708 return scm_make_real (SCM_REAL_VALUE (x
) + SCM_REAL_VALUE (y
));
3709 else if (SCM_COMPLEXP (y
))
3710 return scm_make_complex (SCM_REAL_VALUE (x
) + SCM_COMPLEX_REAL (y
),
3711 SCM_COMPLEX_IMAG (y
));
3712 else if (SCM_FRACTIONP (y
))
3713 return scm_make_real (SCM_REAL_VALUE (x
) + scm_i_fraction2double (y
));
3715 SCM_WTA_DISPATCH_2 (g_sum
, x
, y
, SCM_ARGn
, s_sum
);
3717 else if (SCM_COMPLEXP (x
))
3720 return scm_make_complex (SCM_COMPLEX_REAL (x
) + SCM_INUM (y
),
3721 SCM_COMPLEX_IMAG (x
));
3722 else if (SCM_BIGP (y
))
3724 double real_part
= (mpz_get_d (SCM_I_BIG_MPZ (y
))
3725 + SCM_COMPLEX_REAL (x
));
3726 scm_remember_upto_here_1 (y
);
3727 return scm_make_complex (real_part
, SCM_COMPLEX_IMAG (x
));
3729 else if (SCM_REALP (y
))
3730 return scm_make_complex (SCM_COMPLEX_REAL (x
) + SCM_REAL_VALUE (y
),
3731 SCM_COMPLEX_IMAG (x
));
3732 else if (SCM_COMPLEXP (y
))
3733 return scm_make_complex (SCM_COMPLEX_REAL (x
) + SCM_COMPLEX_REAL (y
),
3734 SCM_COMPLEX_IMAG (x
) + SCM_COMPLEX_IMAG (y
));
3735 else if (SCM_FRACTIONP (y
))
3736 return scm_make_complex (SCM_COMPLEX_REAL (x
) + scm_i_fraction2double (y
),
3737 SCM_COMPLEX_IMAG (x
));
3739 SCM_WTA_DISPATCH_2 (g_sum
, x
, y
, SCM_ARGn
, s_sum
);
3741 else if (SCM_FRACTIONP (x
))
3744 return scm_make_ratio (scm_sum (SCM_FRACTION_NUMERATOR (x
),
3745 scm_product (y
, SCM_FRACTION_DENOMINATOR (x
))),
3746 SCM_FRACTION_DENOMINATOR (x
));
3747 else if (SCM_BIGP (y
))
3748 return scm_make_ratio (scm_sum (SCM_FRACTION_NUMERATOR (x
),
3749 scm_product (y
, SCM_FRACTION_DENOMINATOR (x
))),
3750 SCM_FRACTION_DENOMINATOR (x
));
3751 else if (SCM_REALP (y
))
3752 return scm_make_real (SCM_REAL_VALUE (y
) + scm_i_fraction2double (x
));
3753 else if (SCM_COMPLEXP (y
))
3754 return scm_make_complex (SCM_COMPLEX_REAL (y
) + scm_i_fraction2double (x
),
3755 SCM_COMPLEX_IMAG (y
));
3756 else if (SCM_FRACTIONP (y
))
3757 /* a/b + c/d = (ad + bc) / bd */
3758 return scm_make_ratio (scm_sum (scm_product (SCM_FRACTION_NUMERATOR (x
), SCM_FRACTION_DENOMINATOR (y
)),
3759 scm_product (SCM_FRACTION_NUMERATOR (y
), SCM_FRACTION_DENOMINATOR (x
))),
3760 scm_product (SCM_FRACTION_DENOMINATOR (x
), SCM_FRACTION_DENOMINATOR (y
)));
3762 SCM_WTA_DISPATCH_2 (g_sum
, x
, y
, SCM_ARGn
, s_sum
);
3765 SCM_WTA_DISPATCH_2 (g_sum
, x
, y
, SCM_ARG1
, s_sum
);
3769 SCM_GPROC1 (s_difference
, "-", scm_tc7_asubr
, scm_difference
, g_difference
);
3770 /* If called with one argument @var{z1}, -@var{z1} returned. Otherwise
3771 * the sum of all but the first argument are subtracted from the first
3773 #define FUNC_NAME s_difference
3775 scm_difference (SCM x
, SCM y
)
3780 SCM_WTA_DISPATCH_0 (g_difference
, s_difference
);
3784 long xx
= -SCM_INUM (x
);
3785 if (SCM_FIXABLE (xx
))
3786 return SCM_MAKINUM (xx
);
3788 return scm_i_long2big (xx
);
3790 else if (SCM_BIGP (x
))
3791 /* FIXME: do we really need to normalize here? */
3792 return scm_i_normbig (scm_i_clonebig (x
, 0));
3793 else if (SCM_REALP (x
))
3794 return scm_make_real (-SCM_REAL_VALUE (x
));
3795 else if (SCM_COMPLEXP (x
))
3796 return scm_make_complex (-SCM_COMPLEX_REAL (x
),
3797 -SCM_COMPLEX_IMAG (x
));
3798 else if (SCM_FRACTIONP (x
))
3799 return scm_make_ratio (scm_difference (SCM_FRACTION_NUMERATOR (x
), SCM_UNDEFINED
),
3800 SCM_FRACTION_DENOMINATOR (x
));
3802 SCM_WTA_DISPATCH_1 (g_difference
, x
, SCM_ARG1
, s_difference
);
3809 long int xx
= SCM_INUM (x
);
3810 long int yy
= SCM_INUM (y
);
3811 long int z
= xx
- yy
;
3812 if (SCM_FIXABLE (z
))
3813 return SCM_MAKINUM (z
);
3815 return scm_i_long2big (z
);
3817 else if (SCM_BIGP (y
))
3819 /* inum-x - big-y */
3820 long xx
= SCM_INUM (x
);
3823 return scm_i_clonebig (y
, 0);
3826 int sgn_y
= mpz_sgn (SCM_I_BIG_MPZ (y
));
3827 SCM result
= scm_i_mkbig ();
3830 mpz_ui_sub (SCM_I_BIG_MPZ (result
), xx
, SCM_I_BIG_MPZ (y
));
3833 /* x - y == -(y + -x) */
3834 mpz_add_ui (SCM_I_BIG_MPZ (result
), SCM_I_BIG_MPZ (y
), -xx
);
3835 mpz_neg (SCM_I_BIG_MPZ (result
), SCM_I_BIG_MPZ (result
));
3837 scm_remember_upto_here_1 (y
);
3839 if ((xx
< 0 && (sgn_y
> 0)) || ((xx
> 0) && sgn_y
< 0))
3840 /* we know the result will have to be a bignum */
3843 return scm_i_normbig (result
);
3846 else if (SCM_REALP (y
))
3848 long int xx
= SCM_INUM (x
);
3849 return scm_make_real (xx
- SCM_REAL_VALUE (y
));
3851 else if (SCM_COMPLEXP (y
))
3853 long int xx
= SCM_INUM (x
);
3854 return scm_make_complex (xx
- SCM_COMPLEX_REAL (y
),
3855 - SCM_COMPLEX_IMAG (y
));
3857 else if (SCM_FRACTIONP (y
))
3858 /* a - b/c = (ac - b) / c */
3859 return scm_make_ratio (scm_difference (scm_product (x
, SCM_FRACTION_DENOMINATOR (y
)),
3860 SCM_FRACTION_NUMERATOR (y
)),
3861 SCM_FRACTION_DENOMINATOR (y
));
3863 SCM_WTA_DISPATCH_2 (g_difference
, x
, y
, SCM_ARGn
, s_difference
);
3865 else if (SCM_BIGP (x
))
3869 /* big-x - inum-y */
3870 long yy
= SCM_INUM (y
);
3871 int sgn_x
= mpz_sgn (SCM_I_BIG_MPZ (x
));
3873 scm_remember_upto_here_1 (x
);
3875 return SCM_FIXABLE (-yy
) ? SCM_MAKINUM (-yy
) : scm_long2num (-yy
);
3878 SCM result
= scm_i_mkbig ();
3881 mpz_sub_ui (SCM_I_BIG_MPZ (result
), SCM_I_BIG_MPZ (x
), yy
);
3883 mpz_add_ui (SCM_I_BIG_MPZ (result
), SCM_I_BIG_MPZ (x
), -yy
);
3884 scm_remember_upto_here_1 (x
);
3886 if ((sgn_x
< 0 && (yy
> 0)) || ((sgn_x
> 0) && yy
< 0))
3887 /* we know the result will have to be a bignum */
3890 return scm_i_normbig (result
);
3893 else if (SCM_BIGP (y
))
3895 int sgn_x
= mpz_sgn (SCM_I_BIG_MPZ (x
));
3896 int sgn_y
= mpz_sgn (SCM_I_BIG_MPZ (y
));
3897 SCM result
= scm_i_mkbig ();
3898 mpz_sub (SCM_I_BIG_MPZ (result
),
3901 scm_remember_upto_here_2 (x
, y
);
3902 /* we know the result will have to be a bignum */
3903 if ((sgn_x
== 1) && (sgn_y
== -1))
3905 if ((sgn_x
== -1) && (sgn_y
== 1))
3907 return scm_i_normbig (result
);
3909 else if (SCM_REALP (y
))
3911 double result
= mpz_get_d (SCM_I_BIG_MPZ (x
)) - SCM_REAL_VALUE (y
);
3912 scm_remember_upto_here_1 (x
);
3913 return scm_make_real (result
);
3915 else if (SCM_COMPLEXP (y
))
3917 double real_part
= (mpz_get_d (SCM_I_BIG_MPZ (x
))
3918 - SCM_COMPLEX_REAL (y
));
3919 scm_remember_upto_here_1 (x
);
3920 return scm_make_complex (real_part
, - SCM_COMPLEX_IMAG (y
));
3922 else if (SCM_FRACTIONP (y
))
3923 return scm_make_ratio (scm_difference (scm_product (x
, SCM_FRACTION_DENOMINATOR (y
)),
3924 SCM_FRACTION_NUMERATOR (y
)),
3925 SCM_FRACTION_DENOMINATOR (y
));
3926 else SCM_WTA_DISPATCH_2 (g_difference
, x
, y
, SCM_ARGn
, s_difference
);
3928 else if (SCM_REALP (x
))
3931 return scm_make_real (SCM_REAL_VALUE (x
) - SCM_INUM (y
));
3932 else if (SCM_BIGP (y
))
3934 double result
= SCM_REAL_VALUE (x
) - mpz_get_d (SCM_I_BIG_MPZ (y
));
3935 scm_remember_upto_here_1 (x
);
3936 return scm_make_real (result
);
3938 else if (SCM_REALP (y
))
3939 return scm_make_real (SCM_REAL_VALUE (x
) - SCM_REAL_VALUE (y
));
3940 else if (SCM_COMPLEXP (y
))
3941 return scm_make_complex (SCM_REAL_VALUE (x
) - SCM_COMPLEX_REAL (y
),
3942 -SCM_COMPLEX_IMAG (y
));
3943 else if (SCM_FRACTIONP (y
))
3944 return scm_make_real (SCM_REAL_VALUE (x
) - scm_i_fraction2double (y
));
3946 SCM_WTA_DISPATCH_2 (g_difference
, x
, y
, SCM_ARGn
, s_difference
);
3948 else if (SCM_COMPLEXP (x
))
3951 return scm_make_complex (SCM_COMPLEX_REAL (x
) - SCM_INUM (y
),
3952 SCM_COMPLEX_IMAG (x
));
3953 else if (SCM_BIGP (y
))
3955 double real_part
= (SCM_COMPLEX_REAL (x
)
3956 - mpz_get_d (SCM_I_BIG_MPZ (y
)));
3957 scm_remember_upto_here_1 (x
);
3958 return scm_make_complex (real_part
, SCM_COMPLEX_IMAG (y
));
3960 else if (SCM_REALP (y
))
3961 return scm_make_complex (SCM_COMPLEX_REAL (x
) - SCM_REAL_VALUE (y
),
3962 SCM_COMPLEX_IMAG (x
));
3963 else if (SCM_COMPLEXP (y
))
3964 return scm_make_complex (SCM_COMPLEX_REAL (x
) - SCM_COMPLEX_REAL (y
),
3965 SCM_COMPLEX_IMAG (x
) - SCM_COMPLEX_IMAG (y
));
3966 else if (SCM_FRACTIONP (y
))
3967 return scm_make_complex (SCM_COMPLEX_REAL (x
) - scm_i_fraction2double (y
),
3968 SCM_COMPLEX_IMAG (x
));
3970 SCM_WTA_DISPATCH_2 (g_difference
, x
, y
, SCM_ARGn
, s_difference
);
3972 else if (SCM_FRACTIONP (x
))
3975 /* a/b - c = (a - cb) / b */
3976 return scm_make_ratio (scm_difference (SCM_FRACTION_NUMERATOR (x
),
3977 scm_product(y
, SCM_FRACTION_DENOMINATOR (x
))),
3978 SCM_FRACTION_DENOMINATOR (x
));
3979 else if (SCM_BIGP (y
))
3980 return scm_make_ratio (scm_difference (SCM_FRACTION_NUMERATOR (x
),
3981 scm_product(y
, SCM_FRACTION_DENOMINATOR (x
))),
3982 SCM_FRACTION_DENOMINATOR (x
));
3983 else if (SCM_REALP (y
))
3984 return scm_make_real (scm_i_fraction2double (x
) - SCM_REAL_VALUE (y
));
3985 else if (SCM_COMPLEXP (y
))
3986 return scm_make_complex (scm_i_fraction2double (x
) - SCM_COMPLEX_REAL (y
),
3987 -SCM_COMPLEX_IMAG (y
));
3988 else if (SCM_FRACTIONP (y
))
3989 /* a/b - c/d = (ad - bc) / bd */
3990 return scm_make_ratio (scm_difference (scm_product (SCM_FRACTION_NUMERATOR (x
), SCM_FRACTION_DENOMINATOR (y
)),
3991 scm_product (SCM_FRACTION_NUMERATOR (y
), SCM_FRACTION_DENOMINATOR (x
))),
3992 scm_product (SCM_FRACTION_DENOMINATOR (x
), SCM_FRACTION_DENOMINATOR (y
)));
3994 SCM_WTA_DISPATCH_2 (g_difference
, x
, y
, SCM_ARGn
, s_difference
);
3997 SCM_WTA_DISPATCH_2 (g_difference
, x
, y
, SCM_ARG1
, s_difference
);
4002 SCM_GPROC1 (s_product
, "*", scm_tc7_asubr
, scm_product
, g_product
);
4003 /* "Return the product of all arguments. If called without arguments,\n"
4007 scm_product (SCM x
, SCM y
)
4012 return SCM_MAKINUM (1L);
4013 else if (SCM_NUMBERP (x
))
4016 SCM_WTA_DISPATCH_1 (g_product
, x
, SCM_ARG1
, s_product
);
4028 case 0: return x
; break;
4029 case 1: return y
; break;
4034 long yy
= SCM_INUM (y
);
4036 SCM k
= SCM_MAKINUM (kk
);
4037 if ((kk
== SCM_INUM (k
)) && (kk
/ xx
== yy
))
4041 SCM result
= scm_i_long2big (xx
);
4042 mpz_mul_si (SCM_I_BIG_MPZ (result
), SCM_I_BIG_MPZ (result
), yy
);
4043 return scm_i_normbig (result
);
4046 else if (SCM_BIGP (y
))
4048 SCM result
= scm_i_mkbig ();
4049 mpz_mul_si (SCM_I_BIG_MPZ (result
), SCM_I_BIG_MPZ (y
), xx
);
4050 scm_remember_upto_here_1 (y
);
4053 else if (SCM_REALP (y
))
4054 return scm_make_real (xx
* SCM_REAL_VALUE (y
));
4055 else if (SCM_COMPLEXP (y
))
4056 return scm_make_complex (xx
* SCM_COMPLEX_REAL (y
),
4057 xx
* SCM_COMPLEX_IMAG (y
));
4058 else if (SCM_FRACTIONP (y
))
4059 return scm_make_ratio (scm_product (x
, SCM_FRACTION_NUMERATOR (y
)),
4060 SCM_FRACTION_DENOMINATOR (y
));
4062 SCM_WTA_DISPATCH_2 (g_product
, x
, y
, SCM_ARGn
, s_product
);
4064 else if (SCM_BIGP (x
))
4071 else if (SCM_BIGP (y
))
4073 SCM result
= scm_i_mkbig ();
4074 mpz_mul (SCM_I_BIG_MPZ (result
),
4077 scm_remember_upto_here_2 (x
, y
);
4080 else if (SCM_REALP (y
))
4082 double result
= mpz_get_d (SCM_I_BIG_MPZ (x
)) * SCM_REAL_VALUE (y
);
4083 scm_remember_upto_here_1 (x
);
4084 return scm_make_real (result
);
4086 else if (SCM_COMPLEXP (y
))
4088 double z
= mpz_get_d (SCM_I_BIG_MPZ (x
));
4089 scm_remember_upto_here_1 (x
);
4090 return scm_make_complex (z
* SCM_COMPLEX_REAL (y
),
4091 z
* SCM_COMPLEX_IMAG (y
));
4093 else if (SCM_FRACTIONP (y
))
4094 return scm_make_ratio (scm_product (x
, SCM_FRACTION_NUMERATOR (y
)),
4095 SCM_FRACTION_DENOMINATOR (y
));
4097 SCM_WTA_DISPATCH_2 (g_product
, x
, y
, SCM_ARGn
, s_product
);
4099 else if (SCM_REALP (x
))
4102 return scm_make_real (SCM_INUM (y
) * SCM_REAL_VALUE (x
));
4103 else if (SCM_BIGP (y
))
4105 double result
= mpz_get_d (SCM_I_BIG_MPZ (y
)) * SCM_REAL_VALUE (x
);
4106 scm_remember_upto_here_1 (y
);
4107 return scm_make_real (result
);
4109 else if (SCM_REALP (y
))
4110 return scm_make_real (SCM_REAL_VALUE (x
) * SCM_REAL_VALUE (y
));
4111 else if (SCM_COMPLEXP (y
))
4112 return scm_make_complex (SCM_REAL_VALUE (x
) * SCM_COMPLEX_REAL (y
),
4113 SCM_REAL_VALUE (x
) * SCM_COMPLEX_IMAG (y
));
4114 else if (SCM_FRACTIONP (y
))
4115 return scm_make_real (SCM_REAL_VALUE (x
) * scm_i_fraction2double (y
));
4117 SCM_WTA_DISPATCH_2 (g_product
, x
, y
, SCM_ARGn
, s_product
);
4119 else if (SCM_COMPLEXP (x
))
4122 return scm_make_complex (SCM_INUM (y
) * SCM_COMPLEX_REAL (x
),
4123 SCM_INUM (y
) * SCM_COMPLEX_IMAG (x
));
4124 else if (SCM_BIGP (y
))
4126 double z
= mpz_get_d (SCM_I_BIG_MPZ (y
));
4127 scm_remember_upto_here_1 (y
);
4128 return scm_make_complex (z
* SCM_COMPLEX_REAL (x
),
4129 z
* SCM_COMPLEX_IMAG (x
));
4131 else if (SCM_REALP (y
))
4132 return scm_make_complex (SCM_REAL_VALUE (y
) * SCM_COMPLEX_REAL (x
),
4133 SCM_REAL_VALUE (y
) * SCM_COMPLEX_IMAG (x
));
4134 else if (SCM_COMPLEXP (y
))
4136 return scm_make_complex (SCM_COMPLEX_REAL (x
) * SCM_COMPLEX_REAL (y
)
4137 - SCM_COMPLEX_IMAG (x
) * SCM_COMPLEX_IMAG (y
),
4138 SCM_COMPLEX_REAL (x
) * SCM_COMPLEX_IMAG (y
)
4139 + SCM_COMPLEX_IMAG (x
) * SCM_COMPLEX_REAL (y
));
4141 else if (SCM_FRACTIONP (y
))
4143 double yy
= scm_i_fraction2double (y
);
4144 return scm_make_complex (yy
* SCM_COMPLEX_REAL (x
),
4145 yy
* SCM_COMPLEX_IMAG (x
));
4148 SCM_WTA_DISPATCH_2 (g_product
, x
, y
, SCM_ARGn
, s_product
);
4150 else if (SCM_FRACTIONP (x
))
4153 return scm_make_ratio (scm_product (y
, SCM_FRACTION_NUMERATOR (x
)),
4154 SCM_FRACTION_DENOMINATOR (x
));
4155 else if (SCM_BIGP (y
))
4156 return scm_make_ratio (scm_product (y
, SCM_FRACTION_NUMERATOR (x
)),
4157 SCM_FRACTION_DENOMINATOR (x
));
4158 else if (SCM_REALP (y
))
4159 return scm_make_real (scm_i_fraction2double (x
) * SCM_REAL_VALUE (y
));
4160 else if (SCM_COMPLEXP (y
))
4162 double xx
= scm_i_fraction2double (x
);
4163 return scm_make_complex (xx
* SCM_COMPLEX_REAL (y
),
4164 xx
* SCM_COMPLEX_IMAG (y
));
4166 else if (SCM_FRACTIONP (y
))
4167 /* a/b * c/d = ac / bd */
4168 return scm_make_ratio (scm_product (SCM_FRACTION_NUMERATOR (x
),
4169 SCM_FRACTION_NUMERATOR (y
)),
4170 scm_product (SCM_FRACTION_DENOMINATOR (x
),
4171 SCM_FRACTION_DENOMINATOR (y
)));
4173 SCM_WTA_DISPATCH_2 (g_product
, x
, y
, SCM_ARGn
, s_product
);
4176 SCM_WTA_DISPATCH_2 (g_product
, x
, y
, SCM_ARG1
, s_product
);
4180 scm_num2dbl (SCM a
, const char *why
)
4181 #define FUNC_NAME why
4184 return (double) SCM_INUM (a
);
4185 else if (SCM_BIGP (a
))
4187 double result
= mpz_get_d (SCM_I_BIG_MPZ (a
));
4188 scm_remember_upto_here_1 (a
);
4191 else if (SCM_REALP (a
))
4192 return (SCM_REAL_VALUE (a
));
4193 else if (SCM_FRACTIONP (a
))
4194 return scm_i_fraction2double (a
);
4196 SCM_WRONG_TYPE_ARG (SCM_ARGn
, a
);
4200 #if ((defined (HAVE_ISINF) && defined (HAVE_ISNAN)) \
4201 || (defined (HAVE_FINITE) && defined (HAVE_ISNAN)))
4202 #define ALLOW_DIVIDE_BY_ZERO
4203 /* #define ALLOW_DIVIDE_BY_EXACT_ZERO */
4206 /* The code below for complex division is adapted from the GNU
4207 libstdc++, which adapted it from f2c's libF77, and is subject to
4210 /****************************************************************
4211 Copyright 1990, 1991, 1992, 1993 by AT&T Bell Laboratories and Bellcore.
4213 Permission to use, copy, modify, and distribute this software
4214 and its documentation for any purpose and without fee is hereby
4215 granted, provided that the above copyright notice appear in all
4216 copies and that both that the copyright notice and this
4217 permission notice and warranty disclaimer appear in supporting
4218 documentation, and that the names of AT&T Bell Laboratories or
4219 Bellcore or any of their entities not be used in advertising or
4220 publicity pertaining to distribution of the software without
4221 specific, written prior permission.
4223 AT&T and Bellcore disclaim all warranties with regard to this
4224 software, including all implied warranties of merchantability
4225 and fitness. In no event shall AT&T or Bellcore be liable for
4226 any special, indirect or consequential damages or any damages
4227 whatsoever resulting from loss of use, data or profits, whether
4228 in an action of contract, negligence or other tortious action,
4229 arising out of or in connection with the use or performance of
4231 ****************************************************************/
4233 SCM_GPROC1 (s_divide
, "/", scm_tc7_asubr
, scm_divide
, g_divide
);
4234 /* Divide the first argument by the product of the remaining
4235 arguments. If called with one argument @var{z1}, 1/@var{z1} is
4237 #define FUNC_NAME s_divide
4239 scm_i_divide (SCM x
, SCM y
, int inexact
)
4246 SCM_WTA_DISPATCH_0 (g_divide
, s_divide
);
4247 else if (SCM_INUMP (x
))
4249 long xx
= SCM_INUM (x
);
4250 if (xx
== 1 || xx
== -1)
4252 #ifndef ALLOW_DIVIDE_BY_EXACT_ZERO
4254 scm_num_overflow (s_divide
);
4259 return scm_make_real (1.0 / (double) xx
);
4260 else return scm_make_ratio (SCM_MAKINUM(1), x
);
4263 else if (SCM_BIGP (x
))
4266 return scm_make_real (1.0 / scm_i_big2dbl (x
));
4267 else return scm_make_ratio (SCM_MAKINUM(1), x
);
4269 else if (SCM_REALP (x
))
4271 double xx
= SCM_REAL_VALUE (x
);
4272 #ifndef ALLOW_DIVIDE_BY_ZERO
4274 scm_num_overflow (s_divide
);
4277 return scm_make_real (1.0 / xx
);
4279 else if (SCM_COMPLEXP (x
))
4281 double r
= SCM_COMPLEX_REAL (x
);
4282 double i
= SCM_COMPLEX_IMAG (x
);
4286 double d
= i
* (1.0 + t
* t
);
4287 return scm_make_complex (t
/ d
, -1.0 / d
);
4292 double d
= r
* (1.0 + t
* t
);
4293 return scm_make_complex (1.0 / d
, -t
/ d
);
4296 else if (SCM_FRACTIONP (x
))
4297 return scm_make_ratio (SCM_FRACTION_DENOMINATOR (x
),
4298 SCM_FRACTION_NUMERATOR (x
));
4300 SCM_WTA_DISPATCH_1 (g_divide
, x
, SCM_ARG1
, s_divide
);
4305 long xx
= SCM_INUM (x
);
4308 long yy
= SCM_INUM (y
);
4311 #ifndef ALLOW_DIVIDE_BY_EXACT_ZERO
4312 scm_num_overflow (s_divide
);
4314 return scm_make_real ((double) xx
/ (double) yy
);
4317 else if (xx
% yy
!= 0)
4320 return scm_make_real ((double) xx
/ (double) yy
);
4321 else return scm_make_ratio (x
, y
);
4326 if (SCM_FIXABLE (z
))
4327 return SCM_MAKINUM (z
);
4329 return scm_i_long2big (z
);
4332 else if (SCM_BIGP (y
))
4335 return scm_make_real ((double) xx
/ scm_i_big2dbl (y
));
4336 else return scm_make_ratio (x
, y
);
4338 else if (SCM_REALP (y
))
4340 double yy
= SCM_REAL_VALUE (y
);
4341 #ifndef ALLOW_DIVIDE_BY_ZERO
4343 scm_num_overflow (s_divide
);
4346 return scm_make_real ((double) xx
/ yy
);
4348 else if (SCM_COMPLEXP (y
))
4351 complex_div
: /* y _must_ be a complex number */
4353 double r
= SCM_COMPLEX_REAL (y
);
4354 double i
= SCM_COMPLEX_IMAG (y
);
4358 double d
= i
* (1.0 + t
* t
);
4359 return scm_make_complex ((a
* t
) / d
, -a
/ d
);
4364 double d
= r
* (1.0 + t
* t
);
4365 return scm_make_complex (a
/ d
, -(a
* t
) / d
);
4369 else if (SCM_FRACTIONP (y
))
4370 /* a / b/c = ac / b */
4371 return scm_make_ratio (scm_product (x
, SCM_FRACTION_DENOMINATOR (y
)),
4372 SCM_FRACTION_NUMERATOR (y
));
4374 SCM_WTA_DISPATCH_2 (g_divide
, x
, y
, SCM_ARGn
, s_divide
);
4376 else if (SCM_BIGP (x
))
4380 long int yy
= SCM_INUM (y
);
4383 #ifndef ALLOW_DIVIDE_BY_EXACT_ZERO
4384 scm_num_overflow (s_divide
);
4386 int sgn
= mpz_sgn (SCM_I_BIG_MPZ (x
));
4387 scm_remember_upto_here_1 (x
);
4388 return (sgn
== 0) ? scm_nan () : scm_inf ();
4395 /* FIXME: HMM, what are the relative performance issues here?
4396 We need to test. Is it faster on average to test
4397 divisible_p, then perform whichever operation, or is it
4398 faster to perform the integer div opportunistically and
4399 switch to real if there's a remainder? For now we take the
4400 middle ground: test, then if divisible, use the faster div
4403 long abs_yy
= yy
< 0 ? -yy
: yy
;
4404 int divisible_p
= mpz_divisible_ui_p (SCM_I_BIG_MPZ (x
), abs_yy
);
4408 SCM result
= scm_i_mkbig ();
4409 mpz_divexact_ui (SCM_I_BIG_MPZ (result
), SCM_I_BIG_MPZ (x
), abs_yy
);
4410 scm_remember_upto_here_1 (x
);
4412 mpz_neg (SCM_I_BIG_MPZ (result
), SCM_I_BIG_MPZ (result
));
4413 return scm_i_normbig (result
);
4418 return scm_make_real (scm_i_big2dbl (x
) / (double) yy
);
4419 else return scm_make_ratio (x
, y
);
4423 else if (SCM_BIGP (y
))
4425 int y_is_zero
= (mpz_sgn (SCM_I_BIG_MPZ (y
)) == 0);
4428 #ifndef ALLOW_DIVIDE_BY_EXACT_ZERO
4429 scm_num_overflow (s_divide
);
4431 int sgn
= mpz_sgn (SCM_I_BIG_MPZ (x
));
4432 scm_remember_upto_here_1 (x
);
4433 return (sgn
== 0) ? scm_nan () : scm_inf ();
4439 int divisible_p
= mpz_divisible_p (SCM_I_BIG_MPZ (x
),
4443 SCM result
= scm_i_mkbig ();
4444 mpz_divexact (SCM_I_BIG_MPZ (result
),
4447 scm_remember_upto_here_2 (x
, y
);
4448 return scm_i_normbig (result
);
4454 double dbx
= mpz_get_d (SCM_I_BIG_MPZ (x
));
4455 double dby
= mpz_get_d (SCM_I_BIG_MPZ (y
));
4456 scm_remember_upto_here_2 (x
, y
);
4457 return scm_make_real (dbx
/ dby
);
4459 else return scm_make_ratio (x
, y
);
4463 else if (SCM_REALP (y
))
4465 double yy
= SCM_REAL_VALUE (y
);
4466 #ifndef ALLOW_DIVIDE_BY_ZERO
4468 scm_num_overflow (s_divide
);
4471 return scm_make_real (scm_i_big2dbl (x
) / yy
);
4473 else if (SCM_COMPLEXP (y
))
4475 a
= scm_i_big2dbl (x
);
4478 else if (SCM_FRACTIONP (y
))
4479 return scm_make_ratio (scm_product (x
, SCM_FRACTION_DENOMINATOR (y
)),
4480 SCM_FRACTION_NUMERATOR (y
));
4482 SCM_WTA_DISPATCH_2 (g_divide
, x
, y
, SCM_ARGn
, s_divide
);
4484 else if (SCM_REALP (x
))
4486 double rx
= SCM_REAL_VALUE (x
);
4489 long int yy
= SCM_INUM (y
);
4490 #ifndef ALLOW_DIVIDE_BY_EXACT_ZERO
4492 scm_num_overflow (s_divide
);
4495 return scm_make_real (rx
/ (double) yy
);
4497 else if (SCM_BIGP (y
))
4499 double dby
= mpz_get_d (SCM_I_BIG_MPZ (y
));
4500 scm_remember_upto_here_1 (y
);
4501 return scm_make_real (rx
/ dby
);
4503 else if (SCM_REALP (y
))
4505 double yy
= SCM_REAL_VALUE (y
);
4506 #ifndef ALLOW_DIVIDE_BY_ZERO
4508 scm_num_overflow (s_divide
);
4511 return scm_make_real (rx
/ yy
);
4513 else if (SCM_COMPLEXP (y
))
4518 else if (SCM_FRACTIONP (y
))
4519 return scm_make_real (rx
/ scm_i_fraction2double (y
));
4521 SCM_WTA_DISPATCH_2 (g_divide
, x
, y
, SCM_ARGn
, s_divide
);
4523 else if (SCM_COMPLEXP (x
))
4525 double rx
= SCM_COMPLEX_REAL (x
);
4526 double ix
= SCM_COMPLEX_IMAG (x
);
4529 long int yy
= SCM_INUM (y
);
4530 #ifndef ALLOW_DIVIDE_BY_EXACT_ZERO
4532 scm_num_overflow (s_divide
);
4537 return scm_make_complex (rx
/ d
, ix
/ d
);
4540 else if (SCM_BIGP (y
))
4542 double dby
= mpz_get_d (SCM_I_BIG_MPZ (y
));
4543 scm_remember_upto_here_1 (y
);
4544 return scm_make_complex (rx
/ dby
, ix
/ dby
);
4546 else if (SCM_REALP (y
))
4548 double yy
= SCM_REAL_VALUE (y
);
4549 #ifndef ALLOW_DIVIDE_BY_ZERO
4551 scm_num_overflow (s_divide
);
4554 return scm_make_complex (rx
/ yy
, ix
/ yy
);
4556 else if (SCM_COMPLEXP (y
))
4558 double ry
= SCM_COMPLEX_REAL (y
);
4559 double iy
= SCM_COMPLEX_IMAG (y
);
4563 double d
= iy
* (1.0 + t
* t
);
4564 return scm_make_complex ((rx
* t
+ ix
) / d
, (ix
* t
- rx
) / d
);
4569 double d
= ry
* (1.0 + t
* t
);
4570 return scm_make_complex ((rx
+ ix
* t
) / d
, (ix
- rx
* t
) / d
);
4573 else if (SCM_FRACTIONP (y
))
4575 double yy
= scm_i_fraction2double (y
);
4576 return scm_make_complex (rx
/ yy
, ix
/ yy
);
4579 SCM_WTA_DISPATCH_2 (g_divide
, x
, y
, SCM_ARGn
, s_divide
);
4581 else if (SCM_FRACTIONP (x
))
4585 long int yy
= SCM_INUM (y
);
4586 #ifndef ALLOW_DIVIDE_BY_EXACT_ZERO
4588 scm_num_overflow (s_divide
);
4591 return scm_make_ratio (SCM_FRACTION_NUMERATOR (x
),
4592 scm_product (SCM_FRACTION_DENOMINATOR (x
), y
));
4594 else if (SCM_BIGP (y
))
4596 return scm_make_ratio (SCM_FRACTION_NUMERATOR (x
),
4597 scm_product (SCM_FRACTION_DENOMINATOR (x
), y
));
4599 else if (SCM_REALP (y
))
4601 double yy
= SCM_REAL_VALUE (y
);
4602 #ifndef ALLOW_DIVIDE_BY_ZERO
4604 scm_num_overflow (s_divide
);
4607 return scm_make_real (scm_i_fraction2double (x
) / yy
);
4609 else if (SCM_COMPLEXP (y
))
4611 a
= scm_i_fraction2double (x
);
4614 else if (SCM_FRACTIONP (y
))
4615 return scm_make_ratio (scm_product (SCM_FRACTION_NUMERATOR (x
), SCM_FRACTION_DENOMINATOR (y
)),
4616 scm_product (SCM_FRACTION_NUMERATOR (y
), SCM_FRACTION_DENOMINATOR (x
)));
4618 SCM_WTA_DISPATCH_2 (g_divide
, x
, y
, SCM_ARGn
, s_divide
);
4621 SCM_WTA_DISPATCH_2 (g_divide
, x
, y
, SCM_ARG1
, s_divide
);
4625 scm_divide (SCM x
, SCM y
)
4627 return scm_i_divide (x
, y
, 0);
4630 static SCM
scm_divide2real (SCM x
, SCM y
)
4632 return scm_i_divide (x
, y
, 1);
4638 scm_asinh (double x
)
4643 #define asinh scm_asinh
4644 return log (x
+ sqrt (x
* x
+ 1));
4647 SCM_GPROC1 (s_asinh
, "$asinh", scm_tc7_dsubr
, (SCM (*)()) asinh
, g_asinh
);
4648 /* "Return the inverse hyperbolic sine of @var{x}."
4653 scm_acosh (double x
)
4658 #define acosh scm_acosh
4659 return log (x
+ sqrt (x
* x
- 1));
4662 SCM_GPROC1 (s_acosh
, "$acosh", scm_tc7_dsubr
, (SCM (*)()) acosh
, g_acosh
);
4663 /* "Return the inverse hyperbolic cosine of @var{x}."
4668 scm_atanh (double x
)
4673 #define atanh scm_atanh
4674 return 0.5 * log ((1 + x
) / (1 - x
));
4677 SCM_GPROC1 (s_atanh
, "$atanh", scm_tc7_dsubr
, (SCM (*)()) atanh
, g_atanh
);
4678 /* "Return the inverse hyperbolic tangent of @var{x}."
4682 /* XXX - eventually, we should remove this definition of scm_round and
4683 rename scm_round_number to scm_round. Likewise for scm_truncate
4684 and scm_truncate_number.
4688 scm_truncate (double x
)
4693 #define trunc scm_truncate
4701 scm_round (double x
)
4703 double plus_half
= x
+ 0.5;
4704 double result
= floor (plus_half
);
4705 /* Adjust so that the scm_round is towards even. */
4706 return ((plus_half
== result
&& plus_half
/ 2 != floor (plus_half
/ 2))
4711 SCM_DEFINE (scm_truncate_number
, "truncate", 1, 0, 0,
4713 "Round the number @var{x} towards zero.")
4714 #define FUNC_NAME s_scm_truncate_number
4716 if (SCM_FALSEP (scm_negative_p (x
)))
4717 return scm_floor (x
);
4719 return scm_ceiling (x
);
4723 static SCM exactly_one_half
;
4725 SCM_DEFINE (scm_round_number
, "round", 1, 0, 0,
4727 "Round the number @var{x} towards the nearest integer. "
4728 "When it is exactly halfway between two integers, "
4729 "round towards the even one.")
4730 #define FUNC_NAME s_scm_round_number
4732 SCM plus_half
= scm_sum (x
, exactly_one_half
);
4733 SCM result
= scm_floor (plus_half
);
4734 /* Adjust so that the scm_round is towards even. */
4735 if (!SCM_FALSEP (scm_num_eq_p (plus_half
, result
))
4736 && !SCM_FALSEP (scm_odd_p (result
)))
4737 return scm_difference (result
, SCM_MAKINUM (1));
4743 SCM_PRIMITIVE_GENERIC (scm_floor
, "floor", 1, 0, 0,
4745 "Round the number @var{x} towards minus infinity.")
4746 #define FUNC_NAME s_scm_floor
4748 if (SCM_INUMP (x
) || SCM_BIGP (x
))
4750 else if (SCM_REALP (x
))
4751 return scm_make_real (floor (SCM_REAL_VALUE (x
)));
4752 else if (SCM_FRACTIONP (x
))
4754 SCM q
= scm_quotient (SCM_FRACTION_NUMERATOR (x
),
4755 SCM_FRACTION_DENOMINATOR (x
));
4756 if (SCM_FALSEP (scm_negative_p (x
)))
4758 /* For positive x, rounding towards zero is correct. */
4763 /* For negative x, we need to return q-1 unless x is an
4764 integer. But fractions are never integer, per our
4766 return scm_difference (q
, SCM_MAKINUM (1));
4770 SCM_WTA_DISPATCH_1 (g_scm_floor
, x
, 1, s_scm_floor
);
4774 SCM_PRIMITIVE_GENERIC (scm_ceiling
, "ceiling", 1, 0, 0,
4776 "Round the number @var{x} towards infinity.")
4777 #define FUNC_NAME s_scm_ceiling
4779 if (SCM_INUMP (x
) || SCM_BIGP (x
))
4781 else if (SCM_REALP (x
))
4782 return scm_make_real (ceil (SCM_REAL_VALUE (x
)));
4783 else if (SCM_FRACTIONP (x
))
4785 SCM q
= scm_quotient (SCM_FRACTION_NUMERATOR (x
),
4786 SCM_FRACTION_DENOMINATOR (x
));
4787 if (SCM_FALSEP (scm_positive_p (x
)))
4789 /* For negative x, rounding towards zero is correct. */
4794 /* For positive x, we need to return q+1 unless x is an
4795 integer. But fractions are never integer, per our
4797 return scm_sum (q
, SCM_MAKINUM (1));
4801 SCM_WTA_DISPATCH_1 (g_scm_ceiling
, x
, 1, s_scm_ceiling
);
4805 SCM_GPROC1 (s_i_sqrt
, "$sqrt", scm_tc7_dsubr
, (SCM (*)()) sqrt
, g_i_sqrt
);
4806 /* "Return the square root of the real number @var{x}."
4808 SCM_GPROC1 (s_i_abs
, "$abs", scm_tc7_dsubr
, (SCM (*)()) fabs
, g_i_abs
);
4809 /* "Return the absolute value of the real number @var{x}."
4811 SCM_GPROC1 (s_i_exp
, "$exp", scm_tc7_dsubr
, (SCM (*)()) exp
, g_i_exp
);
4812 /* "Return the @var{x}th power of e."
4814 SCM_GPROC1 (s_i_log
, "$log", scm_tc7_dsubr
, (SCM (*)()) log
, g_i_log
);
4815 /* "Return the natural logarithm of the real number @var{x}."
4817 SCM_GPROC1 (s_i_sin
, "$sin", scm_tc7_dsubr
, (SCM (*)()) sin
, g_i_sin
);
4818 /* "Return the sine of the real number @var{x}."
4820 SCM_GPROC1 (s_i_cos
, "$cos", scm_tc7_dsubr
, (SCM (*)()) cos
, g_i_cos
);
4821 /* "Return the cosine of the real number @var{x}."
4823 SCM_GPROC1 (s_i_tan
, "$tan", scm_tc7_dsubr
, (SCM (*)()) tan
, g_i_tan
);
4824 /* "Return the tangent of the real number @var{x}."
4826 SCM_GPROC1 (s_i_asin
, "$asin", scm_tc7_dsubr
, (SCM (*)()) asin
, g_i_asin
);
4827 /* "Return the arc sine of the real number @var{x}."
4829 SCM_GPROC1 (s_i_acos
, "$acos", scm_tc7_dsubr
, (SCM (*)()) acos
, g_i_acos
);
4830 /* "Return the arc cosine of the real number @var{x}."
4832 SCM_GPROC1 (s_i_atan
, "$atan", scm_tc7_dsubr
, (SCM (*)()) atan
, g_i_atan
);
4833 /* "Return the arc tangent of the real number @var{x}."
4835 SCM_GPROC1 (s_i_sinh
, "$sinh", scm_tc7_dsubr
, (SCM (*)()) sinh
, g_i_sinh
);
4836 /* "Return the hyperbolic sine of the real number @var{x}."
4838 SCM_GPROC1 (s_i_cosh
, "$cosh", scm_tc7_dsubr
, (SCM (*)()) cosh
, g_i_cosh
);
4839 /* "Return the hyperbolic cosine of the real number @var{x}."
4841 SCM_GPROC1 (s_i_tanh
, "$tanh", scm_tc7_dsubr
, (SCM (*)()) tanh
, g_i_tanh
);
4842 /* "Return the hyperbolic tangent of the real number @var{x}."
4850 static void scm_two_doubles (SCM x
,
4852 const char *sstring
,
4856 scm_two_doubles (SCM x
, SCM y
, const char *sstring
, struct dpair
*xy
)
4859 xy
->x
= SCM_INUM (x
);
4860 else if (SCM_BIGP (x
))
4861 xy
->x
= scm_i_big2dbl (x
);
4862 else if (SCM_REALP (x
))
4863 xy
->x
= SCM_REAL_VALUE (x
);
4864 else if (SCM_FRACTIONP (x
))
4865 xy
->x
= scm_i_fraction2double (x
);
4867 scm_wrong_type_arg (sstring
, SCM_ARG1
, x
);
4870 xy
->y
= SCM_INUM (y
);
4871 else if (SCM_BIGP (y
))
4872 xy
->y
= scm_i_big2dbl (y
);
4873 else if (SCM_REALP (y
))
4874 xy
->y
= SCM_REAL_VALUE (y
);
4875 else if (SCM_FRACTIONP (y
))
4876 xy
->y
= scm_i_fraction2double (y
);
4878 scm_wrong_type_arg (sstring
, SCM_ARG2
, y
);
4882 SCM_DEFINE (scm_sys_expt
, "$expt", 2, 0, 0,
4884 "Return @var{x} raised to the power of @var{y}. This\n"
4885 "procedure does not accept complex arguments.")
4886 #define FUNC_NAME s_scm_sys_expt
4889 scm_two_doubles (x
, y
, FUNC_NAME
, &xy
);
4890 return scm_make_real (pow (xy
.x
, xy
.y
));
4895 SCM_DEFINE (scm_sys_atan2
, "$atan2", 2, 0, 0,
4897 "Return the arc tangent of the two arguments @var{x} and\n"
4898 "@var{y}. This is similar to calculating the arc tangent of\n"
4899 "@var{x} / @var{y}, except that the signs of both arguments\n"
4900 "are used to determine the quadrant of the result. This\n"
4901 "procedure does not accept complex arguments.")
4902 #define FUNC_NAME s_scm_sys_atan2
4905 scm_two_doubles (x
, y
, FUNC_NAME
, &xy
);
4906 return scm_make_real (atan2 (xy
.x
, xy
.y
));
4911 SCM_DEFINE (scm_make_rectangular
, "make-rectangular", 2, 0, 0,
4912 (SCM real
, SCM imaginary
),
4913 "Return a complex number constructed of the given @var{real} and\n"
4914 "@var{imaginary} parts.")
4915 #define FUNC_NAME s_scm_make_rectangular
4918 scm_two_doubles (real
, imaginary
, FUNC_NAME
, &xy
);
4919 return scm_make_complex (xy
.x
, xy
.y
);
4925 SCM_DEFINE (scm_make_polar
, "make-polar", 2, 0, 0,
4927 "Return the complex number @var{x} * e^(i * @var{y}).")
4928 #define FUNC_NAME s_scm_make_polar
4932 scm_two_doubles (x
, y
, FUNC_NAME
, &xy
);
4934 sincos (xy
.y
, &s
, &c
);
4939 return scm_make_complex (xy
.x
* c
, xy
.x
* s
);
4944 SCM_GPROC (s_real_part
, "real-part", 1, 0, 0, scm_real_part
, g_real_part
);
4945 /* "Return the real part of the number @var{z}."
4948 scm_real_part (SCM z
)
4952 else if (SCM_BIGP (z
))
4954 else if (SCM_REALP (z
))
4956 else if (SCM_COMPLEXP (z
))
4957 return scm_make_real (SCM_COMPLEX_REAL (z
));
4958 else if (SCM_FRACTIONP (z
))
4959 return scm_make_real (scm_i_fraction2double (z
));
4961 SCM_WTA_DISPATCH_1 (g_real_part
, z
, SCM_ARG1
, s_real_part
);
4965 SCM_GPROC (s_imag_part
, "imag-part", 1, 0, 0, scm_imag_part
, g_imag_part
);
4966 /* "Return the imaginary part of the number @var{z}."
4969 scm_imag_part (SCM z
)
4973 else if (SCM_BIGP (z
))
4975 else if (SCM_REALP (z
))
4977 else if (SCM_COMPLEXP (z
))
4978 return scm_make_real (SCM_COMPLEX_IMAG (z
));
4979 else if (SCM_FRACTIONP (z
))
4982 SCM_WTA_DISPATCH_1 (g_imag_part
, z
, SCM_ARG1
, s_imag_part
);
4985 SCM_GPROC (s_numerator
, "numerator", 1, 0, 0, scm_numerator
, g_numerator
);
4986 /* "Return the numerator of the number @var{z}."
4989 scm_numerator (SCM z
)
4993 else if (SCM_BIGP (z
))
4995 else if (SCM_FRACTIONP (z
))
4997 scm_i_fraction_reduce (z
);
4998 return SCM_FRACTION_NUMERATOR (z
);
5000 else if (SCM_REALP (z
))
5001 return scm_exact_to_inexact (scm_numerator (scm_inexact_to_exact (z
)));
5003 SCM_WTA_DISPATCH_1 (g_numerator
, z
, SCM_ARG1
, s_numerator
);
5007 SCM_GPROC (s_denominator
, "denominator", 1, 0, 0, scm_denominator
, g_denominator
);
5008 /* "Return the denominator of the number @var{z}."
5011 scm_denominator (SCM z
)
5014 return SCM_MAKINUM (1);
5015 else if (SCM_BIGP (z
))
5016 return SCM_MAKINUM (1);
5017 else if (SCM_FRACTIONP (z
))
5019 scm_i_fraction_reduce (z
);
5020 return SCM_FRACTION_DENOMINATOR (z
);
5022 else if (SCM_REALP (z
))
5023 return scm_exact_to_inexact (scm_denominator (scm_inexact_to_exact (z
)));
5025 SCM_WTA_DISPATCH_1 (g_denominator
, z
, SCM_ARG1
, s_denominator
);
5028 SCM_GPROC (s_magnitude
, "magnitude", 1, 0, 0, scm_magnitude
, g_magnitude
);
5029 /* "Return the magnitude of the number @var{z}. This is the same as\n"
5030 * "@code{abs} for real arguments, but also allows complex numbers."
5033 scm_magnitude (SCM z
)
5037 long int zz
= SCM_INUM (z
);
5040 else if (SCM_POSFIXABLE (-zz
))
5041 return SCM_MAKINUM (-zz
);
5043 return scm_i_long2big (-zz
);
5045 else if (SCM_BIGP (z
))
5047 int sgn
= mpz_sgn (SCM_I_BIG_MPZ (z
));
5048 scm_remember_upto_here_1 (z
);
5050 return scm_i_clonebig (z
, 0);
5054 else if (SCM_REALP (z
))
5055 return scm_make_real (fabs (SCM_REAL_VALUE (z
)));
5056 else if (SCM_COMPLEXP (z
))
5057 return scm_make_real (hypot (SCM_COMPLEX_REAL (z
), SCM_COMPLEX_IMAG (z
)));
5058 else if (SCM_FRACTIONP (z
))
5060 if (SCM_FALSEP (scm_negative_p (SCM_FRACTION_NUMERATOR (z
))))
5062 return scm_make_ratio (scm_difference (SCM_FRACTION_NUMERATOR (z
), SCM_UNDEFINED
),
5063 SCM_FRACTION_DENOMINATOR (z
));
5066 SCM_WTA_DISPATCH_1 (g_magnitude
, z
, SCM_ARG1
, s_magnitude
);
5070 SCM_GPROC (s_angle
, "angle", 1, 0, 0, scm_angle
, g_angle
);
5071 /* "Return the angle of the complex number @var{z}."
5076 /* atan(0,-1) is pi and it'd be possible to have that as a constant like
5077 scm_flo0 to save allocating a new flonum with scm_make_real each time.
5078 But if atan2 follows the floating point rounding mode, then the value
5079 is not a constant. Maybe it'd be close enough though. */
5082 if (SCM_INUM (z
) >= 0)
5085 return scm_make_real (atan2 (0.0, -1.0));
5087 else if (SCM_BIGP (z
))
5089 int sgn
= mpz_sgn (SCM_I_BIG_MPZ (z
));
5090 scm_remember_upto_here_1 (z
);
5092 return scm_make_real (atan2 (0.0, -1.0));
5096 else if (SCM_REALP (z
))
5098 if (SCM_REAL_VALUE (z
) >= 0)
5101 return scm_make_real (atan2 (0.0, -1.0));
5103 else if (SCM_COMPLEXP (z
))
5104 return scm_make_real (atan2 (SCM_COMPLEX_IMAG (z
), SCM_COMPLEX_REAL (z
)));
5105 else if (SCM_FRACTIONP (z
))
5107 if (SCM_FALSEP (scm_negative_p (SCM_FRACTION_NUMERATOR (z
))))
5109 else return scm_make_real (atan2 (0.0, -1.0));
5112 SCM_WTA_DISPATCH_1 (g_angle
, z
, SCM_ARG1
, s_angle
);
5116 SCM_GPROC (s_exact_to_inexact
, "exact->inexact", 1, 0, 0, scm_exact_to_inexact
, g_exact_to_inexact
);
5117 /* Convert the number @var{x} to its inexact representation.\n"
5120 scm_exact_to_inexact (SCM z
)
5123 return scm_make_real ((double) SCM_INUM (z
));
5124 else if (SCM_BIGP (z
))
5125 return scm_make_real (scm_i_big2dbl (z
));
5126 else if (SCM_FRACTIONP (z
))
5127 return scm_make_real (scm_i_fraction2double (z
));
5128 else if (SCM_INEXACTP (z
))
5131 SCM_WTA_DISPATCH_1 (g_exact_to_inexact
, z
, 1, s_exact_to_inexact
);
5135 SCM_DEFINE (scm_inexact_to_exact
, "inexact->exact", 1, 0, 0,
5137 "Return an exact number that is numerically closest to @var{z}.")
5138 #define FUNC_NAME s_scm_inexact_to_exact
5142 else if (SCM_BIGP (z
))
5144 else if (SCM_REALP (z
))
5146 if (xisinf (SCM_REAL_VALUE (z
)) || xisnan (SCM_REAL_VALUE (z
)))
5147 SCM_OUT_OF_RANGE (1, z
);
5154 mpq_set_d (frac
, SCM_REAL_VALUE (z
));
5155 q
= scm_make_ratio (scm_i_mpz2num (mpq_numref (frac
)),
5156 scm_i_mpz2num (mpq_denref (frac
)));
5158 /* When scm_make_ratio throws, we leak the memory allocated
5165 else if (SCM_FRACTIONP (z
))
5168 SCM_WRONG_TYPE_ARG (1, z
);
5172 SCM_DEFINE (scm_rationalize
, "rationalize", 2, 0, 0,
5174 "Return an exact number that is within @var{err} of @var{x}.")
5175 #define FUNC_NAME s_scm_rationalize
5179 else if (SCM_BIGP (x
))
5181 else if ((SCM_REALP (x
)) || SCM_FRACTIONP (x
))
5183 /* Use continued fractions to find closest ratio. All
5184 arithmetic is done with exact numbers.
5187 SCM ex
= scm_inexact_to_exact (x
);
5188 SCM int_part
= scm_floor (ex
);
5189 SCM tt
= SCM_MAKINUM (1);
5190 SCM a1
= SCM_MAKINUM (0), a2
= SCM_MAKINUM (1), a
= SCM_MAKINUM (0);
5191 SCM b1
= SCM_MAKINUM (1), b2
= SCM_MAKINUM (0), b
= SCM_MAKINUM (0);
5195 if (!SCM_FALSEP (scm_num_eq_p (ex
, int_part
)))
5198 ex
= scm_difference (ex
, int_part
); /* x = x-int_part */
5199 rx
= scm_divide (ex
, SCM_UNDEFINED
); /* rx = 1/x */
5201 /* We stop after a million iterations just to be absolutely sure
5202 that we don't go into an infinite loop. The process normally
5203 converges after less than a dozen iterations.
5206 err
= scm_abs (err
);
5207 while (++i
< 1000000)
5209 a
= scm_sum (scm_product (a1
, tt
), a2
); /* a = a1*tt + a2 */
5210 b
= scm_sum (scm_product (b1
, tt
), b2
); /* b = b1*tt + b2 */
5211 if (SCM_FALSEP (scm_zero_p (b
)) && /* b != 0 */
5213 (scm_gr_p (scm_abs (scm_difference (ex
, scm_divide (a
, b
))),
5214 err
))) /* abs(x-a/b) <= err */
5215 return scm_sum (int_part
, scm_divide (a
, b
)); /* int_part+a/b */
5216 rx
= scm_divide (scm_difference (rx
, tt
), /* rx = 1/(rx - tt) */
5218 tt
= scm_floor (rx
); /* tt = floor (rx) */
5224 scm_num_overflow (s_scm_rationalize
);
5227 SCM_WRONG_TYPE_ARG (1, x
);
5231 /* if you need to change this, change test-num2integral.c as well */
5232 #if SCM_SIZEOF_LONG_LONG != 0
5234 # define ULLONG_MAX ((unsigned long long) (-1))
5235 # define LLONG_MAX ((long long) (ULLONG_MAX >> 1))
5236 # define LLONG_MIN (~LLONG_MAX)
5240 /* Parameters for creating integer conversion routines.
5242 Define the following preprocessor macros before including
5243 "libguile/num2integral.i.c":
5245 NUM2INTEGRAL - the name of the function for converting from a
5246 Scheme object to the integral type. This function will be
5247 defined when including "num2integral.i.c".
5249 INTEGRAL2NUM - the name of the function for converting from the
5250 integral type to a Scheme object. This function will be defined.
5252 INTEGRAL2BIG - the name of an internal function that createas a
5253 bignum from the integral type. This function will be defined.
5254 The name should start with "scm_i_".
5256 ITYPE - the name of the integral type.
5258 UNSIGNED - Define this to 1 when ITYPE is an unsigned type. Define
5261 UNSIGNED_ITYPE - the name of the the unsigned variant of the
5262 integral type. If you don't define this, it defaults to
5263 "unsigned ITYPE" for signed types and simply "ITYPE" for unsigned
5266 SIZEOF_ITYPE - an expression giving the size of the integral type
5267 in bytes. This expression must be computable by the
5268 preprocessor. (SIZEOF_FOO values are calculated by configure.in
5273 #define NUM2INTEGRAL scm_num2short
5274 #define INTEGRAL2NUM scm_short2num
5275 #define INTEGRAL2BIG scm_i_short2big
5278 #define SIZEOF_ITYPE SIZEOF_SHORT
5279 #include "libguile/num2integral.i.c"
5281 #define NUM2INTEGRAL scm_num2ushort
5282 #define INTEGRAL2NUM scm_ushort2num
5283 #define INTEGRAL2BIG scm_i_ushort2big
5285 #define ITYPE unsigned short
5286 #define SIZEOF_ITYPE SIZEOF_UNSIGNED_SHORT
5287 #include "libguile/num2integral.i.c"
5289 #define NUM2INTEGRAL scm_num2int
5290 #define INTEGRAL2NUM scm_int2num
5291 #define INTEGRAL2BIG scm_i_int2big
5294 #define SIZEOF_ITYPE SIZEOF_INT
5295 #include "libguile/num2integral.i.c"
5297 #define NUM2INTEGRAL scm_num2uint
5298 #define INTEGRAL2NUM scm_uint2num
5299 #define INTEGRAL2BIG scm_i_uint2big
5301 #define ITYPE unsigned int
5302 #define SIZEOF_ITYPE SIZEOF_UNSIGNED_INT
5303 #include "libguile/num2integral.i.c"
5305 #define NUM2INTEGRAL scm_num2long
5306 #define INTEGRAL2NUM scm_long2num
5307 #define INTEGRAL2BIG scm_i_long2big
5310 #define SIZEOF_ITYPE SIZEOF_LONG
5311 #include "libguile/num2integral.i.c"
5313 #define NUM2INTEGRAL scm_num2ulong
5314 #define INTEGRAL2NUM scm_ulong2num
5315 #define INTEGRAL2BIG scm_i_ulong2big
5317 #define ITYPE unsigned long
5318 #define SIZEOF_ITYPE SIZEOF_UNSIGNED_LONG
5319 #include "libguile/num2integral.i.c"
5321 #define NUM2INTEGRAL scm_num2ptrdiff
5322 #define INTEGRAL2NUM scm_ptrdiff2num
5323 #define INTEGRAL2BIG scm_i_ptrdiff2big
5325 #define ITYPE scm_t_ptrdiff
5326 #define UNSIGNED_ITYPE size_t
5327 #define SIZEOF_ITYPE SCM_SIZEOF_SCM_T_PTRDIFF
5328 #include "libguile/num2integral.i.c"
5330 #define NUM2INTEGRAL scm_num2size
5331 #define INTEGRAL2NUM scm_size2num
5332 #define INTEGRAL2BIG scm_i_size2big
5334 #define ITYPE size_t
5335 #define SIZEOF_ITYPE SIZEOF_SIZE_T
5336 #include "libguile/num2integral.i.c"
5338 #if SCM_SIZEOF_LONG_LONG != 0
5340 #ifndef ULONG_LONG_MAX
5341 #define ULONG_LONG_MAX (~0ULL)
5344 #define NUM2INTEGRAL scm_num2long_long
5345 #define INTEGRAL2NUM scm_long_long2num
5346 #define INTEGRAL2BIG scm_i_long_long2big
5348 #define ITYPE long long
5349 #define SIZEOF_ITYPE SIZEOF_LONG_LONG
5350 #include "libguile/num2integral.i.c"
5352 #define NUM2INTEGRAL scm_num2ulong_long
5353 #define INTEGRAL2NUM scm_ulong_long2num
5354 #define INTEGRAL2BIG scm_i_ulong_long2big
5356 #define ITYPE unsigned long long
5357 #define SIZEOF_ITYPE SIZEOF_UNSIGNED_LONG_LONG
5358 #include "libguile/num2integral.i.c"
5360 #endif /* SCM_SIZEOF_LONG_LONG != 0 */
5362 #define NUM2FLOAT scm_num2float
5363 #define FLOAT2NUM scm_float2num
5365 #include "libguile/num2float.i.c"
5367 #define NUM2FLOAT scm_num2double
5368 #define FLOAT2NUM scm_double2num
5369 #define FTYPE double
5370 #include "libguile/num2float.i.c"
5375 #define SIZE_MAX ((size_t) (-1))
5378 #define PTRDIFF_MIN \
5379 ((scm_t_ptrdiff) ((scm_t_ptrdiff) 1 \
5380 << ((sizeof (scm_t_ptrdiff) * SCM_CHAR_BIT) - 1)))
5383 #define PTRDIFF_MAX (~ PTRDIFF_MIN)
5386 #define CHECK(type, v) \
5389 if ((v) != scm_num2##type (scm_##type##2num (v), 1, "check_sanity")) \
5409 CHECK (ptrdiff
, -1);
5411 CHECK (short, SHRT_MAX
);
5412 CHECK (short, SHRT_MIN
);
5413 CHECK (ushort
, USHRT_MAX
);
5414 CHECK (int, INT_MAX
);
5415 CHECK (int, INT_MIN
);
5416 CHECK (uint
, UINT_MAX
);
5417 CHECK (long, LONG_MAX
);
5418 CHECK (long, LONG_MIN
);
5419 CHECK (ulong
, ULONG_MAX
);
5420 CHECK (size
, SIZE_MAX
);
5421 CHECK (ptrdiff
, PTRDIFF_MAX
);
5422 CHECK (ptrdiff
, PTRDIFF_MIN
);
5424 #if SCM_SIZEOF_LONG_LONG != 0
5425 CHECK (long_long
, 0LL);
5426 CHECK (ulong_long
, 0ULL);
5427 CHECK (long_long
, -1LL);
5428 CHECK (long_long
, LLONG_MAX
);
5429 CHECK (long_long
, LLONG_MIN
);
5430 CHECK (ulong_long
, ULLONG_MAX
);
5437 scm_internal_catch (SCM_BOOL_T, check_body, &data, check_handler, &data); \
5438 if (!SCM_FALSEP (data)) abort();
5441 check_body (void *data
)
5443 SCM num
= *(SCM
*) data
;
5444 scm_num2ulong (num
, 1, NULL
);
5446 return SCM_UNSPECIFIED
;
5450 check_handler (void *data
, SCM tag
, SCM throw_args
)
5452 SCM
*num
= (SCM
*) data
;
5455 return SCM_UNSPECIFIED
;
5458 SCM_DEFINE (scm_sys_check_number_conversions
, "%check-number-conversions", 0, 0, 0,
5460 "Number conversion sanity checking.")
5461 #define FUNC_NAME s_scm_sys_check_number_conversions
5463 SCM data
= SCM_MAKINUM (-1);
5465 data
= scm_int2num (INT_MIN
);
5467 data
= scm_ulong2num (ULONG_MAX
);
5468 data
= scm_difference (SCM_INUM0
, data
);
5470 data
= scm_ulong2num (ULONG_MAX
);
5471 data
= scm_sum (SCM_MAKINUM (1), data
); data
= scm_difference (SCM_INUM0
, data
);
5473 data
= scm_int2num (-10000); data
= scm_product (data
, data
); data
= scm_product (data
, data
);
5476 return SCM_UNSPECIFIED
;
5485 abs_most_negative_fixnum
= scm_i_long2big (- SCM_MOST_NEGATIVE_FIXNUM
);
5486 scm_permanent_object (abs_most_negative_fixnum
);
5488 mpz_init_set_si (z_negative_one
, -1);
5490 /* It may be possible to tune the performance of some algorithms by using
5491 * the following constants to avoid the creation of bignums. Please, before
5492 * using these values, remember the two rules of program optimization:
5493 * 1st Rule: Don't do it. 2nd Rule (experts only): Don't do it yet. */
5494 scm_c_define ("most-positive-fixnum",
5495 SCM_MAKINUM (SCM_MOST_POSITIVE_FIXNUM
));
5496 scm_c_define ("most-negative-fixnum",
5497 SCM_MAKINUM (SCM_MOST_NEGATIVE_FIXNUM
));
5499 scm_add_feature ("complex");
5500 scm_add_feature ("inexact");
5501 scm_flo0
= scm_make_real (0.0);
5503 scm_dblprec
= (DBL_DIG
> 20) ? 20 : DBL_DIG
;
5505 { /* determine floating point precision */
5507 double fsum
= 1.0 + f
;
5510 if (++scm_dblprec
> 20)
5518 scm_dblprec
= scm_dblprec
- 1;
5520 #endif /* DBL_DIG */
5526 exactly_one_half
= scm_permanent_object (scm_divide (SCM_MAKINUM (1),
5528 #include "libguile/numbers.x"