HCoop / bpt / guile.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
*** empty log message ***
[bpt/guile.git] / libguile / numbers.c
CommitLineData
238ebcef 1/* Copyright (C) 1995,1996,1997,1998,1999,2000,2001,2002,2003,2004 Free Software Foundation, Inc.
ba74ef4e
MV		 2 *
3 * Portions Copyright 1990, 1991, 1992, 1993 by AT&T Bell Laboratories
4 * and Bellcore. See scm_divide.
5 *
f81e080b 6 *
73be1d9e
MV		 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.
0f2d19dd 11 *
73be1d9e
MV		 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.
0f2d19dd 16 *
73be1d9e
MV		 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
20 */
1bbd0b84 21
0f2d19dd 22\f
ca46fb90
RB		 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.
f92e85f7 28 * All objects satisfying SCM_FRACTIONP are never an integer.
ca46fb90
RB		 29 */
30
31/* TODO:
32
33 - see if special casing bignums and reals in integer-exponent when
34 possible (to use mpz_pow and mpf_pow_ui) is faster.
35
36 - look in to better short-circuiting of common cases in
37 integer-expt and elsewhere.
38
39 - see if direct mpz operations can help in ash and elsewhere.
40
41 */
0f2d19dd 42
fa605590
KR		 43/* tell glibc (2.3) to give prototype for C99 trunc() */
44#define _GNU_SOURCE
45
ee33d62a
RB		 46#if HAVE_CONFIG_H
47# include <config.h>
48#endif
49
0f2d19dd 50#include <math.h>
3c9a524f 51#include <ctype.h>
fc194577 52#include <string.h>
ca46fb90 53#include <gmp.h>
f92e85f7 54
a0599745 55#include "libguile/_scm.h"
a0599745
MD		 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"
a0599745
MD		 61
62#include "libguile/validate.h"
63#include "libguile/numbers.h"
1be6b49c 64#include "libguile/deprecation.h"
f4c627b3 65
f92e85f7
MV		 66#include "libguile/eq.h"
67
0f2d19dd 68\f
f4c627b3 69
ca46fb90
RB		 70/*
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...
74
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 \
534c55a9 83 : (((0xfcff & SCM_CELL_TYPE (x)) == scm_tc7_number) ? SCM_TYP16(x) \
ca46fb90
RB		 84 : SCM_I_NUMTAG_NOTNUM)))
85*/
f92e85f7 86/* the macro above will not work as is with fractions */
f4c627b3
DH		 87
88
34d19ef6 89#define SCM_SWAP(x, y) do { SCM __t = x; x = y; y = __t; } while (0)
09fb7599 90
56e55ac7 91/* FLOBUFLEN is the maximum number of characters neccessary for the
3a9809df
DH		 92 * printed or scm_string representation of an inexact number.
93 */
56e55ac7 94#define FLOBUFLEN (10+2*(sizeof(double)/sizeof(char)*SCM_CHAR_BIT*3+9)/10)
3a9809df 95
7351e207
MV		 96#if defined (SCO)
97#if ! defined (HAVE_ISNAN)
98#define HAVE_ISNAN
99static int
100isnan (double x)
101{
102 return (IsNANorINF (x) && NaN (x) && ! IsINF (x)) ? 1 : 0;
103}
0f2d19dd 104#endif
7351e207
MV		 105#if ! defined (HAVE_ISINF)
106#define HAVE_ISINF
107static int
108isinf (double x)
109{
110 return (IsNANorINF (x) && IsINF (x)) ? 1 : 0;
111}
0f2d19dd 112
7351e207 113#endif
e6f3ef58
MD		 114#endif
115
b127c712
KR		 116
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))
123#else
124#define xmpz_cmp_d(z, d) mpz_cmp_d (z, d)
125#endif
126
a98ce907
KR		 127/* For reference, sparc solaris 7 has infinities (IEEE) but doesn't have
128 isinf. It does have finite and isnan though, hence the use of those.
129 fpclass would be a possibility on that system too. */
f92e85f7
MV		 130static int
131xisinf (double x)
132{
133#if defined (HAVE_ISINF)
134 return isinf (x);
135#elif defined (HAVE_FINITE) && defined (HAVE_ISNAN)
136 return (! (finite (x) || isnan (x)));
137#else
138 return 0;
139#endif
140}
141
142static int
143xisnan (double x)
144{
145#if defined (HAVE_ISNAN)
146 return isnan (x);
147#else
148 return 0;
149#endif
150}
151
0f2d19dd
JB		 152\f
153
713a4259 154static mpz_t z_negative_one;
ac0c002c
DH		 155
156\f
157
570b6821 158SCM_C_INLINE_KEYWORD SCM
ca46fb90
RB		 159scm_i_mkbig ()
160{
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));
164 return z;
165}
166
570b6821 167SCM_C_INLINE_KEYWORD static SCM
ca46fb90
RB		 168scm_i_clonebig (SCM src_big, int same_sign_p)
169{
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));
0aacf84e
MD		 173 if (!same_sign_p)
174 mpz_neg (SCM_I_BIG_MPZ (z), SCM_I_BIG_MPZ (z));
ca46fb90
RB		 175 return z;
176}
177
570b6821 178SCM_C_INLINE_KEYWORD int
ca46fb90
RB		 179scm_i_bigcmp (SCM x, SCM y)
180{
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);
185 return result;
186}
187
570b6821 188SCM_C_INLINE_KEYWORD SCM
ca46fb90
RB		 189scm_i_dbl2big (double d)
190{
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);
194 return z;
195}
196
f92e85f7
MV		 197/* Convert a integer in double representation to a SCM number. */
198
199SCM_C_INLINE_KEYWORD SCM
200scm_i_dbl2num (double u)
201{
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".
206
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. */
212
213 /* XXX - what happens when SCM_MOST_POSITIVE_FIXNUM etc is not
214 representable as a double? */
215
216 if (u < (double) (SCM_MOST_POSITIVE_FIXNUM+1)
217 && u >= (double) SCM_MOST_NEGATIVE_FIXNUM)
218 return SCM_MAKINUM ((long) u);
219 else
220 return scm_i_dbl2big (u);
221}
222
089c9a59
KR		 223/* scm_i_big2dbl() rounds to the closest representable double, in accordance
224 with R5RS exact->inexact.
225
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.
229
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.
234
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
237 twos complement.
238
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. */
244
245double
ca46fb90
RB		 246scm_i_big2dbl (SCM b)
247{
089c9a59
KR		 248 double result;
249 size_t bits;
250
251 bits = mpz_sizeinbase (SCM_I_BIG_MPZ (b), 2);
252
253#if __GNU_MP_VERSION < 4 \
254 || (__GNU_MP_VERSION == 4 && __GNU_MP_VERSION_MINOR < 2)
255 {
256 /* GMP prior to 4.2, force truncate towards zero */
257 mpz_t tmp;
258 if (bits > DBL_MANT_DIG)
259 {
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);
264 mpz_clear (tmp);
265 }
266 else
267 {
268 result = mpz_get_d (SCM_I_BIG_MPZ (b));
269 }
270 }
271#else
272 /* GMP 4.2 and up */
273 result = mpz_get_d (SCM_I_BIG_MPZ (b));
274#endif
275
276 if (bits > DBL_MANT_DIG)
277 {
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)))
282 {
283 result += ldexp ((double) mpz_sgn (SCM_I_BIG_MPZ (b)), pos + 1);
284 }
285 }
286
ca46fb90
RB		 287 scm_remember_upto_here_1 (b);
288 return result;
289}
290
570b6821 291SCM_C_INLINE_KEYWORD SCM
ca46fb90
RB		 292scm_i_normbig (SCM b)
293{
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)))
297 {
298 long val = mpz_get_si (SCM_I_BIG_MPZ (b));
299 if (SCM_FIXABLE (val))
300 b = SCM_MAKINUM (val);
301 }
302 return b;
303}
f872b822 304
f92e85f7
MV		 305static SCM_C_INLINE_KEYWORD SCM
306scm_i_mpz2num (mpz_t b)
307{
308 /* convert a mpz number to a SCM number. */
309 if (mpz_fits_slong_p (b))
310 {
311 long val = mpz_get_si (b);
312 if (SCM_FIXABLE (val))
313 return SCM_MAKINUM (val);
314 }
315
316 {
317 SCM z = scm_double_cell (scm_tc16_big, 0, 0, 0);
318 mpz_init_set (SCM_I_BIG_MPZ (z), b);
319 return z;
320 }
321}
322
323/* this is needed when we want scm_divide to make a float, not a ratio, even if passed two ints */
324static SCM scm_divide2real (SCM x, SCM y);
325
326SCM
327scm_make_ratio (SCM numerator, SCM denominator)
c60e130c 328#define FUNC_NAME "make-ratio"
f92e85f7 329{
c60e130c
MV		 330 /* First make sure the arguments are proper.
331 */
f92e85f7
MV		 332 if (SCM_INUMP (denominator))
333 {
334 if (SCM_EQ_P (denominator, SCM_INUM0))
335 scm_num_overflow ("make-ratio");
336 if (SCM_EQ_P (denominator, SCM_MAKINUM(1)))
337 return numerator;
338 }
339 else
340 {
341 if (!(SCM_BIGP(denominator)))
342 SCM_WRONG_TYPE_ARG (2, denominator);
343 }
c60e130c
MV		 344 if (!SCM_INUMP (numerator) && !SCM_BIGP (numerator))
345 SCM_WRONG_TYPE_ARG (1, numerator);
346
347 /* Then flip signs so that the denominator is positive.
348 */
349 if (SCM_NFALSEP (scm_negative_p (denominator)))
350 {
351 numerator = scm_difference (numerator, SCM_UNDEFINED);
352 denominator = scm_difference (denominator, SCM_UNDEFINED);
353 }
354
355 /* Now consider for each of the four fixnum/bignum combinations
356 whether the rational number is really an integer.
357 */
f92e85f7
MV		 358 if (SCM_INUMP (numerator))
359 {
dd5130ca 360 long x = SCM_INUM (numerator);
f92e85f7
MV		 361 if (SCM_EQ_P (numerator, SCM_INUM0))
362 return SCM_INUM0;
363 if (SCM_INUMP (denominator))
364 {
dd5130ca 365 long y;
f92e85f7
MV		 366 y = SCM_INUM (denominator);
367 if (x == y)
368 return SCM_MAKINUM(1);
369 if ((x % y) == 0)
370 return SCM_MAKINUM (x / y);
f92e85f7 371 }
dd5130ca
KR		 372 else
373 {
374 /* When x == SCM_MOST_NEGATIVE_FIXNUM we could have the negative
3271a325
KR		 375 of that value for the denominator, as a bignum. Apart from
376 that case, abs(bignum) > abs(inum) so inum/bignum is not an
377 integer. */
378 if (x == SCM_MOST_NEGATIVE_FIXNUM
379 && mpz_cmp_ui (SCM_I_BIG_MPZ (denominator),
380 - SCM_MOST_NEGATIVE_FIXNUM) == 0)
dd5130ca
KR		 381 return SCM_MAKINUM(-1);
382 }
f92e85f7 383 }
c60e130c 384 else if (SCM_BIGP (numerator))
f92e85f7 385 {
c60e130c
MV		 386 if (SCM_INUMP (denominator))
387 {
388 long yy = SCM_INUM (denominator);
389 if (mpz_divisible_ui_p (SCM_I_BIG_MPZ (numerator), yy))
390 return scm_divide (numerator, denominator);
391 }
392 else
f92e85f7 393 {
c60e130c
MV		 394 if (SCM_EQ_P (numerator, denominator))
395 return SCM_MAKINUM(1);
396 if (mpz_divisible_p (SCM_I_BIG_MPZ (numerator),
397 SCM_I_BIG_MPZ (denominator)))
398 return scm_divide(numerator, denominator);
f92e85f7 399 }
f92e85f7 400 }
c60e130c
MV		 401
402 /* No, it's a proper fraction.
403 */
404 return scm_double_cell (scm_tc16_fraction,
405 SCM_UNPACK (numerator),
406 SCM_UNPACK (denominator), 0);
f92e85f7 407}
c60e130c 408#undef FUNC_NAME
f92e85f7
MV		 409
410static void scm_i_fraction_reduce (SCM z)
411{
412 if (!(SCM_FRACTION_REDUCED (z)))
413 {
414 SCM divisor;
415 divisor = scm_gcd (SCM_FRACTION_NUMERATOR (z), SCM_FRACTION_DENOMINATOR (z));
416 if (!(SCM_EQ_P (divisor, SCM_MAKINUM(1))))
417 {
418 /* is this safe? */
419 SCM_FRACTION_SET_NUMERATOR (z, scm_divide (SCM_FRACTION_NUMERATOR (z), divisor));
420 SCM_FRACTION_SET_DENOMINATOR (z, scm_divide (SCM_FRACTION_DENOMINATOR (z), divisor));
421 }
422 SCM_FRACTION_REDUCED_SET (z);
423 }
424}
425
426double
427scm_i_fraction2double (SCM z)
428{
429 return scm_num2dbl (scm_divide2real (SCM_FRACTION_NUMERATOR (z),
430 SCM_FRACTION_DENOMINATOR (z)),
431 "fraction2real");
432}
433
a1ec6916 434SCM_DEFINE (scm_exact_p, "exact?", 1, 0, 0,
1bbd0b84 435 (SCM x),
942e5b91
MG		 436 "Return @code{#t} if @var{x} is an exact number, @code{#f}\n"
437 "otherwise.")
1bbd0b84 438#define FUNC_NAME s_scm_exact_p
0f2d19dd 439{
0aacf84e
MD		 440 if (SCM_INUMP (x))
441 return SCM_BOOL_T;
442 if (SCM_BIGP (x))
443 return SCM_BOOL_T;
f92e85f7
MV		 444 if (SCM_FRACTIONP (x))
445 return SCM_BOOL_T;
eb927cb9
MV		 446 if (SCM_NUMBERP (x))
447 return SCM_BOOL_F;
448 SCM_WRONG_TYPE_ARG (1, x);
0f2d19dd 449}
1bbd0b84 450#undef FUNC_NAME
0f2d19dd 451
4219f20d 452
a1ec6916 453SCM_DEFINE (scm_odd_p, "odd?", 1, 0, 0,
1bbd0b84 454 (SCM n),
942e5b91
MG		 455 "Return @code{#t} if @var{n} is an odd number, @code{#f}\n"
456 "otherwise.")
1bbd0b84 457#define FUNC_NAME s_scm_odd_p
0f2d19dd 458{
0aacf84e
MD		 459 if (SCM_INUMP (n))
460 {
461 long val = SCM_INUM (n);
462 return SCM_BOOL ((val & 1L) != 0);
463 }
464 else if (SCM_BIGP (n))
465 {
466 int odd_p = mpz_odd_p (SCM_I_BIG_MPZ (n));
467 scm_remember_upto_here_1 (n);
468 return SCM_BOOL (odd_p);
469 }
470 else if (!SCM_FALSEP (scm_inf_p (n)))
7351e207 471 return SCM_BOOL_T;
f92e85f7
MV		 472 else if (SCM_REALP (n))
473 {
474 double rem = fabs (fmod (SCM_REAL_VALUE(n), 2.0));
475 if (rem == 1.0)
476 return SCM_BOOL_T;
477 else if (rem == 0.0)
478 return SCM_BOOL_F;
479 else
480 SCM_WRONG_TYPE_ARG (1, n);
481 }
0aacf84e 482 else
a1a33b0f 483 SCM_WRONG_TYPE_ARG (1, n);
0f2d19dd 484}
1bbd0b84 485#undef FUNC_NAME
0f2d19dd 486
4219f20d 487
a1ec6916 488SCM_DEFINE (scm_even_p, "even?", 1, 0, 0,
1bbd0b84 489 (SCM n),
942e5b91
MG		 490 "Return @code{#t} if @var{n} is an even number, @code{#f}\n"
491 "otherwise.")
1bbd0b84 492#define FUNC_NAME s_scm_even_p
0f2d19dd 493{
0aacf84e
MD		 494 if (SCM_INUMP (n))
495 {
496 long val = SCM_INUM (n);
497 return SCM_BOOL ((val & 1L) == 0);
498 }
499 else if (SCM_BIGP (n))
500 {
501 int even_p = mpz_even_p (SCM_I_BIG_MPZ (n));
502 scm_remember_upto_here_1 (n);
503 return SCM_BOOL (even_p);
504 }
505 else if (!SCM_FALSEP (scm_inf_p (n)))
7351e207 506 return SCM_BOOL_T;
f92e85f7
MV		 507 else if (SCM_REALP (n))
508 {
509 double rem = fabs (fmod (SCM_REAL_VALUE(n), 2.0));
510 if (rem == 1.0)
511 return SCM_BOOL_F;
512 else if (rem == 0.0)
513 return SCM_BOOL_T;
514 else
515 SCM_WRONG_TYPE_ARG (1, n);
516 }
0aacf84e 517 else
a1a33b0f 518 SCM_WRONG_TYPE_ARG (1, n);
0f2d19dd 519}
1bbd0b84 520#undef FUNC_NAME
0f2d19dd 521
7351e207
MV		 522SCM_DEFINE (scm_inf_p, "inf?", 1, 0, 0,
523 (SCM n),
524 "Return @code{#t} if @var{n} is infinite, @code{#f}\n"
525 "otherwise.")
526#define FUNC_NAME s_scm_inf_p
527{
0aacf84e 528 if (SCM_REALP (n))
7351e207 529 return SCM_BOOL (xisinf (SCM_REAL_VALUE (n)));
0aacf84e 530 else if (SCM_COMPLEXP (n))
7351e207
MV		 531 return SCM_BOOL (xisinf (SCM_COMPLEX_REAL (n))
532 || xisinf (SCM_COMPLEX_IMAG (n)));
0aacf84e 533 else
7351e207 534 return SCM_BOOL_F;
7351e207
MV		 535}
536#undef FUNC_NAME
537
538SCM_DEFINE (scm_nan_p, "nan?", 1, 0, 0,
539 (SCM n),
540 "Return @code{#t} if @var{n} is a NaN, @code{#f}\n"
541 "otherwise.")
542#define FUNC_NAME s_scm_nan_p
543{
0aacf84e 544 if (SCM_REALP (n))
7351e207 545 return SCM_BOOL (xisnan (SCM_REAL_VALUE (n)));
0aacf84e 546 else if (SCM_COMPLEXP (n))
7351e207
MV		 547 return SCM_BOOL (xisnan (SCM_COMPLEX_REAL (n))
548 || xisnan (SCM_COMPLEX_IMAG (n)));
0aacf84e 549 else
7351e207 550 return SCM_BOOL_F;
7351e207
MV		 551}
552#undef FUNC_NAME
553
554/* Guile's idea of infinity. */
555static double guile_Inf;
556
557/* Guile's idea of not a number. */
558static double guile_NaN;
559
560static void
561guile_ieee_init (void)
562{
563#if defined (HAVE_ISINF) || defined (HAVE_FINITE)
564
565/* Some version of gcc on some old version of Linux used to crash when
566 trying to make Inf and NaN. */
567
240a27d2
KR		 568#ifdef INFINITY
569 /* C99 INFINITY, when available.
570 FIXME: The standard allows for INFINITY to be something that overflows
571 at compile time. We ought to have a configure test to check for that
572 before trying to use it. (But in practice we believe this is not a
573 problem on any system guile is likely to target.) */
574 guile_Inf = INFINITY;
575#elif HAVE_DINFINITY
576 /* OSF */
7351e207
MV		 577 extern unsigned int DINFINITY[2];
578 guile_Inf = (*(X_CAST(double *, DINFINITY)));
579#else
580 double tmp = 1e+10;
581 guile_Inf = tmp;
582 for (;;)
583 {
584 guile_Inf *= 1e+10;
585 if (guile_Inf == tmp)
586 break;
587 tmp = guile_Inf;
588 }
589#endif
590
591#endif
592
593#if defined (HAVE_ISNAN)
594
240a27d2
KR		 595#ifdef NAN
596 /* C99 NAN, when available */
597 guile_NaN = NAN;
598#elif HAVE_DQNAN
599 /* OSF */
7351e207
MV		 600 extern unsigned int DQNAN[2];
601 guile_NaN = (*(X_CAST(double *, DQNAN)));
602#else
603 guile_NaN = guile_Inf / guile_Inf;
604#endif
605
606#endif
607}
608
609SCM_DEFINE (scm_inf, "inf", 0, 0, 0,
610 (void),
611 "Return Inf.")
612#define FUNC_NAME s_scm_inf
613{
614 static int initialized = 0;
615 if (! initialized)
616 {
617 guile_ieee_init ();
618 initialized = 1;
619 }
620 return scm_make_real (guile_Inf);
621}
622#undef FUNC_NAME
623
624SCM_DEFINE (scm_nan, "nan", 0, 0, 0,
625 (void),
626 "Return NaN.")
627#define FUNC_NAME s_scm_nan
628{
629 static int initialized = 0;
0aacf84e 630 if (!initialized)
7351e207
MV		 631 {
632 guile_ieee_init ();
633 initialized = 1;
634 }
635 return scm_make_real (guile_NaN);
636}
637#undef FUNC_NAME
638
4219f20d 639
a48d60b1
MD		 640SCM_PRIMITIVE_GENERIC (scm_abs, "abs", 1, 0, 0,
641 (SCM x),
642 "Return the absolute value of @var{x}.")
643#define FUNC_NAME
0f2d19dd 644{
0aacf84e
MD		 645 if (SCM_INUMP (x))
646 {
647 long int xx = SCM_INUM (x);
648 if (xx >= 0)
649 return x;
650 else if (SCM_POSFIXABLE (-xx))
651 return SCM_MAKINUM (-xx);
652 else
653 return scm_i_long2big (-xx);
4219f20d 654 }
0aacf84e
MD		 655 else if (SCM_BIGP (x))
656 {
657 const int sgn = mpz_sgn (SCM_I_BIG_MPZ (x));
658 if (sgn < 0)
659 return scm_i_clonebig (x, 0);
660 else
661 return x;
4219f20d 662 }
0aacf84e 663 else if (SCM_REALP (x))
ae38324d
KR		 664 {
665 /* note that if x is a NaN then xx<0 is false so we return x unchanged */
666 double xx = SCM_REAL_VALUE (x);
667 if (xx < 0.0)
668 return scm_make_real (-xx);
669 else
670 return x;
671 }
f92e85f7
MV		 672 else if (SCM_FRACTIONP (x))
673 {
674 if (SCM_FALSEP (scm_negative_p (SCM_FRACTION_NUMERATOR (x))))
675 return x;
676 return scm_make_ratio (scm_difference (SCM_FRACTION_NUMERATOR (x), SCM_UNDEFINED),
677 SCM_FRACTION_DENOMINATOR (x));
678 }
0aacf84e 679 else
a48d60b1 680 SCM_WTA_DISPATCH_1 (g_scm_abs, x, 1, s_scm_abs);
0f2d19dd 681}
a48d60b1 682#undef FUNC_NAME
0f2d19dd 683
4219f20d 684
9de33deb 685SCM_GPROC (s_quotient, "quotient", 2, 0, 0, scm_quotient, g_quotient);
942e5b91
MG		 686/* "Return the quotient of the numbers @var{x} and @var{y}."
687 */
0f2d19dd 688SCM
6e8d25a6 689scm_quotient (SCM x, SCM y)
0f2d19dd 690{
0aacf84e
MD		 691 if (SCM_INUMP (x))
692 {
693 long xx = SCM_INUM (x);
694 if (SCM_INUMP (y))
695 {
696 long yy = SCM_INUM (y);
697 if (yy == 0)
698 scm_num_overflow (s_quotient);
699 else
700 {
701 long z = xx / yy;
702 if (SCM_FIXABLE (z))
703 return SCM_MAKINUM (z);
704 else
705 return scm_i_long2big (z);
706 }
828865c3 707 }
0aacf84e 708 else if (SCM_BIGP (y))
ac0c002c 709 {
0aacf84e 710 if ((SCM_INUM (x) == SCM_MOST_NEGATIVE_FIXNUM)
4dc09ee4
KR		 711 && (mpz_cmp_ui (SCM_I_BIG_MPZ (y),
712 - SCM_MOST_NEGATIVE_FIXNUM) == 0))
713 {
714 /* Special case: x == fixnum-min && y == abs (fixnum-min) */
715 scm_remember_upto_here_1 (y);
716 return SCM_MAKINUM (-1);
717 }
0aacf84e
MD		 718 else
719 return SCM_MAKINUM (0);
ac0c002c
DH		 720 }
721 else
0aacf84e 722 SCM_WTA_DISPATCH_2 (g_quotient, x, y, SCM_ARG2, s_quotient);
828865c3 723 }
0aacf84e
MD		 724 else if (SCM_BIGP (x))
725 {
726 if (SCM_INUMP (y))
727 {
728 long yy = SCM_INUM (y);
729 if (yy == 0)
730 scm_num_overflow (s_quotient);
731 else if (yy == 1)
732 return x;
733 else
734 {
735 SCM result = scm_i_mkbig ();
736 if (yy < 0)
737 {
738 mpz_tdiv_q_ui (SCM_I_BIG_MPZ (result),
739 SCM_I_BIG_MPZ (x),
740 - yy);
741 mpz_neg (SCM_I_BIG_MPZ (result), SCM_I_BIG_MPZ (result));
742 }
743 else
744 mpz_tdiv_q_ui (SCM_I_BIG_MPZ (result), SCM_I_BIG_MPZ (x), yy);
745 scm_remember_upto_here_1 (x);
746 return scm_i_normbig (result);
747 }
748 }
749 else if (SCM_BIGP (y))
750 {
751 SCM result = scm_i_mkbig ();
752 mpz_tdiv_q (SCM_I_BIG_MPZ (result),
753 SCM_I_BIG_MPZ (x),
754 SCM_I_BIG_MPZ (y));
755 scm_remember_upto_here_2 (x, y);
756 return scm_i_normbig (result);
757 }
758 else
759 SCM_WTA_DISPATCH_2 (g_quotient, x, y, SCM_ARG2, s_quotient);
f872b822 760 }
0aacf84e 761 else
89a7e495 762 SCM_WTA_DISPATCH_2 (g_quotient, x, y, SCM_ARG1, s_quotient);
0f2d19dd
JB		 763}
764
9de33deb 765SCM_GPROC (s_remainder, "remainder", 2, 0, 0, scm_remainder, g_remainder);
942e5b91
MG		 766/* "Return the remainder of the numbers @var{x} and @var{y}.\n"
767 * "@lisp\n"
768 * "(remainder 13 4) @result{} 1\n"
769 * "(remainder -13 4) @result{} -1\n"
770 * "@end lisp"
771 */
0f2d19dd 772SCM
6e8d25a6 773scm_remainder (SCM x, SCM y)
0f2d19dd 774{
0aacf84e
MD		 775 if (SCM_INUMP (x))
776 {
777 if (SCM_INUMP (y))
778 {
779 long yy = SCM_INUM (y);
780 if (yy == 0)
781 scm_num_overflow (s_remainder);
782 else
783 {
784 long z = SCM_INUM (x) % yy;
785 return SCM_MAKINUM (z);
786 }
787 }
788 else if (SCM_BIGP (y))
ac0c002c 789 {
0aacf84e 790 if ((SCM_INUM (x) == SCM_MOST_NEGATIVE_FIXNUM)
4dc09ee4
KR		 791 && (mpz_cmp_ui (SCM_I_BIG_MPZ (y),
792 - SCM_MOST_NEGATIVE_FIXNUM) == 0))
793 {
794 /* Special case: x == fixnum-min && y == abs (fixnum-min) */
795 scm_remember_upto_here_1 (y);
796 return SCM_MAKINUM (0);
797 }
0aacf84e
MD		 798 else
799 return x;
ac0c002c
DH		 800 }
801 else
0aacf84e 802 SCM_WTA_DISPATCH_2 (g_remainder, x, y, SCM_ARG2, s_remainder);
89a7e495 803 }
0aacf84e
MD		 804 else if (SCM_BIGP (x))
805 {
806 if (SCM_INUMP (y))
807 {
808 long yy = SCM_INUM (y);
809 if (yy == 0)
810 scm_num_overflow (s_remainder);
811 else
812 {
813 SCM result = scm_i_mkbig ();
814 if (yy < 0)
815 yy = - yy;
816 mpz_tdiv_r_ui (SCM_I_BIG_MPZ (result), SCM_I_BIG_MPZ(x), yy);
817 scm_remember_upto_here_1 (x);
818 return scm_i_normbig (result);
819 }
820 }
821 else if (SCM_BIGP (y))
822 {
823 SCM result = scm_i_mkbig ();
824 mpz_tdiv_r (SCM_I_BIG_MPZ (result),
825 SCM_I_BIG_MPZ (x),
826 SCM_I_BIG_MPZ (y));
827 scm_remember_upto_here_2 (x, y);
828 return scm_i_normbig (result);
829 }
830 else
831 SCM_WTA_DISPATCH_2 (g_remainder, x, y, SCM_ARG2, s_remainder);
f872b822 832 }
0aacf84e 833 else
89a7e495 834 SCM_WTA_DISPATCH_2 (g_remainder, x, y, SCM_ARG1, s_remainder);
0f2d19dd
JB		 835}
836
89a7e495 837
9de33deb 838SCM_GPROC (s_modulo, "modulo", 2, 0, 0, scm_modulo, g_modulo);
942e5b91
MG		 839/* "Return the modulo of the numbers @var{x} and @var{y}.\n"
840 * "@lisp\n"
841 * "(modulo 13 4) @result{} 1\n"
842 * "(modulo -13 4) @result{} 3\n"
843 * "@end lisp"
844 */
0f2d19dd 845SCM
6e8d25a6 846scm_modulo (SCM x, SCM y)
0f2d19dd 847{
0aacf84e
MD		 848 if (SCM_INUMP (x))
849 {
850 long xx = SCM_INUM (x);
851 if (SCM_INUMP (y))
852 {
853 long yy = SCM_INUM (y);
854 if (yy == 0)
855 scm_num_overflow (s_modulo);
856 else
857 {
858 /* FIXME: I think this may be a bug on some arches -- results
859 of % with negative second arg are undefined... */
860 long z = xx % yy;
861 long result;
862
863 if (yy < 0)
864 {
865 if (z > 0)
866 result = z + yy;
867 else
868 result = z;
869 }
870 else
871 {
872 if (z < 0)
873 result = z + yy;
874 else
875 result = z;
876 }
877 return SCM_MAKINUM (result);
878 }
879 }
880 else if (SCM_BIGP (y))
881 {
882 int sgn_y = mpz_sgn (SCM_I_BIG_MPZ (y));
0aacf84e
MD		 883 {
884 mpz_t z_x;
885 SCM result;
886
887 if (sgn_y < 0)
888 {
889 SCM pos_y = scm_i_clonebig (y, 0);
890 /* do this after the last scm_op */
891 mpz_init_set_si (z_x, xx);
892 result = pos_y; /* re-use this bignum */
893 mpz_mod (SCM_I_BIG_MPZ (result),
894 z_x,
895 SCM_I_BIG_MPZ (pos_y));
896 scm_remember_upto_here_1 (pos_y);
897 }
898 else
899 {
900 result = scm_i_mkbig ();
901 /* do this after the last scm_op */
902 mpz_init_set_si (z_x, xx);
903 mpz_mod (SCM_I_BIG_MPZ (result),
904 z_x,
905 SCM_I_BIG_MPZ (y));
906 scm_remember_upto_here_1 (y);
907 }
ca46fb90 908
0aacf84e
MD		 909 if ((sgn_y < 0) && mpz_sgn (SCM_I_BIG_MPZ (result)) != 0)
910 mpz_add (SCM_I_BIG_MPZ (result),
911 SCM_I_BIG_MPZ (y),
912 SCM_I_BIG_MPZ (result));
913 scm_remember_upto_here_1 (y);
914 /* and do this before the next one */
915 mpz_clear (z_x);
916 return scm_i_normbig (result);
917 }
918 }
919 else
920 SCM_WTA_DISPATCH_2 (g_modulo, x, y, SCM_ARG2, s_modulo);
f872b822 921 }
0aacf84e
MD		 922 else if (SCM_BIGP (x))
923 {
924 if (SCM_INUMP (y))
925 {
926 long yy = SCM_INUM (y);
927 if (yy == 0)
928 scm_num_overflow (s_modulo);
929 else
930 {
931 SCM result = scm_i_mkbig ();
932 mpz_mod_ui (SCM_I_BIG_MPZ (result),
933 SCM_I_BIG_MPZ (x),
934 (yy < 0) ? - yy : yy);
935 scm_remember_upto_here_1 (x);
936 if ((yy < 0) && (mpz_sgn (SCM_I_BIG_MPZ (result)) != 0))
937 mpz_sub_ui (SCM_I_BIG_MPZ (result),
938 SCM_I_BIG_MPZ (result),
939 - yy);
940 return scm_i_normbig (result);
941 }
942 }
943 else if (SCM_BIGP (y))
944 {
0aacf84e
MD		 945 {
946 SCM result = scm_i_mkbig ();
947 int y_sgn = mpz_sgn (SCM_I_BIG_MPZ (y));
948 SCM pos_y = scm_i_clonebig (y, y_sgn >= 0);
949 mpz_mod (SCM_I_BIG_MPZ (result),
950 SCM_I_BIG_MPZ (x),
951 SCM_I_BIG_MPZ (pos_y));
ca46fb90 952
0aacf84e
MD		 953 scm_remember_upto_here_1 (x);
954 if ((y_sgn < 0) && (mpz_sgn (SCM_I_BIG_MPZ (result)) != 0))
955 mpz_add (SCM_I_BIG_MPZ (result),
956 SCM_I_BIG_MPZ (y),
957 SCM_I_BIG_MPZ (result));
958 scm_remember_upto_here_2 (y, pos_y);
959 return scm_i_normbig (result);
960 }
961 }
962 else
963 SCM_WTA_DISPATCH_2 (g_modulo, x, y, SCM_ARG2, s_modulo);
828865c3 964 }
0aacf84e 965 else
09fb7599 966 SCM_WTA_DISPATCH_2 (g_modulo, x, y, SCM_ARG1, s_modulo);
0f2d19dd
JB		 967}
968
9de33deb 969SCM_GPROC1 (s_gcd, "gcd", scm_tc7_asubr, scm_gcd, g_gcd);
942e5b91
MG		 970/* "Return the greatest common divisor of all arguments.\n"
971 * "If called without arguments, 0 is returned."
972 */
0f2d19dd 973SCM
6e8d25a6 974scm_gcd (SCM x, SCM y)
0f2d19dd 975{
ca46fb90 976 if (SCM_UNBNDP (y))
0aacf84e 977 return SCM_UNBNDP (x) ? SCM_INUM0 : x;
ca46fb90
RB		 978
979 if (SCM_INUMP (x))
980 {
981 if (SCM_INUMP (y))
982 {
983 long xx = SCM_INUM (x);
984 long yy = SCM_INUM (y);
985 long u = xx < 0 ? -xx : xx;
986 long v = yy < 0 ? -yy : yy;
987 long result;
0aacf84e
MD		 988 if (xx == 0)
989 result = v;
990 else if (yy == 0)
991 result = u;
992 else
993 {
994 long k = 1;
995 long t;
996 /* Determine a common factor 2^k */
997 while (!(1 & (u | v)))
998 {
999 k <<= 1;
1000 u >>= 1;
1001 v >>= 1;
1002 }
1003 /* Now, any factor 2^n can be eliminated */
1004 if (u & 1)
1005 t = -v;
1006 else
1007 {
1008 t = u;
1009 b3:
1010 t = SCM_SRS (t, 1);
1011 }
1012 if (!(1 & t))
1013 goto b3;
1014 if (t > 0)
1015 u = t;
1016 else
1017 v = -t;
1018 t = u - v;
1019 if (t != 0)
1020 goto b3;
1021 result = u * k;
1022 }
1023 return (SCM_POSFIXABLE (result)
1024 ? SCM_MAKINUM (result)
1025 : scm_i_long2big (result));
ca46fb90
RB		 1026 }
1027 else if (SCM_BIGP (y))
1028 {
0bff4dce
KR		 1029 SCM_SWAP (x, y);
1030 goto big_inum;
ca46fb90
RB		 1031 }
1032 else
1033 SCM_WTA_DISPATCH_2 (g_gcd, x, y, SCM_ARG2, s_gcd);
f872b822 1034 }
ca46fb90
RB		 1035 else if (SCM_BIGP (x))
1036 {
1037 if (SCM_INUMP (y))
1038 {
1039 unsigned long result;
0bff4dce
KR		 1040 long yy;
1041 big_inum:
1042 yy = SCM_INUM (y);
8c5b0afc
KR		 1043 if (yy == 0)
1044 return scm_abs (x);
0aacf84e
MD		 1045 if (yy < 0)
1046 yy = -yy;
ca46fb90
RB		 1047 result = mpz_gcd_ui (NULL, SCM_I_BIG_MPZ (x), yy);
1048 scm_remember_upto_here_1 (x);
0aacf84e
MD		 1049 return (SCM_POSFIXABLE (result)
1050 ? SCM_MAKINUM (result)
1051 : scm_ulong2num (result));
ca46fb90
RB		 1052 }
1053 else if (SCM_BIGP (y))
1054 {
1055 SCM result = scm_i_mkbig ();
0aacf84e
MD		 1056 mpz_gcd (SCM_I_BIG_MPZ (result),
1057 SCM_I_BIG_MPZ (x),
1058 SCM_I_BIG_MPZ (y));
1059 scm_remember_upto_here_2 (x, y);
ca46fb90
RB		 1060 return scm_i_normbig (result);
1061 }
1062 else
1063 SCM_WTA_DISPATCH_2 (g_gcd, x, y, SCM_ARG2, s_gcd);
09fb7599 1064 }
ca46fb90 1065 else
09fb7599 1066 SCM_WTA_DISPATCH_2 (g_gcd, x, y, SCM_ARG1, s_gcd);
0f2d19dd
JB		 1067}
1068
9de33deb 1069SCM_GPROC1 (s_lcm, "lcm", scm_tc7_asubr, scm_lcm, g_lcm);
942e5b91
MG		 1070/* "Return the least common multiple of the arguments.\n"
1071 * "If called without arguments, 1 is returned."
1072 */
0f2d19dd 1073SCM
6e8d25a6 1074scm_lcm (SCM n1, SCM n2)
0f2d19dd 1075{
ca46fb90
RB		 1076 if (SCM_UNBNDP (n2))
1077 {
1078 if (SCM_UNBNDP (n1))
1079 return SCM_MAKINUM (1L);
09fb7599
DH		 1080 n2 = SCM_MAKINUM (1L);
1081 }
09fb7599 1082
09fb7599 1083 SCM_GASSERT2 (SCM_INUMP (n1) || SCM_BIGP (n1),
ca46fb90 1084 g_lcm, n1, n2, SCM_ARG1, s_lcm);
09fb7599 1085 SCM_GASSERT2 (SCM_INUMP (n2) || SCM_BIGP (n2),
ca46fb90 1086 g_lcm, n1, n2, SCM_ARGn, s_lcm);
09fb7599 1087
ca46fb90
RB		 1088 if (SCM_INUMP (n1))
1089 {
1090 if (SCM_INUMP (n2))
1091 {
1092 SCM d = scm_gcd (n1, n2);
1093 if (SCM_EQ_P (d, SCM_INUM0))
1094 return d;
1095 else
1096 return scm_abs (scm_product (n1, scm_quotient (n2, d)));
1097 }
1098 else
1099 {
1100 /* inum n1, big n2 */
1101 inumbig:
1102 {
1103 SCM result = scm_i_mkbig ();
1104 long nn1 = SCM_INUM (n1);
1105 if (nn1 == 0) return SCM_INUM0;
1106 if (nn1 < 0) nn1 = - nn1;
1107 mpz_lcm_ui (SCM_I_BIG_MPZ (result), SCM_I_BIG_MPZ (n2), nn1);
1108 scm_remember_upto_here_1 (n2);
1109 return result;
1110 }
1111 }
1112 }
1113 else
1114 {
1115 /* big n1 */
1116 if (SCM_INUMP (n2))
1117 {
1118 SCM_SWAP (n1, n2);
1119 goto inumbig;
1120 }
1121 else
1122 {
1123 SCM result = scm_i_mkbig ();
1124 mpz_lcm(SCM_I_BIG_MPZ (result),
1125 SCM_I_BIG_MPZ (n1),
1126 SCM_I_BIG_MPZ (n2));
1127 scm_remember_upto_here_2(n1, n2);
1128 /* shouldn't need to normalize b/c lcm of 2 bigs should be big */
1129 return result;
1130 }
f872b822 1131 }
0f2d19dd
JB		 1132}
1133
0f2d19dd 1134#ifndef scm_long2num
c1bfcf60
GB		 1135#define SCM_LOGOP_RETURN(x) scm_ulong2num(x)
1136#else
1137#define SCM_LOGOP_RETURN(x) SCM_MAKINUM(x)
1138#endif
1139
8a525303
GB		 1140/* Emulating 2's complement bignums with sign magnitude arithmetic:
1141
1142 Logand:
1143 X Y Result Method:
1144 (len)
1145 + + + x (map digit:logand X Y)
1146 + - + x (map digit:logand X (lognot (+ -1 Y)))
1147 - + + y (map digit:logand (lognot (+ -1 X)) Y)
1148 - - - (+ 1 (map digit:logior (+ -1 X) (+ -1 Y)))
1149
1150 Logior:
1151 X Y Result Method:
1152
1153 + + + (map digit:logior X Y)
1154 + - - y (+ 1 (map digit:logand (lognot X) (+ -1 Y)))
1155 - + - x (+ 1 (map digit:logand (+ -1 X) (lognot Y)))
1156 - - - x (+ 1 (map digit:logand (+ -1 X) (+ -1 Y)))
1157
1158 Logxor:
1159 X Y Result Method:
1160
1161 + + + (map digit:logxor X Y)
1162 + - - (+ 1 (map digit:logxor X (+ -1 Y)))
1163 - + - (+ 1 (map digit:logxor (+ -1 X) Y))
1164 - - + (map digit:logxor (+ -1 X) (+ -1 Y))
1165
1166 Logtest:
1167 X Y Result
1168
1169 + + (any digit:logand X Y)
1170 + - (any digit:logand X (lognot (+ -1 Y)))
1171 - + (any digit:logand (lognot (+ -1 X)) Y)
1172 - - #t
1173
1174*/
1175
c3ee7520 1176SCM_DEFINE1 (scm_logand, "logand", scm_tc7_asubr,
1bbd0b84 1177 (SCM n1, SCM n2),
3c3db128
GH		 1178 "Return the bitwise AND of the integer arguments.\n\n"
1179 "@lisp\n"
1180 "(logand) @result{} -1\n"
1181 "(logand 7) @result{} 7\n"
535f2a51 1182 "(logand #b111 #b011 #b001) @result{} 1\n"
3c3db128 1183 "@end lisp")
1bbd0b84 1184#define FUNC_NAME s_scm_logand
0f2d19dd 1185{
9a00c9fc
DH		 1186 long int nn1;
1187
0aacf84e
MD		 1188 if (SCM_UNBNDP (n2))
1189 {
1190 if (SCM_UNBNDP (n1))
1191 return SCM_MAKINUM (-1);
1192 else if (!SCM_NUMBERP (n1))
1193 SCM_WRONG_TYPE_ARG (SCM_ARG1, n1);
1194 else if (SCM_NUMBERP (n1))
1195 return n1;
1196 else
1197 SCM_WRONG_TYPE_ARG (SCM_ARG1, n1);
d28da049 1198 }
09fb7599 1199
0aacf84e
MD		 1200 if (SCM_INUMP (n1))
1201 {
9a00c9fc 1202 nn1 = SCM_INUM (n1);
0aacf84e
MD		 1203 if (SCM_INUMP (n2))
1204 {
1205 long nn2 = SCM_INUM (n2);
1206 return SCM_MAKINUM (nn1 & nn2);
1207 }
1208 else if SCM_BIGP (n2)
1209 {
1210 intbig:
1211 if (n1 == 0)
1212 return SCM_INUM0;
1213 {
1214 SCM result_z = scm_i_mkbig ();
1215 mpz_t nn1_z;
1216 mpz_init_set_si (nn1_z, nn1);
1217 mpz_and (SCM_I_BIG_MPZ (result_z), nn1_z, SCM_I_BIG_MPZ (n2));
1218 scm_remember_upto_here_1 (n2);
1219 mpz_clear (nn1_z);
1220 return scm_i_normbig (result_z);
1221 }
1222 }
1223 else
1224 SCM_WRONG_TYPE_ARG (SCM_ARG2, n2);
1225 }
1226 else if (SCM_BIGP (n1))
1227 {
1228 if (SCM_INUMP (n2))
1229 {
1230 SCM_SWAP (n1, n2);
1231 nn1 = SCM_INUM (n1);
1232 goto intbig;
1233 }
1234 else if (SCM_BIGP (n2))
1235 {
1236 SCM result_z = scm_i_mkbig ();
1237 mpz_and (SCM_I_BIG_MPZ (result_z),
1238 SCM_I_BIG_MPZ (n1),
1239 SCM_I_BIG_MPZ (n2));
1240 scm_remember_upto_here_2 (n1, n2);
1241 return scm_i_normbig (result_z);
1242 }
1243 else
1244 SCM_WRONG_TYPE_ARG (SCM_ARG2, n2);
09fb7599 1245 }
0aacf84e 1246 else
09fb7599 1247 SCM_WRONG_TYPE_ARG (SCM_ARG1, n1);
0f2d19dd 1248}
1bbd0b84 1249#undef FUNC_NAME
0f2d19dd 1250
09fb7599 1251
c3ee7520 1252SCM_DEFINE1 (scm_logior, "logior", scm_tc7_asubr,
1bbd0b84 1253 (SCM n1, SCM n2),
3c3db128
GH		 1254 "Return the bitwise OR of the integer arguments.\n\n"
1255 "@lisp\n"
1256 "(logior) @result{} 0\n"
1257 "(logior 7) @result{} 7\n"
1258 "(logior #b000 #b001 #b011) @result{} 3\n"
1e6808ea 1259 "@end lisp")
1bbd0b84 1260#define FUNC_NAME s_scm_logior
0f2d19dd 1261{
9a00c9fc
DH		 1262 long int nn1;
1263
0aacf84e
MD		 1264 if (SCM_UNBNDP (n2))
1265 {
1266 if (SCM_UNBNDP (n1))
1267 return SCM_INUM0;
1268 else if (SCM_NUMBERP (n1))
1269 return n1;
1270 else
1271 SCM_WRONG_TYPE_ARG (SCM_ARG1, n1);
d28da049 1272 }
09fb7599 1273
0aacf84e
MD		 1274 if (SCM_INUMP (n1))
1275 {
9a00c9fc 1276 nn1 = SCM_INUM (n1);
0aacf84e
MD		 1277 if (SCM_INUMP (n2))
1278 {
1279 long nn2 = SCM_INUM (n2);
1280 return SCM_MAKINUM (nn1 | nn2);
1281 }
1282 else if (SCM_BIGP (n2))
1283 {
1284 intbig:
1285 if (nn1 == 0)
1286 return n2;
1287 {
1288 SCM result_z = scm_i_mkbig ();
1289 mpz_t nn1_z;
1290 mpz_init_set_si (nn1_z, nn1);
1291 mpz_ior (SCM_I_BIG_MPZ (result_z), nn1_z, SCM_I_BIG_MPZ (n2));
1292 scm_remember_upto_here_1 (n2);
1293 mpz_clear (nn1_z);
1294 return result_z;
1295 }
1296 }
1297 else
1298 SCM_WRONG_TYPE_ARG (SCM_ARG2, n2);
1299 }
1300 else if (SCM_BIGP (n1))
1301 {
1302 if (SCM_INUMP (n2))
1303 {
1304 SCM_SWAP (n1, n2);
1305 nn1 = SCM_INUM (n1);
1306 goto intbig;
1307 }
1308 else if (SCM_BIGP (n2))
1309 {
1310 SCM result_z = scm_i_mkbig ();
1311 mpz_ior (SCM_I_BIG_MPZ (result_z),
1312 SCM_I_BIG_MPZ (n1),
1313 SCM_I_BIG_MPZ (n2));
1314 scm_remember_upto_here_2 (n1, n2);
1315 return result_z;
1316 }
1317 else
1318 SCM_WRONG_TYPE_ARG (SCM_ARG2, n2);
09fb7599 1319 }
0aacf84e 1320 else
09fb7599 1321 SCM_WRONG_TYPE_ARG (SCM_ARG1, n1);
0f2d19dd 1322}
1bbd0b84 1323#undef FUNC_NAME
0f2d19dd 1324
09fb7599 1325
c3ee7520 1326SCM_DEFINE1 (scm_logxor, "logxor", scm_tc7_asubr,
1bbd0b84 1327 (SCM n1, SCM n2),
3c3db128
GH		 1328 "Return the bitwise XOR of the integer arguments. A bit is\n"
1329 "set in the result if it is set in an odd number of arguments.\n"
1330 "@lisp\n"
1331 "(logxor) @result{} 0\n"
1332 "(logxor 7) @result{} 7\n"
1333 "(logxor #b000 #b001 #b011) @result{} 2\n"
1334 "(logxor #b000 #b001 #b011 #b011) @result{} 1\n"
1e6808ea 1335 "@end lisp")
1bbd0b84 1336#define FUNC_NAME s_scm_logxor
0f2d19dd 1337{
9a00c9fc
DH		 1338 long int nn1;
1339
0aacf84e
MD		 1340 if (SCM_UNBNDP (n2))
1341 {
1342 if (SCM_UNBNDP (n1))
1343 return SCM_INUM0;
1344 else if (SCM_NUMBERP (n1))
1345 return n1;
1346 else
1347 SCM_WRONG_TYPE_ARG (SCM_ARG1, n1);
d28da049 1348 }
09fb7599 1349
0aacf84e
MD		 1350 if (SCM_INUMP (n1))
1351 {
9a00c9fc 1352 nn1 = SCM_INUM (n1);
0aacf84e
MD		 1353 if (SCM_INUMP (n2))
1354 {
1355 long nn2 = SCM_INUM (n2);
1356 return SCM_MAKINUM (nn1 ^ nn2);
1357 }
1358 else if (SCM_BIGP (n2))
1359 {
1360 intbig:
1361 {
1362 SCM result_z = scm_i_mkbig ();
1363 mpz_t nn1_z;
1364 mpz_init_set_si (nn1_z, nn1);
1365 mpz_xor (SCM_I_BIG_MPZ (result_z), nn1_z, SCM_I_BIG_MPZ (n2));
1366 scm_remember_upto_here_1 (n2);
1367 mpz_clear (nn1_z);
1368 return scm_i_normbig (result_z);
1369 }
1370 }
1371 else
1372 SCM_WRONG_TYPE_ARG (SCM_ARG2, n2);
1373 }
1374 else if (SCM_BIGP (n1))
1375 {
1376 if (SCM_INUMP (n2))
1377 {
1378 SCM_SWAP (n1, n2);
1379 nn1 = SCM_INUM (n1);
1380 goto intbig;
1381 }
1382 else if (SCM_BIGP (n2))
1383 {
1384 SCM result_z = scm_i_mkbig ();
1385 mpz_xor (SCM_I_BIG_MPZ (result_z),
1386 SCM_I_BIG_MPZ (n1),
1387 SCM_I_BIG_MPZ (n2));
1388 scm_remember_upto_here_2 (n1, n2);
1389 return scm_i_normbig (result_z);
1390 }
1391 else
1392 SCM_WRONG_TYPE_ARG (SCM_ARG2, n2);
09fb7599 1393 }
0aacf84e 1394 else
09fb7599 1395 SCM_WRONG_TYPE_ARG (SCM_ARG1, n1);
0f2d19dd 1396}
1bbd0b84 1397#undef FUNC_NAME
0f2d19dd 1398
09fb7599 1399
a1ec6916 1400SCM_DEFINE (scm_logtest, "logtest", 2, 0, 0,
1e6808ea
MG		 1401 (SCM j, SCM k),
1402 "@lisp\n"
b380b885
MD		 1403 "(logtest j k) @equiv{} (not (zero? (logand j k)))\n\n"
1404 "(logtest #b0100 #b1011) @result{} #f\n"
1405 "(logtest #b0100 #b0111) @result{} #t\n"
1e6808ea 1406 "@end lisp")
1bbd0b84 1407#define FUNC_NAME s_scm_logtest
0f2d19dd 1408{
1e6808ea 1409 long int nj;
9a00c9fc 1410
0aacf84e
MD		 1411 if (SCM_INUMP (j))
1412 {
1e6808ea 1413 nj = SCM_INUM (j);
0aacf84e
MD		 1414 if (SCM_INUMP (k))
1415 {
1416 long nk = SCM_INUM (k);
1417 return SCM_BOOL (nj & nk);
1418 }
1419 else if (SCM_BIGP (k))
1420 {
1421 intbig:
1422 if (nj == 0)
1423 return SCM_BOOL_F;
1424 {
1425 SCM result;
1426 mpz_t nj_z;
1427 mpz_init_set_si (nj_z, nj);
1428 mpz_and (nj_z, nj_z, SCM_I_BIG_MPZ (k));
1429 scm_remember_upto_here_1 (k);
1430 result = SCM_BOOL (mpz_sgn (nj_z) != 0);
1431 mpz_clear (nj_z);
1432 return result;
1433 }
1434 }
1435 else
1436 SCM_WRONG_TYPE_ARG (SCM_ARG2, k);
1437 }
1438 else if (SCM_BIGP (j))
1439 {
1440 if (SCM_INUMP (k))
1441 {
1442 SCM_SWAP (j, k);
1443 nj = SCM_INUM (j);
1444 goto intbig;
1445 }
1446 else if (SCM_BIGP (k))
1447 {
1448 SCM result;
1449 mpz_t result_z;
1450 mpz_init (result_z);
1451 mpz_and (result_z,
1452 SCM_I_BIG_MPZ (j),
1453 SCM_I_BIG_MPZ (k));
1454 scm_remember_upto_here_2 (j, k);
1455 result = SCM_BOOL (mpz_sgn (result_z) != 0);
1456 mpz_clear (result_z);
1457 return result;
1458 }
1459 else
1460 SCM_WRONG_TYPE_ARG (SCM_ARG2, k);
1461 }
1462 else
1463 SCM_WRONG_TYPE_ARG (SCM_ARG1, j);
0f2d19dd 1464}
1bbd0b84 1465#undef FUNC_NAME
0f2d19dd 1466
c1bfcf60 1467
a1ec6916 1468SCM_DEFINE (scm_logbit_p, "logbit?", 2, 0, 0,
2cd04b42 1469 (SCM index, SCM j),
1e6808ea 1470 "@lisp\n"
b380b885
MD		 1471 "(logbit? index j) @equiv{} (logtest (integer-expt 2 index) j)\n\n"
1472 "(logbit? 0 #b1101) @result{} #t\n"
1473 "(logbit? 1 #b1101) @result{} #f\n"
1474 "(logbit? 2 #b1101) @result{} #t\n"
1475 "(logbit? 3 #b1101) @result{} #t\n"
1476 "(logbit? 4 #b1101) @result{} #f\n"
1e6808ea 1477 "@end lisp")
1bbd0b84 1478#define FUNC_NAME s_scm_logbit_p
0f2d19dd 1479{
78166ad5
DH		 1480 unsigned long int iindex;
1481
1482 SCM_VALIDATE_INUM_MIN (SCM_ARG1, index, 0);
1483 iindex = (unsigned long int) SCM_INUM (index);
1484
0aacf84e 1485 if (SCM_INUMP (j))
78166ad5 1486 return SCM_BOOL ((1L << iindex) & SCM_INUM (j));
0aacf84e
MD		 1487 else if (SCM_BIGP (j))
1488 {
1489 int val = mpz_tstbit (SCM_I_BIG_MPZ (j), iindex);
1490 scm_remember_upto_here_1 (j);
1491 return SCM_BOOL (val);
1492 }
1493 else
78166ad5 1494 SCM_WRONG_TYPE_ARG (SCM_ARG2, j);
0f2d19dd 1495}
1bbd0b84 1496#undef FUNC_NAME
0f2d19dd 1497
78166ad5 1498
a1ec6916 1499SCM_DEFINE (scm_lognot, "lognot", 1, 0, 0,
1bbd0b84 1500 (SCM n),
4d814788 1501 "Return the integer which is the ones-complement of the integer\n"
1e6808ea
MG		 1502 "argument.\n"
1503 "\n"
b380b885
MD		 1504 "@lisp\n"
1505 "(number->string (lognot #b10000000) 2)\n"
1506 " @result{} \"-10000001\"\n"
1507 "(number->string (lognot #b0) 2)\n"
1508 " @result{} \"-1\"\n"
1e6808ea 1509 "@end lisp")
1bbd0b84 1510#define FUNC_NAME s_scm_lognot
0f2d19dd 1511{
f9811f9f
KR		 1512 if (SCM_INUMP (n)) {
1513 /* No overflow here, just need to toggle all the bits making up the inum.
1514 Enhancement: No need to strip the tag and add it back, could just xor
1515 a block of 1 bits, if that worked with the various debug versions of
1516 the SCM typedef. */
1517 return SCM_MAKINUM (~ SCM_INUM (n));
1518
1519 } else if (SCM_BIGP (n)) {
1520 SCM result = scm_i_mkbig ();
1521 mpz_com (SCM_I_BIG_MPZ (result), SCM_I_BIG_MPZ (n));
1522 scm_remember_upto_here_1 (n);
1523 return result;
1524
1525 } else {
1526 SCM_WRONG_TYPE_ARG (SCM_ARG1, n);
1527 }
0f2d19dd 1528}
1bbd0b84 1529#undef FUNC_NAME
0f2d19dd 1530
518b7508
KR		 1531/* returns 0 if IN is not an integer. OUT must already be
1532 initialized. */
1533static int
1534coerce_to_big (SCM in, mpz_t out)
1535{
1536 if (SCM_BIGP (in))
1537 mpz_set (out, SCM_I_BIG_MPZ (in));
1538 else if (SCM_INUMP (in))
1539 mpz_set_si (out, SCM_INUM (in));
1540 else
1541 return 0;
1542
1543 return 1;
1544}
1545
d885e204 1546SCM_DEFINE (scm_modulo_expt, "modulo-expt", 3, 0, 0,
518b7508
KR		 1547 (SCM n, SCM k, SCM m),
1548 "Return @var{n} raised to the integer exponent\n"
1549 "@var{k}, modulo @var{m}.\n"
1550 "\n"
1551 "@lisp\n"
1552 "(modulo-expt 2 3 5)\n"
1553 " @result{} 3\n"
1554 "@end lisp")
d885e204 1555#define FUNC_NAME s_scm_modulo_expt
518b7508
KR		 1556{
1557 mpz_t n_tmp;
1558 mpz_t k_tmp;
1559 mpz_t m_tmp;
1560
1561 /* There are two classes of error we might encounter --
1562 1) Math errors, which we'll report by calling scm_num_overflow,
1563 and
1564 2) wrong-type errors, which of course we'll report by calling
1565 SCM_WRONG_TYPE_ARG.
1566 We don't report those errors immediately, however; instead we do
1567 some cleanup first. These variables tell us which error (if
1568 any) we should report after cleaning up.
1569 */
1570 int report_overflow = 0;
1571
1572 int position_of_wrong_type = 0;
1573 SCM value_of_wrong_type = SCM_INUM0;
1574
1575 SCM result = SCM_UNDEFINED;
1576
1577 mpz_init (n_tmp);
1578 mpz_init (k_tmp);
1579 mpz_init (m_tmp);
1580
1581 if (SCM_EQ_P (m, SCM_INUM0))
1582 {
1583 report_overflow = 1;
1584 goto cleanup;
1585 }
1586
1587 if (!coerce_to_big (n, n_tmp))
1588 {
1589 value_of_wrong_type = n;
1590 position_of_wrong_type = 1;
1591 goto cleanup;
1592 }
1593
1594 if (!coerce_to_big (k, k_tmp))
1595 {
1596 value_of_wrong_type = k;
1597 position_of_wrong_type = 2;
1598 goto cleanup;
1599 }
1600
1601 if (!coerce_to_big (m, m_tmp))
1602 {
1603 value_of_wrong_type = m;
1604 position_of_wrong_type = 3;
1605 goto cleanup;
1606 }
1607
1608 /* if the exponent K is negative, and we simply call mpz_powm, we
1609 will get a divide-by-zero exception when an inverse 1/n mod m
1610 doesn't exist (or is not unique). Since exceptions are hard to
1611 handle, we'll attempt the inversion "by hand" -- that way, we get
1612 a simple failure code, which is easy to handle. */
1613
1614 if (-1 == mpz_sgn (k_tmp))
1615 {
1616 if (!mpz_invert (n_tmp, n_tmp, m_tmp))
1617 {
1618 report_overflow = 1;
1619 goto cleanup;
1620 }
1621 mpz_neg (k_tmp, k_tmp);
1622 }
1623
1624 result = scm_i_mkbig ();
1625 mpz_powm (SCM_I_BIG_MPZ (result),
1626 n_tmp,
1627 k_tmp,
1628 m_tmp);
b7b8c575
KR		 1629
1630 if (mpz_sgn (m_tmp) < 0 && mpz_sgn (SCM_I_BIG_MPZ (result)) != 0)
1631 mpz_add (SCM_I_BIG_MPZ (result), SCM_I_BIG_MPZ (result), m_tmp);
1632
518b7508
KR		 1633 cleanup:
1634 mpz_clear (m_tmp);
1635 mpz_clear (k_tmp);
1636 mpz_clear (n_tmp);
1637
1638 if (report_overflow)
1639 scm_num_overflow (FUNC_NAME);
1640
1641 if (position_of_wrong_type)
1642 SCM_WRONG_TYPE_ARG (position_of_wrong_type,
1643 value_of_wrong_type);
1644
1645 return scm_i_normbig (result);
1646}
1647#undef FUNC_NAME
1648
a1ec6916 1649SCM_DEFINE (scm_integer_expt, "integer-expt", 2, 0, 0,
2cd04b42 1650 (SCM n, SCM k),
1e6808ea
MG		 1651 "Return @var{n} raised to the non-negative integer exponent\n"
1652 "@var{k}.\n"
1653 "\n"
b380b885
MD		 1654 "@lisp\n"
1655 "(integer-expt 2 5)\n"
1656 " @result{} 32\n"
1657 "(integer-expt -3 3)\n"
1658 " @result{} -27\n"
1659 "@end lisp")
1bbd0b84 1660#define FUNC_NAME s_scm_integer_expt
0f2d19dd 1661{
1c35cb19
RB		 1662 long i2 = 0;
1663 SCM z_i2 = SCM_BOOL_F;
1664 int i2_is_big = 0;
f872b822 1665 SCM acc = SCM_MAKINUM (1L);
ca46fb90 1666
d57ed702 1667 /* 0^0 == 1 according to R5RS */
4260a7fc 1668 if (SCM_EQ_P (n, SCM_INUM0) || SCM_EQ_P (n, acc))
7b3381f4 1669 return SCM_FALSEP (scm_zero_p(k)) ? n : acc;
4260a7fc
DH		 1670 else if (SCM_EQ_P (n, SCM_MAKINUM (-1L)))
1671 return SCM_FALSEP (scm_even_p (k)) ? n : acc;
ca46fb90 1672
ca46fb90
RB		 1673 if (SCM_INUMP (k))
1674 i2 = SCM_INUM (k);
1675 else if (SCM_BIGP (k))
1676 {
1677 z_i2 = scm_i_clonebig (k, 1);
ca46fb90
RB		 1678 scm_remember_upto_here_1 (k);
1679 i2_is_big = 1;
1680 }
1681 else if (SCM_REALP (k))
2830fd91
MD		 1682 {
1683 double r = SCM_REAL_VALUE (k);
ca46fb90
RB		 1684 if (floor (r) != r)
1685 SCM_WRONG_TYPE_ARG (2, k);
1686 if ((r > SCM_MOST_POSITIVE_FIXNUM) || (r < SCM_MOST_NEGATIVE_FIXNUM))
1687 {
1688 z_i2 = scm_i_mkbig ();
753ac1e7 1689 mpz_set_d (SCM_I_BIG_MPZ (z_i2), r);
ca46fb90
RB		 1690 i2_is_big = 1;
1691 }
1692 else
1693 {
1694 i2 = r;
1695 }
2830fd91
MD		 1696 }
1697 else
ca46fb90
RB		 1698 SCM_WRONG_TYPE_ARG (2, k);
1699
1700 if (i2_is_big)
f872b822 1701 {
ca46fb90
RB		 1702 if (mpz_sgn(SCM_I_BIG_MPZ (z_i2)) == -1)
1703 {
1704 mpz_neg (SCM_I_BIG_MPZ (z_i2), SCM_I_BIG_MPZ (z_i2));
1705 n = scm_divide (n, SCM_UNDEFINED);
1706 }
1707 while (1)
1708 {
1709 if (mpz_sgn(SCM_I_BIG_MPZ (z_i2)) == 0)
1710 {
ca46fb90
RB		 1711 return acc;
1712 }
1713 if (mpz_cmp_ui(SCM_I_BIG_MPZ (z_i2), 1) == 0)
1714 {
ca46fb90
RB		 1715 return scm_product (acc, n);
1716 }
1717 if (mpz_tstbit(SCM_I_BIG_MPZ (z_i2), 0))
1718 acc = scm_product (acc, n);
1719 n = scm_product (n, n);
1720 mpz_fdiv_q_2exp (SCM_I_BIG_MPZ (z_i2), SCM_I_BIG_MPZ (z_i2), 1);
1721 }
f872b822 1722 }
ca46fb90 1723 else
f872b822 1724 {
ca46fb90
RB		 1725 if (i2 < 0)
1726 {
1727 i2 = -i2;
1728 n = scm_divide (n, SCM_UNDEFINED);
1729 }
1730 while (1)
1731 {
1732 if (0 == i2)
1733 return acc;
1734 if (1 == i2)
1735 return scm_product (acc, n);
1736 if (i2 & 1)
1737 acc = scm_product (acc, n);
1738 n = scm_product (n, n);
1739 i2 >>= 1;
1740 }
f872b822 1741 }
0f2d19dd 1742}
1bbd0b84 1743#undef FUNC_NAME
0f2d19dd 1744
a1ec6916 1745SCM_DEFINE (scm_ash, "ash", 2, 0, 0,
1bbd0b84 1746 (SCM n, SCM cnt),
32f19569
KR		 1747 "Return @var{n} shifted left by @var{cnt} bits, or shifted right\n"
1748 "if @var{cnt} is negative. This is an ``arithmetic'' shift.\n"
1e6808ea 1749 "\n"
e7644cb2 1750 "This is effectively a multiplication by 2^@var{cnt}, and when\n"
32f19569
KR		 1751 "@var{cnt} is negative it's a division, rounded towards negative\n"
1752 "infinity. (Note that this is not the same rounding as\n"
1753 "@code{quotient} does.)\n"
1754 "\n"
1755 "With @var{n} viewed as an infinite precision twos complement,\n"
1756 "@code{ash} means a left shift introducing zero bits, or a right\n"
1757 "shift dropping bits.\n"
1e6808ea 1758 "\n"
b380b885 1759 "@lisp\n"
1e6808ea
MG		 1760 "(number->string (ash #b1 3) 2) @result{} \"1000\"\n"
1761 "(number->string (ash #b1010 -1) 2) @result{} \"101\"\n"
32f19569
KR		 1762 "\n"
1763 ";; -23 is bits ...11101001, -6 is bits ...111010\n"
1764 "(ash -23 -2) @result{} -6\n"
a3c8b9fc 1765 "@end lisp")
1bbd0b84 1766#define FUNC_NAME s_scm_ash
0f2d19dd 1767{
3ab9f56e
DH		 1768 long bits_to_shift;
1769
3ab9f56e
DH		 1770 SCM_VALIDATE_INUM (2, cnt);
1771
1772 bits_to_shift = SCM_INUM (cnt);
ca46fb90
RB		 1773
1774 if (bits_to_shift < 0)
1775 {
1776 /* Shift right by abs(cnt) bits. This is realized as a division
1777 by div:=2^abs(cnt). However, to guarantee the floor
1778 rounding, negative values require some special treatment.
1779 */
1780 SCM div = scm_integer_expt (SCM_MAKINUM (2),
1781 SCM_MAKINUM (-bits_to_shift));
f92e85f7
MV		 1782
1783 /* scm_quotient assumes its arguments are integers, but it's legal to (ash 1/2 -1) */
ca46fb90
RB		 1784 if (SCM_FALSEP (scm_negative_p (n)))
1785 return scm_quotient (n, div);
1786 else
1787 return scm_sum (SCM_MAKINUM (-1L),
1788 scm_quotient (scm_sum (SCM_MAKINUM (1L), n), div));
1789 }
1790 else
3ab9f56e 1791 /* Shift left is done by multiplication with 2^CNT */
f872b822 1792 return scm_product (n, scm_integer_expt (SCM_MAKINUM (2), cnt));
0f2d19dd 1793}
1bbd0b84 1794#undef FUNC_NAME
0f2d19dd 1795
3c9f20f8 1796
a1ec6916 1797SCM_DEFINE (scm_bit_extract, "bit-extract", 3, 0, 0,
1bbd0b84 1798 (SCM n, SCM start, SCM end),
1e6808ea
MG		 1799 "Return the integer composed of the @var{start} (inclusive)\n"
1800 "through @var{end} (exclusive) bits of @var{n}. The\n"
1801 "@var{start}th bit becomes the 0-th bit in the result.\n"
1802 "\n"
b380b885
MD		 1803 "@lisp\n"
1804 "(number->string (bit-extract #b1101101010 0 4) 2)\n"
1805 " @result{} \"1010\"\n"
1806 "(number->string (bit-extract #b1101101010 4 9) 2)\n"
1807 " @result{} \"10110\"\n"
1808 "@end lisp")
1bbd0b84 1809#define FUNC_NAME s_scm_bit_extract
0f2d19dd 1810{
7f848242 1811 unsigned long int istart, iend, bits;
34d19ef6 1812 SCM_VALIDATE_INUM_MIN_COPY (2, start,0, istart);
c1bfcf60
GB		 1813 SCM_VALIDATE_INUM_MIN_COPY (3, end, 0, iend);
1814 SCM_ASSERT_RANGE (3, end, (iend >= istart));
78166ad5 1815
7f848242
KR		 1816 /* how many bits to keep */
1817 bits = iend - istart;
1818
0aacf84e
MD		 1819 if (SCM_INUMP (n))
1820 {
1821 long int in = SCM_INUM (n);
7f848242
KR		 1822
1823 /* When istart>=SCM_I_FIXNUM_BIT we can just limit the shift to
d77ad560 1824 SCM_I_FIXNUM_BIT-1 to get either 0 or -1 per the sign of "in". */
857ae6af 1825 in = SCM_SRS (in, min (istart, SCM_I_FIXNUM_BIT-1));
ac0c002c 1826
0aacf84e
MD		 1827 if (in < 0 && bits >= SCM_I_FIXNUM_BIT)
1828 {
1829 /* Since we emulate two's complement encoded numbers, this
1830 * special case requires us to produce a result that has
7f848242 1831 * more bits than can be stored in a fixnum.
0aacf84e 1832 */
7f848242
KR		 1833 SCM result = scm_i_long2big (in);
1834 mpz_fdiv_r_2exp (SCM_I_BIG_MPZ (result), SCM_I_BIG_MPZ (result),
1835 bits);
1836 return result;
0aacf84e 1837 }
ac0c002c 1838
7f848242 1839 /* mask down to requisite bits */
857ae6af 1840 bits = min (bits, SCM_I_FIXNUM_BIT);
7f848242 1841 return SCM_MAKINUM (in & ((1L << bits) - 1));
0aacf84e
MD		 1842 }
1843 else if (SCM_BIGP (n))
ac0c002c 1844 {
7f848242
KR		 1845 SCM result;
1846 if (bits == 1)
1847 {
1848 result = SCM_MAKINUM (mpz_tstbit (SCM_I_BIG_MPZ (n), istart));
1849 }
1850 else
1851 {
1852 /* ENHANCE-ME: It'd be nice not to allocate a new bignum when
1853 bits<SCM_I_FIXNUM_BIT. Would want some help from GMP to get
1854 such bits into a ulong. */
1855 result = scm_i_mkbig ();
1856 mpz_fdiv_q_2exp (SCM_I_BIG_MPZ(result), SCM_I_BIG_MPZ(n), istart);
1857 mpz_fdiv_r_2exp (SCM_I_BIG_MPZ(result), SCM_I_BIG_MPZ(result), bits);
1858 result = scm_i_normbig (result);
1859 }
1860 scm_remember_upto_here_1 (n);
1861 return result;
ac0c002c 1862 }
0aacf84e 1863 else
78166ad5 1864 SCM_WRONG_TYPE_ARG (SCM_ARG1, n);
0f2d19dd 1865}
1bbd0b84 1866#undef FUNC_NAME
0f2d19dd 1867
7f848242 1868
e4755e5c
JB		 1869static const char scm_logtab[] = {
1870 0, 1, 1, 2, 1, 2, 2, 3, 1, 2, 2, 3, 2, 3, 3, 4
1871};
1cc91f1b 1872
a1ec6916 1873SCM_DEFINE (scm_logcount, "logcount", 1, 0, 0,
1bbd0b84 1874 (SCM n),
1e6808ea
MG		 1875 "Return the number of bits in integer @var{n}. If integer is\n"
1876 "positive, the 1-bits in its binary representation are counted.\n"
1877 "If negative, the 0-bits in its two's-complement binary\n"
1878 "representation are counted. If 0, 0 is returned.\n"
1879 "\n"
b380b885
MD		 1880 "@lisp\n"
1881 "(logcount #b10101010)\n"
ca46fb90
RB		 1882 " @result{} 4\n"
1883 "(logcount 0)\n"
1884 " @result{} 0\n"
1885 "(logcount -2)\n"
1886 " @result{} 1\n"
1887 "@end lisp")
1888#define FUNC_NAME s_scm_logcount
1889{
1890 if (SCM_INUMP (n))
f872b822 1891 {
ca46fb90
RB		 1892 unsigned long int c = 0;
1893 long int nn = SCM_INUM (n);
1894 if (nn < 0)
1895 nn = -1 - nn;
1896 while (nn)
1897 {
1898 c += scm_logtab[15 & nn];
1899 nn >>= 4;
1900 }
1901 return SCM_MAKINUM (c);
f872b822 1902 }
ca46fb90 1903 else if (SCM_BIGP (n))
f872b822 1904 {
ca46fb90 1905 unsigned long count;
713a4259
KR		 1906 if (mpz_sgn (SCM_I_BIG_MPZ (n)) >= 0)
1907 count = mpz_popcount (SCM_I_BIG_MPZ (n));
ca46fb90 1908 else
713a4259
KR		 1909 count = mpz_hamdist (SCM_I_BIG_MPZ (n), z_negative_one);
1910 scm_remember_upto_here_1 (n);
ca46fb90 1911 return SCM_MAKINUM (count);
f872b822 1912 }
ca46fb90
RB		 1913 else
1914 SCM_WRONG_TYPE_ARG (SCM_ARG1, n);
0f2d19dd 1915}
ca46fb90 1916#undef FUNC_NAME
0f2d19dd
JB		 1917
1918
ca46fb90
RB		 1919static const char scm_ilentab[] = {
1920 0, 1, 2, 2, 3, 3, 3, 3, 4, 4, 4, 4, 4, 4, 4, 4
1921};
1922
0f2d19dd 1923
ca46fb90
RB		 1924SCM_DEFINE (scm_integer_length, "integer-length", 1, 0, 0,
1925 (SCM n),
1926 "Return the number of bits necessary to represent @var{n}.\n"
1927 "\n"
1928 "@lisp\n"
1929 "(integer-length #b10101010)\n"
1930 " @result{} 8\n"
1931 "(integer-length 0)\n"
1932 " @result{} 0\n"
1933 "(integer-length #b1111)\n"
1934 " @result{} 4\n"
1935 "@end lisp")
1936#define FUNC_NAME s_scm_integer_length
1937{
0aacf84e
MD		 1938 if (SCM_INUMP (n))
1939 {
1940 unsigned long int c = 0;
1941 unsigned int l = 4;
1942 long int nn = SCM_INUM (n);
1943 if (nn < 0)
1944 nn = -1 - nn;
1945 while (nn)
1946 {
1947 c += 4;
1948 l = scm_ilentab [15 & nn];
1949 nn >>= 4;
1950 }
1951 return SCM_MAKINUM (c - 4 + l);
1952 }
1953 else if (SCM_BIGP (n))
1954 {
1955 /* mpz_sizeinbase looks at the absolute value of negatives, whereas we
1956 want a ones-complement. If n is ...111100..00 then mpz_sizeinbase is
1957 1 too big, so check for that and adjust. */
1958 size_t size = mpz_sizeinbase (SCM_I_BIG_MPZ (n), 2);
1959 if (mpz_sgn (SCM_I_BIG_MPZ (n)) < 0
1960 && mpz_scan0 (SCM_I_BIG_MPZ (n), /* no 0 bits above the lowest 1 */
1961 mpz_scan1 (SCM_I_BIG_MPZ (n), 0)) == ULONG_MAX)
1962 size--;
1963 scm_remember_upto_here_1 (n);
1964 return SCM_MAKINUM (size);
1965 }
1966 else
ca46fb90 1967 SCM_WRONG_TYPE_ARG (SCM_ARG1, n);
ca46fb90
RB		 1968}
1969#undef FUNC_NAME
0f2d19dd
JB		 1970
1971/*** NUMBERS -> STRINGS ***/
0f2d19dd 1972int scm_dblprec;
e4755e5c 1973static const double fx[] =
f872b822
MD		 1974{ 0.0, 5e-1, 5e-2, 5e-3, 5e-4, 5e-5,
1975 5e-6, 5e-7, 5e-8, 5e-9, 5e-10,
1976 5e-11, 5e-12, 5e-13, 5e-14, 5e-15,
1977 5e-16, 5e-17, 5e-18, 5e-19, 5e-20};
0f2d19dd 1978
1be6b49c 1979static size_t
1bbd0b84 1980idbl2str (double f, char *a)
0f2d19dd
JB		 1981{
1982 int efmt, dpt, d, i, wp = scm_dblprec;
1be6b49c 1983 size_t ch = 0;
0f2d19dd
JB		 1984 int exp = 0;
1985
f872b822 1986 if (f == 0.0)
abb7e44d
MV		 1987 {
1988#ifdef HAVE_COPYSIGN
1989 double sgn = copysign (1.0, f);
1990
1991 if (sgn < 0.0)
1992 a[ch++] = '-';
1993#endif
1994
1995 goto zero; /*{a[0]='0'; a[1]='.'; a[2]='0'; return 3;} */
1996 }
7351e207
MV		 1997
1998 if (xisinf (f))
1999 {
2000 if (f < 0)
2001 strcpy (a, "-inf.0");
2002 else
2003 strcpy (a, "+inf.0");
2004 return ch+6;
2005 }
2006 else if (xisnan (f))
2007 {
2008 strcpy (a, "+nan.0");
2009 return ch+6;
2010 }
2011
f872b822
MD		 2012 if (f < 0.0)
2013 {
2014 f = -f;
2015 a[ch++] = '-';
2016 }
7351e207 2017
f872b822
MD		 2018#ifdef DBL_MIN_10_EXP /* Prevent unnormalized values, as from
2019 make-uniform-vector, from causing infinite loops. */
2020 while (f < 1.0)
2021 {
2022 f *= 10.0;
2023 if (exp-- < DBL_MIN_10_EXP)
7351e207
MV		 2024 {
2025 a[ch++] = '#';
2026 a[ch++] = '.';
2027 a[ch++] = '#';
2028 return ch;
2029 }
f872b822
MD		 2030 }
2031 while (f > 10.0)
2032 {
2033 f *= 0.10;
2034 if (exp++ > DBL_MAX_10_EXP)
7351e207
MV		 2035 {
2036 a[ch++] = '#';
2037 a[ch++] = '.';
2038 a[ch++] = '#';
2039 return ch;
2040 }
f872b822
MD		 2041 }
2042#else
2043 while (f < 1.0)
2044 {
2045 f *= 10.0;
2046 exp--;
2047 }
2048 while (f > 10.0)
2049 {
2050 f /= 10.0;
2051 exp++;
2052 }
2053#endif
2054 if (f + fx[wp] >= 10.0)
2055 {
2056 f = 1.0;
2057 exp++;
2058 }
0f2d19dd 2059 zero:
f872b822
MD		 2060#ifdef ENGNOT
2061 dpt = (exp + 9999) % 3;
0f2d19dd
JB		 2062 exp -= dpt++;
2063 efmt = 1;
f872b822
MD		 2064#else
2065 efmt = (exp < -3) || (exp > wp + 2);
0f2d19dd 2066 if (!efmt)
cda139a7
MD		 2067 {
2068 if (exp < 0)
2069 {
2070 a[ch++] = '0';
2071 a[ch++] = '.';
2072 dpt = exp;
f872b822
MD		 2073 while (++dpt)
2074 a[ch++] = '0';
cda139a7
MD		 2075 }
2076 else
f872b822 2077 dpt = exp + 1;
cda139a7 2078 }
0f2d19dd
JB		 2079 else
2080 dpt = 1;
f872b822
MD		 2081#endif
2082
2083 do
2084 {
2085 d = f;
2086 f -= d;
2087 a[ch++] = d + '0';
2088 if (f < fx[wp])
2089 break;
2090 if (f + fx[wp] >= 1.0)
2091 {
2092 a[ch - 1]++;
2093 break;
2094 }
2095 f *= 10.0;
2096 if (!(--dpt))
2097 a[ch++] = '.';
0f2d19dd 2098 }
f872b822 2099 while (wp--);
0f2d19dd
JB		 2100
2101 if (dpt > 0)
cda139a7 2102 {
f872b822 2103#ifndef ENGNOT
cda139a7
MD		 2104 if ((dpt > 4) && (exp > 6))
2105 {
f872b822 2106 d = (a[0] == '-' ? 2 : 1);
cda139a7 2107 for (i = ch++; i > d; i--)
f872b822 2108 a[i] = a[i - 1];
cda139a7
MD		 2109 a[d] = '.';
2110 efmt = 1;
2111 }
2112 else
f872b822 2113#endif
cda139a7 2114 {
f872b822
MD		 2115 while (--dpt)
2116 a[ch++] = '0';
cda139a7
MD		 2117 a[ch++] = '.';
2118 }
2119 }
f872b822
MD		 2120 if (a[ch - 1] == '.')
2121 a[ch++] = '0'; /* trailing zero */
2122 if (efmt && exp)
2123 {
2124 a[ch++] = 'e';
2125 if (exp < 0)
2126 {
2127 exp = -exp;
2128 a[ch++] = '-';
2129 }
2130 for (i = 10; i <= exp; i *= 10);
2131 for (i /= 10; i; i /= 10)
2132 {
2133 a[ch++] = exp / i + '0';
2134 exp %= i;
2135 }
0f2d19dd 2136 }
0f2d19dd
JB		 2137 return ch;
2138}
2139
1cc91f1b 2140
1be6b49c 2141static size_t
1bbd0b84 2142iflo2str (SCM flt, char *str)
0f2d19dd 2143{
1be6b49c 2144 size_t i;
3c9a524f 2145 if (SCM_REALP (flt))
f3ae5d60 2146 i = idbl2str (SCM_REAL_VALUE (flt), str);
0f2d19dd 2147 else
f872b822 2148 {
f3ae5d60
MD		 2149 i = idbl2str (SCM_COMPLEX_REAL (flt), str);
2150 if (SCM_COMPLEX_IMAG (flt) != 0.0)
2151 {
7351e207
MV		 2152 double imag = SCM_COMPLEX_IMAG (flt);
2153 /* Don't output a '+' for negative numbers or for Inf and
2154 NaN. They will provide their own sign. */
2155 if (0 <= imag && !xisinf (imag) && !xisnan (imag))
f3ae5d60 2156 str[i++] = '+';
7351e207 2157 i += idbl2str (imag, &str[i]);
f3ae5d60
MD		 2158 str[i++] = 'i';
2159 }
f872b822 2160 }
0f2d19dd
JB		 2161 return i;
2162}
0f2d19dd 2163
5c11cc9d 2164/* convert a long to a string (unterminated). returns the number of
1bbd0b84
GB		 2165 characters in the result.
2166 rad is output base
2167 p is destination: worst case (base 2) is SCM_INTBUFLEN */
1be6b49c 2168size_t
1bbd0b84 2169scm_iint2str (long num, int rad, char *p)
0f2d19dd 2170{
1be6b49c
ML		 2171 size_t j = 1;
2172 size_t i;
5c11cc9d
GH		 2173 unsigned long n = (num < 0) ? -num : num;
2174
f872b822 2175 for (n /= rad; n > 0; n /= rad)
5c11cc9d
GH		 2176 j++;
2177
2178 i = j;
2179 if (num < 0)
f872b822 2180 {
f872b822 2181 *p++ = '-';
5c11cc9d
GH		 2182 j++;
2183 n = -num;
f872b822 2184 }
5c11cc9d
GH		 2185 else
2186 n = num;
f872b822
MD		 2187 while (i--)
2188 {
5c11cc9d
GH		 2189 int d = n % rad;
2190
f872b822
MD		 2191 n /= rad;
2192 p[i] = d + ((d < 10) ? '0' : 'a' - 10);
2193 }
0f2d19dd
JB		 2194 return j;
2195}
2196
a1ec6916 2197SCM_DEFINE (scm_number_to_string, "number->string", 1, 1, 0,
bb628794
DH		 2198 (SCM n, SCM radix),
2199 "Return a string holding the external representation of the\n"
942e5b91
MG		 2200 "number @var{n} in the given @var{radix}. If @var{n} is\n"
2201 "inexact, a radix of 10 will be used.")
1bbd0b84 2202#define FUNC_NAME s_scm_number_to_string
0f2d19dd 2203{
1bbd0b84 2204 int base;
98cb6e75 2205
0aacf84e 2206 if (SCM_UNBNDP (radix))
98cb6e75 2207 base = 10;
0aacf84e
MD		 2208 else
2209 {
2210 SCM_VALIDATE_INUM (2, radix);
2211 base = SCM_INUM (radix);
2212 /* FIXME: ask if range limit was OK, and if so, document */
2213 SCM_ASSERT_RANGE (2, radix, (base >= 2) && (base <= 36));
2214 }
98cb6e75 2215
0aacf84e
MD		 2216 if (SCM_INUMP (n))
2217 {
2218 char num_buf [SCM_INTBUFLEN];
2219 size_t length = scm_iint2str (SCM_INUM (n), base, num_buf);
2220 return scm_mem2string (num_buf, length);
2221 }
2222 else if (SCM_BIGP (n))
2223 {
2224 char *str = mpz_get_str (NULL, base, SCM_I_BIG_MPZ (n));
2225 scm_remember_upto_here_1 (n);
2226 return scm_take0str (str);
2227 }
f92e85f7
MV		 2228 else if (SCM_FRACTIONP (n))
2229 {
2230 scm_i_fraction_reduce (n);
2231 return scm_string_append (scm_list_3 (scm_number_to_string (SCM_FRACTION_NUMERATOR (n), radix),
2232 scm_mem2string ("/", 1),
2233 scm_number_to_string (SCM_FRACTION_DENOMINATOR (n), radix)));
2234 }
0aacf84e
MD		 2235 else if (SCM_INEXACTP (n))
2236 {
2237 char num_buf [FLOBUFLEN];
2238 return scm_mem2string (num_buf, iflo2str (n, num_buf));
2239 }
2240 else
bb628794 2241 SCM_WRONG_TYPE_ARG (1, n);
0f2d19dd 2242}
1bbd0b84 2243#undef FUNC_NAME
0f2d19dd
JB		 2244
2245
ca46fb90
RB		 2246/* These print routines used to be stubbed here so that scm_repl.c
2247 wouldn't need SCM_BIGDIG conditionals (pre GMP) */
1cc91f1b 2248
0f2d19dd 2249int
e81d98ec 2250scm_print_real (SCM sexp, SCM port, scm_print_state *pstate SCM_UNUSED)
0f2d19dd 2251{
56e55ac7 2252 char num_buf[FLOBUFLEN];
f872b822 2253 scm_lfwrite (num_buf, iflo2str (sexp, num_buf), port);
0f2d19dd
JB		 2254 return !0;
2255}
2256
f3ae5d60 2257int
e81d98ec 2258scm_print_complex (SCM sexp, SCM port, scm_print_state *pstate SCM_UNUSED)
f92e85f7 2259
f3ae5d60 2260{
56e55ac7 2261 char num_buf[FLOBUFLEN];
f3ae5d60
MD		 2262 scm_lfwrite (num_buf, iflo2str (sexp, num_buf), port);
2263 return !0;
2264}
1cc91f1b 2265
f92e85f7
MV		 2266int
2267scm_i_print_fraction (SCM sexp, SCM port, scm_print_state *pstate SCM_UNUSED)
2268{
2269 SCM str;
2270 scm_i_fraction_reduce (sexp);
2271 str = scm_number_to_string (sexp, SCM_UNDEFINED);
2272 scm_lfwrite (SCM_STRING_CHARS (str), SCM_STRING_LENGTH (str), port);
2273 scm_remember_upto_here_1 (str);
2274 return !0;
2275}
2276
0f2d19dd 2277int
e81d98ec 2278scm_bigprint (SCM exp, SCM port, scm_print_state *pstate SCM_UNUSED)
0f2d19dd 2279{
ca46fb90
RB		 2280 char *str = mpz_get_str (NULL, 10, SCM_I_BIG_MPZ (exp));
2281 scm_remember_upto_here_1 (exp);
2282 scm_lfwrite (str, (size_t) strlen (str), port);
2283 free (str);
0f2d19dd
JB		 2284 return !0;
2285}
2286/*** END nums->strs ***/
2287
3c9a524f 2288
0f2d19dd 2289/*** STRINGS -> NUMBERS ***/
2a8fecee 2290
3c9a524f
DH		 2291/* The following functions implement the conversion from strings to numbers.
2292 * The implementation somehow follows the grammar for numbers as it is given
2293 * in R5RS. Thus, the functions resemble syntactic units (<ureal R>,
2294 * <uinteger R>, ...) that are used to build up numbers in the grammar. Some
2295 * points should be noted about the implementation:
2296 * * Each function keeps a local index variable 'idx' that points at the
2297 * current position within the parsed string. The global index is only
2298 * updated if the function could parse the corresponding syntactic unit
2299 * successfully.
2300 * * Similarly, the functions keep track of indicators of inexactness ('#',
2301 * '.' or exponents) using local variables ('hash_seen', 'x'). Again, the
2302 * global exactness information is only updated after each part has been
2303 * successfully parsed.
2304 * * Sequences of digits are parsed into temporary variables holding fixnums.
2305 * Only if these fixnums would overflow, the result variables are updated
2306 * using the standard functions scm_add, scm_product, scm_divide etc. Then,
2307 * the temporary variables holding the fixnums are cleared, and the process
2308 * starts over again. If for example fixnums were able to store five decimal
2309 * digits, a number 1234567890 would be parsed in two parts 12345 and 67890,
2310 * and the result was computed as 12345 * 100000 + 67890. In other words,
2311 * only every five digits two bignum operations were performed.
2312 */
2313
2314enum t_exactness {NO_EXACTNESS, INEXACT, EXACT};
2315
2316/* R5RS, section 7.1.1, lexical structure of numbers: <uinteger R>. */
2317
2318/* In non ASCII-style encodings the following macro might not work. */
71df73ac
KR		 2319#define XDIGIT2UINT(d) \
2320 (isdigit ((int) (unsigned char) d) \
2321 ? (d) - '0' \
2322 : tolower ((int) (unsigned char) d) - 'a' + 10)
3c9a524f 2323
2a8fecee 2324static SCM
3c9a524f
DH		 2325mem2uinteger (const char* mem, size_t len, unsigned int *p_idx,
2326 unsigned int radix, enum t_exactness *p_exactness)
2a8fecee 2327{
3c9a524f
DH		 2328 unsigned int idx = *p_idx;
2329 unsigned int hash_seen = 0;
2330 scm_t_bits shift = 1;
2331 scm_t_bits add = 0;
2332 unsigned int digit_value;
2333 SCM result;
2334 char c;
2335
2336 if (idx == len)
2337 return SCM_BOOL_F;
2a8fecee 2338
3c9a524f 2339 c = mem[idx];
71df73ac 2340 if (!isxdigit ((int) (unsigned char) c))
3c9a524f
DH		 2341 return SCM_BOOL_F;
2342 digit_value = XDIGIT2UINT (c);
2343 if (digit_value >= radix)
2344 return SCM_BOOL_F;
2345
2346 idx++;
2347 result = SCM_MAKINUM (digit_value);
2348 while (idx != len)
f872b822 2349 {
3c9a524f 2350 char c = mem[idx];
71df73ac 2351 if (isxdigit ((int) (unsigned char) c))
f872b822 2352 {
3c9a524f 2353 if (hash_seen)
1fe5e088 2354 break;
3c9a524f
DH		 2355 digit_value = XDIGIT2UINT (c);
2356 if (digit_value >= radix)
1fe5e088 2357 break;
f872b822 2358 }
3c9a524f
DH		 2359 else if (c == '#')
2360 {
2361 hash_seen = 1;
2362 digit_value = 0;
2363 }
2364 else
2365 break;
2366
2367 idx++;
2368 if (SCM_MOST_POSITIVE_FIXNUM / radix < shift)
2369 {
2370 result = scm_product (result, SCM_MAKINUM (shift));
2371 if (add > 0)
2372 result = scm_sum (result, SCM_MAKINUM (add));
2373
2374 shift = radix;
2375 add = digit_value;
2376 }
2377 else
2378 {
2379 shift = shift * radix;
2380 add = add * radix + digit_value;
2381 }
2382 };
2383
2384 if (shift > 1)
2385 result = scm_product (result, SCM_MAKINUM (shift));
2386 if (add > 0)
2387 result = scm_sum (result, SCM_MAKINUM (add));
2388
2389 *p_idx = idx;
2390 if (hash_seen)
2391 *p_exactness = INEXACT;
2392
2393 return result;
2a8fecee
JB		 2394}
2395
2396
3c9a524f
DH		 2397/* R5RS, section 7.1.1, lexical structure of numbers: <decimal 10>. Only
2398 * covers the parts of the rules that start at a potential point. The value
2399 * of the digits up to the point have been parsed by the caller and are given
79d34f68
DH		 2400 * in variable result. The content of *p_exactness indicates, whether a hash
2401 * has already been seen in the digits before the point.
3c9a524f 2402 */
1cc91f1b 2403
3c9a524f
DH		 2404/* In non ASCII-style encodings the following macro might not work. */
2405#define DIGIT2UINT(d) ((d) - '0')
2406
2407static SCM
79d34f68 2408mem2decimal_from_point (SCM result, const char* mem, size_t len,
3c9a524f 2409 unsigned int *p_idx, enum t_exactness *p_exactness)
0f2d19dd 2410{
3c9a524f
DH		 2411 unsigned int idx = *p_idx;
2412 enum t_exactness x = *p_exactness;
3c9a524f
DH		 2413
2414 if (idx == len)
79d34f68 2415 return result;
3c9a524f
DH		 2416
2417 if (mem[idx] == '.')
2418 {
2419 scm_t_bits shift = 1;
2420 scm_t_bits add = 0;
2421 unsigned int digit_value;
79d34f68 2422 SCM big_shift = SCM_MAKINUM (1);
3c9a524f
DH		 2423
2424 idx++;
2425 while (idx != len)
2426 {
2427 char c = mem[idx];
71df73ac 2428 if (isdigit ((int) (unsigned char) c))
3c9a524f
DH		 2429 {
2430 if (x == INEXACT)
2431 return SCM_BOOL_F;
2432 else
2433 digit_value = DIGIT2UINT (c);
2434 }
2435 else if (c == '#')
2436 {
2437 x = INEXACT;
2438 digit_value = 0;
2439 }
2440 else
2441 break;
2442
2443 idx++;
2444 if (SCM_MOST_POSITIVE_FIXNUM / 10 < shift)
2445 {
2446 big_shift = scm_product (big_shift, SCM_MAKINUM (shift));
79d34f68 2447 result = scm_product (result, SCM_MAKINUM (shift));
3c9a524f 2448 if (add > 0)
79d34f68 2449 result = scm_sum (result, SCM_MAKINUM (add));
3c9a524f
DH		 2450
2451 shift = 10;
2452 add = digit_value;
2453 }
2454 else
2455 {
2456 shift = shift * 10;
2457 add = add * 10 + digit_value;
2458 }
2459 };
2460
2461 if (add > 0)
2462 {
2463 big_shift = scm_product (big_shift, SCM_MAKINUM (shift));
79d34f68
DH		 2464 result = scm_product (result, SCM_MAKINUM (shift));
2465 result = scm_sum (result, SCM_MAKINUM (add));
3c9a524f
DH		 2466 }
2467
d8592269 2468 result = scm_divide (result, big_shift);
79d34f68 2469
3c9a524f
DH		 2470 /* We've seen a decimal point, thus the value is implicitly inexact. */
2471 x = INEXACT;
f872b822 2472 }
3c9a524f 2473
3c9a524f 2474 if (idx != len)
f872b822 2475 {
3c9a524f
DH		 2476 int sign = 1;
2477 unsigned int start;
2478 char c;
2479 int exponent;
2480 SCM e;
2481
2482 /* R5RS, section 7.1.1, lexical structure of numbers: <suffix> */
2483
2484 switch (mem[idx])
f872b822 2485 {
3c9a524f
DH		 2486 case 'd': case 'D':
2487 case 'e': case 'E':
2488 case 'f': case 'F':
2489 case 'l': case 'L':
2490 case 's': case 'S':
2491 idx++;
2492 start = idx;
2493 c = mem[idx];
2494 if (c == '-')
2495 {
2496 idx++;
2497 sign = -1;
2498 c = mem[idx];
2499 }
2500 else if (c == '+')
2501 {
2502 idx++;
2503 sign = 1;
2504 c = mem[idx];
2505 }
2506 else
2507 sign = 1;
2508
71df73ac 2509 if (!isdigit ((int) (unsigned char) c))
3c9a524f
DH		 2510 return SCM_BOOL_F;
2511
2512 idx++;
2513 exponent = DIGIT2UINT (c);
2514 while (idx != len)
f872b822 2515 {
3c9a524f 2516 char c = mem[idx];
71df73ac 2517 if (isdigit ((int) (unsigned char) c))
3c9a524f
DH		 2518 {
2519 idx++;
2520 if (exponent <= SCM_MAXEXP)
2521 exponent = exponent * 10 + DIGIT2UINT (c);
2522 }
2523 else
2524 break;
f872b822 2525 }
3c9a524f
DH		 2526
2527 if (exponent > SCM_MAXEXP)
f872b822 2528 {
3c9a524f
DH		 2529 size_t exp_len = idx - start;
2530 SCM exp_string = scm_mem2string (&mem[start], exp_len);
2531 SCM exp_num = scm_string_to_number (exp_string, SCM_UNDEFINED);
2532 scm_out_of_range ("string->number", exp_num);
f872b822 2533 }
3c9a524f
DH		 2534
2535 e = scm_integer_expt (SCM_MAKINUM (10), SCM_MAKINUM (exponent));
2536 if (sign == 1)
2537 result = scm_product (result, e);
2538 else
f92e85f7 2539 result = scm_divide2real (result, e);
3c9a524f
DH		 2540
2541 /* We've seen an exponent, thus the value is implicitly inexact. */
2542 x = INEXACT;
2543
f872b822 2544 break;
3c9a524f 2545
f872b822 2546 default:
3c9a524f 2547 break;
f872b822 2548 }
0f2d19dd 2549 }
3c9a524f
DH		 2550
2551 *p_idx = idx;
2552 if (x == INEXACT)
2553 *p_exactness = x;
2554
2555 return result;
0f2d19dd 2556}
0f2d19dd 2557
3c9a524f
DH		 2558
2559/* R5RS, section 7.1.1, lexical structure of numbers: <ureal R> */
2560
2561static SCM
2562mem2ureal (const char* mem, size_t len, unsigned int *p_idx,
2563 unsigned int radix, enum t_exactness *p_exactness)
0f2d19dd 2564{
3c9a524f 2565 unsigned int idx = *p_idx;
164d2481 2566 SCM result;
3c9a524f
DH		 2567
2568 if (idx == len)
2569 return SCM_BOOL_F;
2570
7351e207
MV		 2571 if (idx+5 <= len && !strncmp (mem+idx, "inf.0", 5))
2572 {
2573 *p_idx = idx+5;
2574 return scm_inf ();
2575 }
2576
2577 if (idx+4 < len && !strncmp (mem+idx, "nan.", 4))
2578 {
2579 enum t_exactness x = EXACT;
2580
d8592269
MV		 2581 /* Cobble up the fractional part. We might want to set the
2582 NaN's mantissa from it. */
7351e207
MV		 2583 idx += 4;
2584 mem2uinteger (mem, len, &idx, 10, &x);
2585 *p_idx = idx;
2586 return scm_nan ();
2587 }
2588
3c9a524f
DH		 2589 if (mem[idx] == '.')
2590 {
2591 if (radix != 10)
2592 return SCM_BOOL_F;
2593 else if (idx + 1 == len)
2594 return SCM_BOOL_F;
71df73ac 2595 else if (!isdigit ((int) (unsigned char) mem[idx + 1]))
3c9a524f
DH		 2596 return SCM_BOOL_F;
2597 else
164d2481
MV		 2598 result = mem2decimal_from_point (SCM_MAKINUM (0), mem, len,
2599 p_idx, p_exactness);
f872b822 2600 }
3c9a524f
DH		 2601 else
2602 {
2603 enum t_exactness x = EXACT;
2604 SCM uinteger;
3c9a524f
DH		 2605
2606 uinteger = mem2uinteger (mem, len, &idx, radix, &x);
2607 if (SCM_FALSEP (uinteger))
2608 return SCM_BOOL_F;
2609
2610 if (idx == len)
2611 result = uinteger;
2612 else if (mem[idx] == '/')
f872b822 2613 {
3c9a524f
DH		 2614 SCM divisor;
2615
2616 idx++;
2617
2618 divisor = mem2uinteger (mem, len, &idx, radix, &x);
2619 if (SCM_FALSEP (divisor))
2620 return SCM_BOOL_F;
2621
f92e85f7
MV		 2622 /* both are int/big here, I assume */
2623 result = scm_make_ratio (uinteger, divisor);
f872b822 2624 }
3c9a524f
DH		 2625 else if (radix == 10)
2626 {
2627 result = mem2decimal_from_point (uinteger, mem, len, &idx, &x);
2628 if (SCM_FALSEP (result))
2629 return SCM_BOOL_F;
2630 }
2631 else
2632 result = uinteger;
2633
2634 *p_idx = idx;
2635 if (x == INEXACT)
2636 *p_exactness = x;
f872b822 2637 }
164d2481
MV		 2638
2639 /* When returning an inexact zero, make sure it is represented as a
2640 floating point value so that we can change its sign.
2641 */
2642 if (SCM_EQ_P (result, SCM_MAKINUM(0)) && *p_exactness == INEXACT)
2643 result = scm_make_real (0.0);
2644
2645 return result;
3c9a524f 2646}
0f2d19dd 2647
0f2d19dd 2648
3c9a524f 2649/* R5RS, section 7.1.1, lexical structure of numbers: <complex R> */
0f2d19dd 2650
3c9a524f
DH		 2651static SCM
2652mem2complex (const char* mem, size_t len, unsigned int idx,
2653 unsigned int radix, enum t_exactness *p_exactness)
2654{
2655 char c;
2656 int sign = 0;
2657 SCM ureal;
2658
2659 if (idx == len)
2660 return SCM_BOOL_F;
2661
2662 c = mem[idx];
2663 if (c == '+')
2664 {
2665 idx++;
2666 sign = 1;
2667 }
2668 else if (c == '-')
2669 {
2670 idx++;
2671 sign = -1;
0f2d19dd 2672 }
0f2d19dd 2673
3c9a524f
DH		 2674 if (idx == len)
2675 return SCM_BOOL_F;
2676
2677 ureal = mem2ureal (mem, len, &idx, radix, p_exactness);
2678 if (SCM_FALSEP (ureal))
f872b822 2679 {
3c9a524f
DH		 2680 /* input must be either +i or -i */
2681
2682 if (sign == 0)
2683 return SCM_BOOL_F;
2684
2685 if (mem[idx] == 'i' || mem[idx] == 'I')
f872b822 2686 {
3c9a524f
DH		 2687 idx++;
2688 if (idx != len)
2689 return SCM_BOOL_F;
2690
2691 return scm_make_rectangular (SCM_MAKINUM (0), SCM_MAKINUM (sign));
f872b822 2692 }
3c9a524f
DH		 2693 else
2694 return SCM_BOOL_F;
0f2d19dd 2695 }
3c9a524f
DH		 2696 else
2697 {
fc194577 2698 if (sign == -1 && SCM_FALSEP (scm_nan_p (ureal)))
3c9a524f 2699 ureal = scm_difference (ureal, SCM_UNDEFINED);
f872b822 2700
3c9a524f
DH		 2701 if (idx == len)
2702 return ureal;
2703
2704 c = mem[idx];
2705 switch (c)
f872b822 2706 {
3c9a524f
DH		 2707 case 'i': case 'I':
2708 /* either +<ureal>i or -<ureal>i */
2709
2710 idx++;
2711 if (sign == 0)
2712 return SCM_BOOL_F;
2713 if (idx != len)
2714 return SCM_BOOL_F;
2715 return scm_make_rectangular (SCM_MAKINUM (0), ureal);
2716
2717 case '@':
2718 /* polar input: <real>@<real>. */
2719
2720 idx++;
2721 if (idx == len)
2722 return SCM_BOOL_F;
2723 else
f872b822 2724 {
3c9a524f
DH		 2725 int sign;
2726 SCM angle;
2727 SCM result;
2728
2729 c = mem[idx];
2730 if (c == '+')
2731 {
2732 idx++;
2733 sign = 1;
2734 }
2735 else if (c == '-')
2736 {
2737 idx++;
2738 sign = -1;
2739 }
2740 else
2741 sign = 1;
2742
2743 angle = mem2ureal (mem, len, &idx, radix, p_exactness);
2744 if (SCM_FALSEP (angle))
2745 return SCM_BOOL_F;
2746 if (idx != len)
2747 return SCM_BOOL_F;
2748
fc194577 2749 if (sign == -1 && SCM_FALSEP (scm_nan_p (ureal)))
3c9a524f
DH		 2750 angle = scm_difference (angle, SCM_UNDEFINED);
2751
2752 result = scm_make_polar (ureal, angle);
2753 return result;
f872b822 2754 }
3c9a524f
DH		 2755 case '+':
2756 case '-':
2757 /* expecting input matching <real>[+-]<ureal>?i */
0f2d19dd 2758
3c9a524f
DH		 2759 idx++;
2760 if (idx == len)
2761 return SCM_BOOL_F;
2762 else
2763 {
2764 int sign = (c == '+') ? 1 : -1;
2765 SCM imag = mem2ureal (mem, len, &idx, radix, p_exactness);
0f2d19dd 2766
3c9a524f
DH		 2767 if (SCM_FALSEP (imag))
2768 imag = SCM_MAKINUM (sign);
fc194577 2769 else if (sign == -1 && SCM_FALSEP (scm_nan_p (ureal)))
1fe5e088 2770 imag = scm_difference (imag, SCM_UNDEFINED);
0f2d19dd 2771
3c9a524f
DH		 2772 if (idx == len)
2773 return SCM_BOOL_F;
2774 if (mem[idx] != 'i' && mem[idx] != 'I')
2775 return SCM_BOOL_F;
0f2d19dd 2776
3c9a524f
DH		 2777 idx++;
2778 if (idx != len)
2779 return SCM_BOOL_F;
0f2d19dd 2780
1fe5e088 2781 return scm_make_rectangular (ureal, imag);
3c9a524f
DH		 2782 }
2783 default:
2784 return SCM_BOOL_F;
2785 }
2786 }
0f2d19dd 2787}
0f2d19dd
JB		 2788
2789
3c9a524f
DH		 2790/* R5RS, section 7.1.1, lexical structure of numbers: <number> */
2791
2792enum t_radix {NO_RADIX=0, DUAL=2, OCT=8, DEC=10, HEX=16};
1cc91f1b 2793
0f2d19dd 2794SCM
3c9a524f 2795scm_i_mem2number (const char* mem, size_t len, unsigned int default_radix)
0f2d19dd 2796{
3c9a524f
DH		 2797 unsigned int idx = 0;
2798 unsigned int radix = NO_RADIX;
2799 enum t_exactness forced_x = NO_EXACTNESS;
2800 enum t_exactness implicit_x = EXACT;
2801 SCM result;
2802
2803 /* R5RS, section 7.1.1, lexical structure of numbers: <prefix R> */
2804 while (idx + 2 < len && mem[idx] == '#')
2805 {
2806 switch (mem[idx + 1])
2807 {
2808 case 'b': case 'B':
2809 if (radix != NO_RADIX)
2810 return SCM_BOOL_F;
2811 radix = DUAL;
2812 break;
2813 case 'd': case 'D':
2814 if (radix != NO_RADIX)
2815 return SCM_BOOL_F;
2816 radix = DEC;
2817 break;
2818 case 'i': case 'I':
2819 if (forced_x != NO_EXACTNESS)
2820 return SCM_BOOL_F;
2821 forced_x = INEXACT;
2822 break;
2823 case 'e': case 'E':
2824 if (forced_x != NO_EXACTNESS)
2825 return SCM_BOOL_F;
2826 forced_x = EXACT;
2827 break;
2828 case 'o': case 'O':
2829 if (radix != NO_RADIX)
2830 return SCM_BOOL_F;
2831 radix = OCT;
2832 break;
2833 case 'x': case 'X':
2834 if (radix != NO_RADIX)
2835 return SCM_BOOL_F;
2836 radix = HEX;
2837 break;
2838 default:
f872b822 2839 return SCM_BOOL_F;
3c9a524f
DH		 2840 }
2841 idx += 2;
2842 }
2843
2844 /* R5RS, section 7.1.1, lexical structure of numbers: <complex R> */
2845 if (radix == NO_RADIX)
2846 result = mem2complex (mem, len, idx, default_radix, &implicit_x);
2847 else
2848 result = mem2complex (mem, len, idx, (unsigned int) radix, &implicit_x);
2849
2850 if (SCM_FALSEP (result))
2851 return SCM_BOOL_F;
f872b822 2852
3c9a524f 2853 switch (forced_x)
f872b822 2854 {
3c9a524f
DH		 2855 case EXACT:
2856 if (SCM_INEXACTP (result))
3c9a524f
DH		 2857 return scm_inexact_to_exact (result);
2858 else
2859 return result;
2860 case INEXACT:
2861 if (SCM_INEXACTP (result))
2862 return result;
2863 else
2864 return scm_exact_to_inexact (result);
2865 case NO_EXACTNESS:
2866 default:
2867 if (implicit_x == INEXACT)
2868 {
2869 if (SCM_INEXACTP (result))
2870 return result;
2871 else
2872 return scm_exact_to_inexact (result);
2873 }
2874 else
2875 return result;
f872b822 2876 }
0f2d19dd
JB		 2877}
2878
2879
a1ec6916 2880SCM_DEFINE (scm_string_to_number, "string->number", 1, 1, 0,
bb628794 2881 (SCM string, SCM radix),
1e6808ea 2882 "Return a number of the maximally precise representation\n"
942e5b91 2883 "expressed by the given @var{string}. @var{radix} must be an\n"
5352393c
MG		 2884 "exact integer, either 2, 8, 10, or 16. If supplied, @var{radix}\n"
2885 "is a default radix that may be overridden by an explicit radix\n"
2886 "prefix in @var{string} (e.g. \"#o177\"). If @var{radix} is not\n"
2887 "supplied, then the default radix is 10. If string is not a\n"
2888 "syntactically valid notation for a number, then\n"
2889 "@code{string->number} returns @code{#f}.")
1bbd0b84 2890#define FUNC_NAME s_scm_string_to_number
0f2d19dd
JB		 2891{
2892 SCM answer;
1bbd0b84 2893 int base;
a6d9e5ab 2894 SCM_VALIDATE_STRING (1, string);
34d19ef6 2895 SCM_VALIDATE_INUM_MIN_DEF_COPY (2, radix,2,10, base);
3c9a524f 2896 answer = scm_i_mem2number (SCM_STRING_CHARS (string),
d8592269
MV		 2897 SCM_STRING_LENGTH (string),
2898 base);
bb628794 2899 return scm_return_first (answer, string);
0f2d19dd 2900}
1bbd0b84 2901#undef FUNC_NAME
3c9a524f
DH		 2902
2903
0f2d19dd
JB		 2904/*** END strs->nums ***/
2905
5986c47d 2906
0f2d19dd 2907SCM
f3ae5d60 2908scm_make_real (double x)
0f2d19dd 2909{
3553e1d1
GH		 2910 SCM z = scm_double_cell (scm_tc16_real, 0, 0, 0);
2911
3a9809df 2912 SCM_REAL_VALUE (z) = x;
0f2d19dd
JB		 2913 return z;
2914}
0f2d19dd 2915
5986c47d 2916
f3ae5d60
MD		 2917SCM
2918scm_make_complex (double x, double y)
2919{
0aacf84e 2920 if (y == 0.0)
3a9809df 2921 return scm_make_real (x);
0aacf84e
MD		 2922 else
2923 {
2924 SCM z;
29c4382a 2925 SCM_NEWSMOB (z, scm_tc16_complex, scm_gc_malloc (sizeof (scm_t_complex),
0aacf84e
MD		 2926 "complex"));
2927 SCM_COMPLEX_REAL (z) = x;
2928 SCM_COMPLEX_IMAG (z) = y;
2929 return z;
2930 }
f3ae5d60 2931}
1cc91f1b 2932
5986c47d 2933
0f2d19dd 2934SCM
1bbd0b84 2935scm_bigequal (SCM x, SCM y)
0f2d19dd 2936{
47ae1f0e 2937 int result = mpz_cmp (SCM_I_BIG_MPZ (x), SCM_I_BIG_MPZ (y));
ca46fb90
RB		 2938 scm_remember_upto_here_2 (x, y);
2939 return SCM_BOOL (0 == result);
0f2d19dd
JB		 2940}
2941
0f2d19dd 2942SCM
f3ae5d60 2943scm_real_equalp (SCM x, SCM y)
0f2d19dd 2944{
f3ae5d60 2945 return SCM_BOOL (SCM_REAL_VALUE (x) == SCM_REAL_VALUE (y));
0f2d19dd
JB		 2946}
2947
f3ae5d60
MD		 2948SCM
2949scm_complex_equalp (SCM x, SCM y)
2950{
2951 return SCM_BOOL (SCM_COMPLEX_REAL (x) == SCM_COMPLEX_REAL (y)
2952 && SCM_COMPLEX_IMAG (x) == SCM_COMPLEX_IMAG (y));
2953}
0f2d19dd 2954
f92e85f7
MV		 2955SCM
2956scm_i_fraction_equalp (SCM x, SCM y)
2957{
2958 scm_i_fraction_reduce (x);
2959 scm_i_fraction_reduce (y);
02164269
MV		 2960 if (SCM_FALSEP (scm_equal_p (SCM_FRACTION_NUMERATOR (x),
2961 SCM_FRACTION_NUMERATOR (y)))
2962 || SCM_FALSEP (scm_equal_p (SCM_FRACTION_DENOMINATOR (x),
2963 SCM_FRACTION_DENOMINATOR (y))))
2964 return SCM_BOOL_F;
2965 else
2966 return SCM_BOOL_T;
f92e85f7 2967}
0f2d19dd
JB		 2968
2969
1bbd0b84 2970SCM_REGISTER_PROC (s_number_p, "number?", 1, 0, 0, scm_number_p);
942e5b91
MG		 2971/* "Return @code{#t} if @var{x} is a number, @code{#f}\n"
2972 * "else. Note that the sets of complex, real, rational and\n"
2973 * "integer values form subsets of the set of numbers, i. e. the\n"
2974 * "predicate will be fulfilled for any number."
2975 */
a1ec6916 2976SCM_DEFINE (scm_number_p, "complex?", 1, 0, 0,
1bbd0b84 2977 (SCM x),
942e5b91 2978 "Return @code{#t} if @var{x} is a complex number, @code{#f}\n"
bb2c02f2 2979 "otherwise. Note that the sets of real, rational and integer\n"
942e5b91
MG		 2980 "values form subsets of the set of complex numbers, i. e. the\n"
2981 "predicate will also be fulfilled if @var{x} is a real,\n"
2982 "rational or integer number.")
1bbd0b84 2983#define FUNC_NAME s_scm_number_p
0f2d19dd 2984{
bb628794 2985 return SCM_BOOL (SCM_NUMBERP (x));
0f2d19dd 2986}
1bbd0b84 2987#undef FUNC_NAME
0f2d19dd
JB		 2988
2989
f92e85f7
MV		 2990SCM_DEFINE (scm_real_p, "real?", 1, 0, 0,
2991 (SCM x),
2992 "Return @code{#t} if @var{x} is a real number, @code{#f}\n"
2993 "otherwise. Note that the set of integer values forms a subset of\n"
2994 "the set of real numbers, i. e. the predicate will also be\n"
2995 "fulfilled if @var{x} is an integer number.")
2996#define FUNC_NAME s_scm_real_p
2997{
2998 /* we can't represent irrational numbers. */
2999 return scm_rational_p (x);
3000}
3001#undef FUNC_NAME
3002
3003SCM_DEFINE (scm_rational_p, "rational?", 1, 0, 0,
1bbd0b84 3004 (SCM x),
942e5b91 3005 "Return @code{#t} if @var{x} is a rational number, @code{#f}\n"
bb2c02f2 3006 "otherwise. Note that the set of integer values forms a subset of\n"
942e5b91 3007 "the set of rational numbers, i. e. the predicate will also be\n"
f92e85f7
MV		 3008 "fulfilled if @var{x} is an integer number.")
3009#define FUNC_NAME s_scm_rational_p
0f2d19dd 3010{
0aacf84e 3011 if (SCM_INUMP (x))
0f2d19dd 3012 return SCM_BOOL_T;
0aacf84e 3013 else if (SCM_IMP (x))
0f2d19dd 3014 return SCM_BOOL_F;
0aacf84e 3015 else if (SCM_BIGP (x))
0f2d19dd 3016 return SCM_BOOL_T;
f92e85f7
MV		 3017 else if (SCM_FRACTIONP (x))
3018 return SCM_BOOL_T;
3019 else if (SCM_REALP (x))
3020 /* due to their limited precision, all floating point numbers are
3021 rational as well. */
3022 return SCM_BOOL_T;
0aacf84e 3023 else
bb628794 3024 return SCM_BOOL_F;
0f2d19dd 3025}
1bbd0b84 3026#undef FUNC_NAME
0f2d19dd
JB		 3027
3028
a1ec6916 3029SCM_DEFINE (scm_integer_p, "integer?", 1, 0, 0,
1bbd0b84 3030 (SCM x),
942e5b91
MG		 3031 "Return @code{#t} if @var{x} is an integer number, @code{#f}\n"
3032 "else.")
1bbd0b84 3033#define FUNC_NAME s_scm_integer_p
0f2d19dd
JB		 3034{
3035 double r;
f872b822
MD		 3036 if (SCM_INUMP (x))
3037 return SCM_BOOL_T;
3038 if (SCM_IMP (x))
3039 return SCM_BOOL_F;
f872b822
MD		 3040 if (SCM_BIGP (x))
3041 return SCM_BOOL_T;
3c9a524f 3042 if (!SCM_INEXACTP (x))
f872b822 3043 return SCM_BOOL_F;
3c9a524f 3044 if (SCM_COMPLEXP (x))
f872b822 3045 return SCM_BOOL_F;
5986c47d 3046 r = SCM_REAL_VALUE (x);
f872b822
MD		 3047 if (r == floor (r))
3048 return SCM_BOOL_T;
0f2d19dd
JB		 3049 return SCM_BOOL_F;
3050}
1bbd0b84 3051#undef FUNC_NAME
0f2d19dd
JB		 3052
3053
a1ec6916 3054SCM_DEFINE (scm_inexact_p, "inexact?", 1, 0, 0,
1bbd0b84 3055 (SCM x),
942e5b91
MG		 3056 "Return @code{#t} if @var{x} is an inexact number, @code{#f}\n"
3057 "else.")
1bbd0b84 3058#define FUNC_NAME s_scm_inexact_p
0f2d19dd 3059{
eb927cb9
MV		 3060 if (SCM_INEXACTP (x))
3061 return SCM_BOOL_T;
3062 if (SCM_NUMBERP (x))
3063 return SCM_BOOL_F;
3064 SCM_WRONG_TYPE_ARG (1, x);
0f2d19dd 3065}
1bbd0b84 3066#undef FUNC_NAME
0f2d19dd
JB		 3067
3068
152f82bf 3069SCM_GPROC1 (s_eq_p, "=", scm_tc7_rpsubr, scm_num_eq_p, g_eq_p);
942e5b91 3070/* "Return @code{#t} if all parameters are numerically equal." */
0f2d19dd 3071SCM
6e8d25a6 3072scm_num_eq_p (SCM x, SCM y)
0f2d19dd 3073{
d8b95e27 3074 again:
0aacf84e
MD		 3075 if (SCM_INUMP (x))
3076 {
3077 long xx = SCM_INUM (x);
3078 if (SCM_INUMP (y))
3079 {
3080 long yy = SCM_INUM (y);
3081 return SCM_BOOL (xx == yy);
3082 }
3083 else if (SCM_BIGP (y))
3084 return SCM_BOOL_F;
3085 else if (SCM_REALP (y))
3086 return SCM_BOOL ((double) xx == SCM_REAL_VALUE (y));
3087 else if (SCM_COMPLEXP (y))
3088 return SCM_BOOL (((double) xx == SCM_COMPLEX_REAL (y))
3089 && (0.0 == SCM_COMPLEX_IMAG (y)));
f92e85f7
MV		 3090 else if (SCM_FRACTIONP (y))
3091 return SCM_BOOL_F;
0aacf84e
MD		 3092 else
3093 SCM_WTA_DISPATCH_2 (g_eq_p, x, y, SCM_ARGn, s_eq_p);
f872b822 3094 }
0aacf84e
MD		 3095 else if (SCM_BIGP (x))
3096 {
3097 if (SCM_INUMP (y))
3098 return SCM_BOOL_F;
3099 else if (SCM_BIGP (y))
3100 {
3101 int cmp = mpz_cmp (SCM_I_BIG_MPZ (x), SCM_I_BIG_MPZ (y));
3102 scm_remember_upto_here_2 (x, y);
3103 return SCM_BOOL (0 == cmp);
3104 }
3105 else if (SCM_REALP (y))
3106 {
3107 int cmp;
3108 if (xisnan (SCM_REAL_VALUE (y)))
3109 return SCM_BOOL_F;
3110 cmp = xmpz_cmp_d (SCM_I_BIG_MPZ (x), SCM_REAL_VALUE (y));
3111 scm_remember_upto_here_1 (x);
3112 return SCM_BOOL (0 == cmp);
3113 }
3114 else if (SCM_COMPLEXP (y))
3115 {
3116 int cmp;
3117 if (0.0 != SCM_COMPLEX_IMAG (y))
3118 return SCM_BOOL_F;
3119 if (xisnan (SCM_COMPLEX_REAL (y)))
3120 return SCM_BOOL_F;
3121 cmp = xmpz_cmp_d (SCM_I_BIG_MPZ (x), SCM_COMPLEX_REAL (y));
3122 scm_remember_upto_here_1 (x);
3123 return SCM_BOOL (0 == cmp);
3124 }
f92e85f7
MV		 3125 else if (SCM_FRACTIONP (y))
3126 return SCM_BOOL_F;
0aacf84e
MD		 3127 else
3128 SCM_WTA_DISPATCH_2 (g_eq_p, x, y, SCM_ARGn, s_eq_p);
f4c627b3 3129 }
0aacf84e
MD		 3130 else if (SCM_REALP (x))
3131 {
3132 if (SCM_INUMP (y))
3133 return SCM_BOOL (SCM_REAL_VALUE (x) == (double) SCM_INUM (y));
3134 else if (SCM_BIGP (y))
3135 {
3136 int cmp;
3137 if (xisnan (SCM_REAL_VALUE (x)))
3138 return SCM_BOOL_F;
3139 cmp = xmpz_cmp_d (SCM_I_BIG_MPZ (y), SCM_REAL_VALUE (x));
3140 scm_remember_upto_here_1 (y);
3141 return SCM_BOOL (0 == cmp);
3142 }
3143 else if (SCM_REALP (y))
3144 return SCM_BOOL (SCM_REAL_VALUE (x) == SCM_REAL_VALUE (y));
3145 else if (SCM_COMPLEXP (y))
3146 return SCM_BOOL ((SCM_REAL_VALUE (x) == SCM_COMPLEX_REAL (y))
3147 && (0.0 == SCM_COMPLEX_IMAG (y)));
f92e85f7 3148 else if (SCM_FRACTIONP (y))
d8b95e27
KR		 3149 {
3150 double xx = SCM_REAL_VALUE (x);
3151 if (xisnan (xx))
3152 return SCM_BOOL_F;
3153 if (xisinf (xx))
3154 return SCM_BOOL (xx < 0.0);
3155 x = scm_inexact_to_exact (x); /* with x as frac or int */
3156 goto again;
3157 }
0aacf84e
MD		 3158 else
3159 SCM_WTA_DISPATCH_2 (g_eq_p, x, y, SCM_ARGn, s_eq_p);
f872b822 3160 }
0aacf84e
MD		 3161 else if (SCM_COMPLEXP (x))
3162 {
3163 if (SCM_INUMP (y))
3164 return SCM_BOOL ((SCM_COMPLEX_REAL (x) == (double) SCM_INUM (y))
3165 && (SCM_COMPLEX_IMAG (x) == 0.0));
3166 else if (SCM_BIGP (y))
3167 {
3168 int cmp;
3169 if (0.0 != SCM_COMPLEX_IMAG (x))
3170 return SCM_BOOL_F;
3171 if (xisnan (SCM_COMPLEX_REAL (x)))
3172 return SCM_BOOL_F;
3173 cmp = xmpz_cmp_d (SCM_I_BIG_MPZ (y), SCM_COMPLEX_REAL (x));
3174 scm_remember_upto_here_1 (y);
3175 return SCM_BOOL (0 == cmp);
3176 }
3177 else if (SCM_REALP (y))
3178 return SCM_BOOL ((SCM_COMPLEX_REAL (x) == SCM_REAL_VALUE (y))
3179 && (SCM_COMPLEX_IMAG (x) == 0.0));
3180 else if (SCM_COMPLEXP (y))
3181 return SCM_BOOL ((SCM_COMPLEX_REAL (x) == SCM_COMPLEX_REAL (y))
3182 && (SCM_COMPLEX_IMAG (x) == SCM_COMPLEX_IMAG (y)));
f92e85f7 3183 else if (SCM_FRACTIONP (y))
d8b95e27
KR		 3184 {
3185 double xx;
3186 if (SCM_COMPLEX_IMAG (x) != 0.0)
3187 return SCM_BOOL_F;
3188 xx = SCM_COMPLEX_REAL (x);
3189 if (xisnan (xx))
3190 return SCM_BOOL_F;
3191 if (xisinf (xx))
3192 return SCM_BOOL (xx < 0.0);
3193 x = scm_inexact_to_exact (x); /* with x as frac or int */
3194 goto again;
3195 }
f92e85f7
MV		 3196 else
3197 SCM_WTA_DISPATCH_2 (g_eq_p, x, y, SCM_ARGn, s_eq_p);
3198 }
3199 else if (SCM_FRACTIONP (x))
3200 {
3201 if (SCM_INUMP (y))
3202 return SCM_BOOL_F;
3203 else if (SCM_BIGP (y))
3204 return SCM_BOOL_F;
3205 else if (SCM_REALP (y))
d8b95e27
KR		 3206 {
3207 double yy = SCM_REAL_VALUE (y);
3208 if (xisnan (yy))
3209 return SCM_BOOL_F;
3210 if (xisinf (yy))
3211 return SCM_BOOL (0.0 < yy);
3212 y = scm_inexact_to_exact (y); /* with y as frac or int */
3213 goto again;
3214 }
f92e85f7 3215 else if (SCM_COMPLEXP (y))
d8b95e27
KR		 3216 {
3217 double yy;
3218 if (SCM_COMPLEX_IMAG (y) != 0.0)
3219 return SCM_BOOL_F;
3220 yy = SCM_COMPLEX_REAL (y);
3221 if (xisnan (yy))
3222 return SCM_BOOL_F;
3223 if (xisinf (yy))
3224 return SCM_BOOL (0.0 < yy);
3225 y = scm_inexact_to_exact (y); /* with y as frac or int */
3226 goto again;
3227 }
f92e85f7
MV		 3228 else if (SCM_FRACTIONP (y))
3229 return scm_i_fraction_equalp (x, y);
0aacf84e
MD		 3230 else
3231 SCM_WTA_DISPATCH_2 (g_eq_p, x, y, SCM_ARGn, s_eq_p);
f4c627b3 3232 }
0aacf84e 3233 else
f4c627b3 3234 SCM_WTA_DISPATCH_2 (g_eq_p, x, y, SCM_ARG1, s_eq_p);
0f2d19dd
JB		 3235}
3236
3237
a5f0b599
KR		 3238/* OPTIMIZE-ME: For int/frac and frac/frac compares, the multiplications
3239 done are good for inums, but for bignums an answer can almost always be
3240 had by just examining a few high bits of the operands, as done by GMP in
3241 mpq_cmp. flonum/frac compares likewise, but with the slight complication
3242 of the float exponent to take into account. */
3243
152f82bf 3244SCM_GPROC1 (s_less_p, "<", scm_tc7_rpsubr, scm_less_p, g_less_p);
942e5b91
MG		 3245/* "Return @code{#t} if the list of parameters is monotonically\n"
3246 * "increasing."
3247 */
0f2d19dd 3248SCM
6e8d25a6 3249scm_less_p (SCM x, SCM y)
0f2d19dd 3250{
a5f0b599 3251 again:
0aacf84e
MD		 3252 if (SCM_INUMP (x))
3253 {
3254 long xx = SCM_INUM (x);
3255 if (SCM_INUMP (y))
3256 {
3257 long yy = SCM_INUM (y);
3258 return SCM_BOOL (xx < yy);
3259 }
3260 else if (SCM_BIGP (y))
3261 {
3262 int sgn = mpz_sgn (SCM_I_BIG_MPZ (y));
3263 scm_remember_upto_here_1 (y);
3264 return SCM_BOOL (sgn > 0);
3265 }
3266 else if (SCM_REALP (y))
3267 return SCM_BOOL ((double) xx < SCM_REAL_VALUE (y));
f92e85f7 3268 else if (SCM_FRACTIONP (y))
a5f0b599
KR		 3269 {
3270 /* "x < a/b" becomes "x*b < a" */
3271 int_frac:
3272 x = scm_product (x, SCM_FRACTION_DENOMINATOR (y));
3273 y = SCM_FRACTION_NUMERATOR (y);
3274 goto again;
3275 }
0aacf84e
MD		 3276 else
3277 SCM_WTA_DISPATCH_2 (g_less_p, x, y, SCM_ARGn, s_less_p);
f872b822 3278 }
0aacf84e
MD		 3279 else if (SCM_BIGP (x))
3280 {
3281 if (SCM_INUMP (y))
3282 {
3283 int sgn = mpz_sgn (SCM_I_BIG_MPZ (x));
3284 scm_remember_upto_here_1 (x);
3285 return SCM_BOOL (sgn < 0);
3286 }
3287 else if (SCM_BIGP (y))
3288 {
3289 int cmp = mpz_cmp (SCM_I_BIG_MPZ (x), SCM_I_BIG_MPZ (y));
3290 scm_remember_upto_here_2 (x, y);
3291 return SCM_BOOL (cmp < 0);
3292 }
3293 else if (SCM_REALP (y))
3294 {
3295 int cmp;
3296 if (xisnan (SCM_REAL_VALUE (y)))
3297 return SCM_BOOL_F;
3298 cmp = xmpz_cmp_d (SCM_I_BIG_MPZ (x), SCM_REAL_VALUE (y));
3299 scm_remember_upto_here_1 (x);
3300 return SCM_BOOL (cmp < 0);
3301 }
f92e85f7 3302 else if (SCM_FRACTIONP (y))
a5f0b599 3303 goto int_frac;
0aacf84e
MD		 3304 else
3305 SCM_WTA_DISPATCH_2 (g_less_p, x, y, SCM_ARGn, s_less_p);
f4c627b3 3306 }
0aacf84e
MD		 3307 else if (SCM_REALP (x))
3308 {
3309 if (SCM_INUMP (y))
3310 return SCM_BOOL (SCM_REAL_VALUE (x) < (double) SCM_INUM (y));
3311 else if (SCM_BIGP (y))
3312 {
3313 int cmp;
3314 if (xisnan (SCM_REAL_VALUE (x)))
3315 return SCM_BOOL_F;
3316 cmp = xmpz_cmp_d (SCM_I_BIG_MPZ (y), SCM_REAL_VALUE (x));
3317 scm_remember_upto_here_1 (y);
3318 return SCM_BOOL (cmp > 0);
3319 }
3320 else if (SCM_REALP (y))
3321 return SCM_BOOL (SCM_REAL_VALUE (x) < SCM_REAL_VALUE (y));
f92e85f7 3322 else if (SCM_FRACTIONP (y))
a5f0b599
KR		 3323 {
3324 double xx = SCM_REAL_VALUE (x);
3325 if (xisnan (xx))
3326 return SCM_BOOL_F;
3327 if (xisinf (xx))
3328 return SCM_BOOL (xx < 0.0);
3329 x = scm_inexact_to_exact (x); /* with x as frac or int */
3330 goto again;
3331 }
f92e85f7
MV		 3332 else
3333 SCM_WTA_DISPATCH_2 (g_less_p, x, y, SCM_ARGn, s_less_p);
3334 }
3335 else if (SCM_FRACTIONP (x))
3336 {
a5f0b599
KR		 3337 if (SCM_INUMP (y) || SCM_BIGP (y))
3338 {
3339 /* "a/b < y" becomes "a < y*b" */
3340 y = scm_product (y, SCM_FRACTION_DENOMINATOR (x));
3341 x = SCM_FRACTION_NUMERATOR (x);
3342 goto again;
3343 }
f92e85f7 3344 else if (SCM_REALP (y))
a5f0b599
KR		 3345 {
3346 double yy = SCM_REAL_VALUE (y);
3347 if (xisnan (yy))
3348 return SCM_BOOL_F;
3349 if (xisinf (yy))
3350 return SCM_BOOL (0.0 < yy);
3351 y = scm_inexact_to_exact (y); /* with y as frac or int */
3352 goto again;
3353 }
f92e85f7 3354 else if (SCM_FRACTIONP (y))
a5f0b599
KR		 3355 {
3356 /* "a/b < c/d" becomes "a*d < c*b" */
3357 SCM new_x = scm_product (SCM_FRACTION_NUMERATOR (x),
3358 SCM_FRACTION_DENOMINATOR (y));
3359 SCM new_y = scm_product (SCM_FRACTION_NUMERATOR (y),
3360 SCM_FRACTION_DENOMINATOR (x));
3361 x = new_x;
3362 y = new_y;
3363 goto again;
3364 }
0aacf84e
MD		 3365 else
3366 SCM_WTA_DISPATCH_2 (g_less_p, x, y, SCM_ARGn, s_less_p);
f872b822 3367 }
0aacf84e 3368 else
f4c627b3 3369 SCM_WTA_DISPATCH_2 (g_less_p, x, y, SCM_ARG1, s_less_p);
0f2d19dd
JB		 3370}
3371
3372
c76b1eaf 3373SCM_GPROC1 (s_scm_gr_p, ">", scm_tc7_rpsubr, scm_gr_p, g_gr_p);
942e5b91
MG		 3374/* "Return @code{#t} if the list of parameters is monotonically\n"
3375 * "decreasing."
c76b1eaf 3376 */
1bbd0b84 3377#define FUNC_NAME s_scm_gr_p
c76b1eaf
MD		 3378SCM
3379scm_gr_p (SCM x, SCM y)
0f2d19dd 3380{
c76b1eaf
MD		 3381 if (!SCM_NUMBERP (x))
3382 SCM_WTA_DISPATCH_2 (g_gr_p, x, y, SCM_ARG1, FUNC_NAME);
3383 else if (!SCM_NUMBERP (y))
3384 SCM_WTA_DISPATCH_2 (g_gr_p, x, y, SCM_ARG2, FUNC_NAME);
3385 else
3386 return scm_less_p (y, x);
0f2d19dd 3387}
1bbd0b84 3388#undef FUNC_NAME
0f2d19dd
JB		 3389
3390
c76b1eaf 3391SCM_GPROC1 (s_scm_leq_p, "<=", scm_tc7_rpsubr, scm_leq_p, g_leq_p);
942e5b91 3392/* "Return @code{#t} if the list of parameters is monotonically\n"
c76b1eaf
MD		 3393 * "non-decreasing."
3394 */
1bbd0b84 3395#define FUNC_NAME s_scm_leq_p
c76b1eaf
MD		 3396SCM
3397scm_leq_p (SCM x, SCM y)
0f2d19dd 3398{
c76b1eaf
MD		 3399 if (!SCM_NUMBERP (x))
3400 SCM_WTA_DISPATCH_2 (g_leq_p, x, y, SCM_ARG1, FUNC_NAME);
3401 else if (!SCM_NUMBERP (y))
3402 SCM_WTA_DISPATCH_2 (g_leq_p, x, y, SCM_ARG2, FUNC_NAME);
fc194577
MV		 3403 else if (SCM_NFALSEP (scm_nan_p (x)) || SCM_NFALSEP (scm_nan_p (y)))
3404 return SCM_BOOL_F;
c76b1eaf
MD		 3405 else
3406 return SCM_BOOL_NOT (scm_less_p (y, x));
0f2d19dd 3407}
1bbd0b84 3408#undef FUNC_NAME
0f2d19dd
JB		 3409
3410
c76b1eaf 3411SCM_GPROC1 (s_scm_geq_p, ">=", scm_tc7_rpsubr, scm_geq_p, g_geq_p);
942e5b91 3412/* "Return @code{#t} if the list of parameters is monotonically\n"
c76b1eaf
MD		 3413 * "non-increasing."
3414 */
1bbd0b84 3415#define FUNC_NAME s_scm_geq_p
c76b1eaf
MD		 3416SCM
3417scm_geq_p (SCM x, SCM y)
0f2d19dd 3418{
c76b1eaf
MD		 3419 if (!SCM_NUMBERP (x))
3420 SCM_WTA_DISPATCH_2 (g_geq_p, x, y, SCM_ARG1, FUNC_NAME);
3421 else if (!SCM_NUMBERP (y))
3422 SCM_WTA_DISPATCH_2 (g_geq_p, x, y, SCM_ARG2, FUNC_NAME);
fc194577
MV		 3423 else if (SCM_NFALSEP (scm_nan_p (x)) || SCM_NFALSEP (scm_nan_p (y)))
3424 return SCM_BOOL_F;
c76b1eaf 3425 else
fc194577 3426 return SCM_BOOL_NOT (scm_less_p (x, y));
0f2d19dd 3427}
1bbd0b84 3428#undef FUNC_NAME
0f2d19dd
JB		 3429
3430
152f82bf 3431SCM_GPROC (s_zero_p, "zero?", 1, 0, 0, scm_zero_p, g_zero_p);
942e5b91
MG		 3432/* "Return @code{#t} if @var{z} is an exact or inexact number equal to\n"
3433 * "zero."
3434 */
0f2d19dd 3435SCM
6e8d25a6 3436scm_zero_p (SCM z)
0f2d19dd 3437{
0aacf84e 3438 if (SCM_INUMP (z))
c2ff8ab0 3439 return SCM_BOOL (SCM_EQ_P (z, SCM_INUM0));
0aacf84e 3440 else if (SCM_BIGP (z))
c2ff8ab0 3441 return SCM_BOOL_F;
0aacf84e 3442 else if (SCM_REALP (z))
c2ff8ab0 3443 return SCM_BOOL (SCM_REAL_VALUE (z) == 0.0);
0aacf84e 3444 else if (SCM_COMPLEXP (z))
c2ff8ab0
DH		 3445 return SCM_BOOL (SCM_COMPLEX_REAL (z) == 0.0
3446 && SCM_COMPLEX_IMAG (z) == 0.0);
f92e85f7
MV		 3447 else if (SCM_FRACTIONP (z))
3448 return SCM_BOOL_F;
0aacf84e 3449 else
c2ff8ab0 3450 SCM_WTA_DISPATCH_1 (g_zero_p, z, SCM_ARG1, s_zero_p);
0f2d19dd
JB		 3451}
3452
3453
152f82bf 3454SCM_GPROC (s_positive_p, "positive?", 1, 0, 0, scm_positive_p, g_positive_p);
942e5b91
MG		 3455/* "Return @code{#t} if @var{x} is an exact or inexact number greater than\n"
3456 * "zero."
3457 */
0f2d19dd 3458SCM
6e8d25a6 3459scm_positive_p (SCM x)
0f2d19dd 3460{
0aacf84e 3461 if (SCM_INUMP (x))
c2ff8ab0 3462 return SCM_BOOL (SCM_INUM (x) > 0);
0aacf84e
MD		 3463 else if (SCM_BIGP (x))
3464 {
3465 int sgn = mpz_sgn (SCM_I_BIG_MPZ (x));
3466 scm_remember_upto_here_1 (x);
3467 return SCM_BOOL (sgn > 0);
3468 }
3469 else if (SCM_REALP (x))
c2ff8ab0 3470 return SCM_BOOL(SCM_REAL_VALUE (x) > 0.0);
f92e85f7
MV		 3471 else if (SCM_FRACTIONP (x))
3472 return scm_positive_p (SCM_FRACTION_NUMERATOR (x));
0aacf84e 3473 else
c2ff8ab0 3474 SCM_WTA_DISPATCH_1 (g_positive_p, x, SCM_ARG1, s_positive_p);
0f2d19dd
JB		 3475}
3476
3477
152f82bf 3478SCM_GPROC (s_negative_p, "negative?", 1, 0, 0, scm_negative_p, g_negative_p);
942e5b91
MG		 3479/* "Return @code{#t} if @var{x} is an exact or inexact number less than\n"
3480 * "zero."
3481 */
0f2d19dd 3482SCM
6e8d25a6 3483scm_negative_p (SCM x)
0f2d19dd 3484{
0aacf84e 3485 if (SCM_INUMP (x))
c2ff8ab0 3486 return SCM_BOOL (SCM_INUM (x) < 0);
0aacf84e
MD		 3487 else if (SCM_BIGP (x))
3488 {
3489 int sgn = mpz_sgn (SCM_I_BIG_MPZ (x));
3490 scm_remember_upto_here_1 (x);
3491 return SCM_BOOL (sgn < 0);
3492 }
3493 else if (SCM_REALP (x))
c2ff8ab0 3494 return SCM_BOOL(SCM_REAL_VALUE (x) < 0.0);
f92e85f7
MV		 3495 else if (SCM_FRACTIONP (x))
3496 return scm_negative_p (SCM_FRACTION_NUMERATOR (x));
0aacf84e 3497 else
c2ff8ab0 3498 SCM_WTA_DISPATCH_1 (g_negative_p, x, SCM_ARG1, s_negative_p);
0f2d19dd
JB		 3499}
3500
3501
2a06f791
KR		 3502/* scm_min and scm_max return an inexact when either argument is inexact, as
3503 required by r5rs. On that basis, for exact/inexact combinations the
3504 exact is converted to inexact to compare and possibly return. This is
3505 unlike scm_less_p above which takes some trouble to preserve all bits in
3506 its test, such trouble is not required for min and max. */
3507
9de33deb 3508SCM_GPROC1 (s_max, "max", scm_tc7_asubr, scm_max, g_max);
942e5b91
MG		 3509/* "Return the maximum of all parameter values."
3510 */
0f2d19dd 3511SCM
6e8d25a6 3512scm_max (SCM x, SCM y)
0f2d19dd 3513{
0aacf84e
MD		 3514 if (SCM_UNBNDP (y))
3515 {
3516 if (SCM_UNBNDP (x))
3517 SCM_WTA_DISPATCH_0 (g_max, s_max);
dab4e67a 3518 else if (SCM_INUMP(x) || SCM_BIGP(x) || SCM_REALP(x) || SCM_FRACTIONP(x))
0aacf84e
MD		 3519 return x;
3520 else
3521 SCM_WTA_DISPATCH_1 (g_max, x, SCM_ARG1, s_max);
f872b822 3522 }
f4c627b3 3523
0aacf84e
MD		 3524 if (SCM_INUMP (x))
3525 {
3526 long xx = SCM_INUM (x);
3527 if (SCM_INUMP (y))
3528 {
3529 long yy = SCM_INUM (y);
3530 return (xx < yy) ? y : x;
3531 }
3532 else if (SCM_BIGP (y))
3533 {
3534 int sgn = mpz_sgn (SCM_I_BIG_MPZ (y));
3535 scm_remember_upto_here_1 (y);
3536 return (sgn < 0) ? x : y;
3537 }
3538 else if (SCM_REALP (y))
3539 {
3540 double z = xx;
3541 /* if y==NaN then ">" is false and we return NaN */
3542 return (z > SCM_REAL_VALUE (y)) ? scm_make_real (z) : y;
3543 }
f92e85f7
MV		 3544 else if (SCM_FRACTIONP (y))
3545 {
e4bc5d6c
KR		 3546 use_less:
3547 return (SCM_FALSEP (scm_less_p (x, y)) ? x : y);
f92e85f7 3548 }
0aacf84e
MD		 3549 else
3550 SCM_WTA_DISPATCH_2 (g_max, x, y, SCM_ARGn, s_max);
f872b822 3551 }
0aacf84e
MD		 3552 else if (SCM_BIGP (x))
3553 {
3554 if (SCM_INUMP (y))
3555 {
3556 int sgn = mpz_sgn (SCM_I_BIG_MPZ (x));
3557 scm_remember_upto_here_1 (x);
3558 return (sgn < 0) ? y : x;
3559 }
3560 else if (SCM_BIGP (y))
3561 {
3562 int cmp = mpz_cmp (SCM_I_BIG_MPZ (x), SCM_I_BIG_MPZ (y));
3563 scm_remember_upto_here_2 (x, y);
3564 return (cmp > 0) ? x : y;
3565 }
3566 else if (SCM_REALP (y))
3567 {
2a06f791
KR		 3568 /* if y==NaN then xx>yy is false, so we return the NaN y */
3569 double xx, yy;
3570 big_real:
3571 xx = scm_i_big2dbl (x);
3572 yy = SCM_REAL_VALUE (y);
3573 return (xx > yy ? scm_make_real (xx) : y);
0aacf84e 3574 }
f92e85f7
MV		 3575 else if (SCM_FRACTIONP (y))
3576 {
e4bc5d6c 3577 goto use_less;
f92e85f7 3578 }
0aacf84e
MD		 3579 else
3580 SCM_WTA_DISPATCH_2 (g_max, x, y, SCM_ARGn, s_max);
f4c627b3 3581 }
0aacf84e
MD		 3582 else if (SCM_REALP (x))
3583 {
3584 if (SCM_INUMP (y))
3585 {
3586 double z = SCM_INUM (y);
3587 /* if x==NaN then "<" is false and we return NaN */
3588 return (SCM_REAL_VALUE (x) < z) ? scm_make_real (z) : x;
3589 }
3590 else if (SCM_BIGP (y))
3591 {
b6f8f763 3592 SCM_SWAP (x, y);
2a06f791 3593 goto big_real;
0aacf84e
MD		 3594 }
3595 else if (SCM_REALP (y))
3596 {
3597 /* if x==NaN then our explicit check means we return NaN
3598 if y==NaN then ">" is false and we return NaN
3599 calling isnan is unavoidable, since it's the only way to know
3600 which of x or y causes any compares to be false */
3601 double xx = SCM_REAL_VALUE (x);
3602 return (xisnan (xx) || xx > SCM_REAL_VALUE (y)) ? x : y;
3603 }
f92e85f7
MV		 3604 else if (SCM_FRACTIONP (y))
3605 {
3606 double yy = scm_i_fraction2double (y);
3607 double xx = SCM_REAL_VALUE (x);
3608 return (xx < yy) ? scm_make_real (yy) : x;
3609 }
3610 else
3611 SCM_WTA_DISPATCH_2 (g_max, x, y, SCM_ARGn, s_max);
3612 }
3613 else if (SCM_FRACTIONP (x))
3614 {
3615 if (SCM_INUMP (y))
3616 {
e4bc5d6c 3617 goto use_less;
f92e85f7
MV		 3618 }
3619 else if (SCM_BIGP (y))
3620 {
e4bc5d6c 3621 goto use_less;
f92e85f7
MV		 3622 }
3623 else if (SCM_REALP (y))
3624 {
3625 double xx = scm_i_fraction2double (x);
3626 return (xx < SCM_REAL_VALUE (y)) ? y : scm_make_real (xx);
3627 }
3628 else if (SCM_FRACTIONP (y))
3629 {
e4bc5d6c 3630 goto use_less;
f92e85f7 3631 }
0aacf84e
MD		 3632 else
3633 SCM_WTA_DISPATCH_2 (g_max, x, y, SCM_ARGn, s_max);
f872b822 3634 }
0aacf84e 3635 else
f4c627b3 3636 SCM_WTA_DISPATCH_2 (g_max, x, y, SCM_ARG1, s_max);
0f2d19dd
JB		 3637}
3638
3639
9de33deb 3640SCM_GPROC1 (s_min, "min", scm_tc7_asubr, scm_min, g_min);
942e5b91
MG		 3641/* "Return the minium of all parameter values."
3642 */
0f2d19dd 3643SCM
6e8d25a6 3644scm_min (SCM x, SCM y)
0f2d19dd 3645{
0aacf84e
MD		 3646 if (SCM_UNBNDP (y))
3647 {
3648 if (SCM_UNBNDP (x))
3649 SCM_WTA_DISPATCH_0 (g_min, s_min);
dab4e67a 3650 else if (SCM_INUMP(x) || SCM_BIGP(x) || SCM_REALP(x) || SCM_FRACTIONP(x))
0aacf84e
MD		 3651 return x;
3652 else
3653 SCM_WTA_DISPATCH_1 (g_min, x, SCM_ARG1, s_min);
f872b822 3654 }
f4c627b3 3655
0aacf84e
MD		 3656 if (SCM_INUMP (x))
3657 {
3658 long xx = SCM_INUM (x);
3659 if (SCM_INUMP (y))
3660 {
3661 long yy = SCM_INUM (y);
3662 return (xx < yy) ? x : y;
3663 }
3664 else if (SCM_BIGP (y))
3665 {
3666 int sgn = mpz_sgn (SCM_I_BIG_MPZ (y));
3667 scm_remember_upto_here_1 (y);
3668 return (sgn < 0) ? y : x;
3669 }
3670 else if (SCM_REALP (y))
3671 {
3672 double z = xx;
3673 /* if y==NaN then "<" is false and we return NaN */
3674 return (z < SCM_REAL_VALUE (y)) ? scm_make_real (z) : y;
3675 }
f92e85f7
MV		 3676 else if (SCM_FRACTIONP (y))
3677 {
e4bc5d6c
KR		 3678 use_less:
3679 return (SCM_FALSEP (scm_less_p (x, y)) ? y : x);
f92e85f7 3680 }
0aacf84e
MD		 3681 else
3682 SCM_WTA_DISPATCH_2 (g_min, x, y, SCM_ARGn, s_min);
f872b822 3683 }
0aacf84e
MD		 3684 else if (SCM_BIGP (x))
3685 {
3686 if (SCM_INUMP (y))
3687 {
3688 int sgn = mpz_sgn (SCM_I_BIG_MPZ (x));
3689 scm_remember_upto_here_1 (x);
3690 return (sgn < 0) ? x : y;
3691 }
3692 else if (SCM_BIGP (y))
3693 {
3694 int cmp = mpz_cmp (SCM_I_BIG_MPZ (x), SCM_I_BIG_MPZ (y));
3695 scm_remember_upto_here_2 (x, y);
3696 return (cmp > 0) ? y : x;
3697 }
3698 else if (SCM_REALP (y))
3699 {
2a06f791
KR		 3700 /* if y==NaN then xx<yy is false, so we return the NaN y */
3701 double xx, yy;
3702 big_real:
3703 xx = scm_i_big2dbl (x);
3704 yy = SCM_REAL_VALUE (y);
3705 return (xx < yy ? scm_make_real (xx) : y);
0aacf84e 3706 }
f92e85f7
MV		 3707 else if (SCM_FRACTIONP (y))
3708 {
e4bc5d6c 3709 goto use_less;
f92e85f7 3710 }
0aacf84e
MD		 3711 else
3712 SCM_WTA_DISPATCH_2 (g_min, x, y, SCM_ARGn, s_min);
f4c627b3 3713 }
0aacf84e
MD		 3714 else if (SCM_REALP (x))
3715 {
3716 if (SCM_INUMP (y))
3717 {
3718 double z = SCM_INUM (y);
3719 /* if x==NaN then "<" is false and we return NaN */
3720 return (z < SCM_REAL_VALUE (x)) ? scm_make_real (z) : x;
3721 }
3722 else if (SCM_BIGP (y))
3723 {
b6f8f763 3724 SCM_SWAP (x, y);
2a06f791 3725 goto big_real;
0aacf84e
MD		 3726 }
3727 else if (SCM_REALP (y))
3728 {
3729 /* if x==NaN then our explicit check means we return NaN
3730 if y==NaN then "<" is false and we return NaN
3731 calling isnan is unavoidable, since it's the only way to know
3732 which of x or y causes any compares to be false */
3733 double xx = SCM_REAL_VALUE (x);
3734 return (xisnan (xx) || xx < SCM_REAL_VALUE (y)) ? x : y;
3735 }
f92e85f7
MV		 3736 else if (SCM_FRACTIONP (y))
3737 {
3738 double yy = scm_i_fraction2double (y);
3739 double xx = SCM_REAL_VALUE (x);
3740 return (yy < xx) ? scm_make_real (yy) : x;
3741 }
0aacf84e
MD		 3742 else
3743 SCM_WTA_DISPATCH_2 (g_min, x, y, SCM_ARGn, s_min);
f872b822 3744 }
f92e85f7
MV		 3745 else if (SCM_FRACTIONP (x))
3746 {
3747 if (SCM_INUMP (y))
3748 {
e4bc5d6c 3749 goto use_less;
f92e85f7
MV		 3750 }
3751 else if (SCM_BIGP (y))
3752 {
e4bc5d6c 3753 goto use_less;
f92e85f7
MV		 3754 }
3755 else if (SCM_REALP (y))
3756 {
3757 double xx = scm_i_fraction2double (x);
3758 return (SCM_REAL_VALUE (y) < xx) ? y : scm_make_real (xx);
3759 }
3760 else if (SCM_FRACTIONP (y))
3761 {
e4bc5d6c 3762 goto use_less;
f92e85f7
MV		 3763 }
3764 else
3765 SCM_WTA_DISPATCH_2 (g_max, x, y, SCM_ARGn, s_max);
3766 }
0aacf84e 3767 else
f4c627b3 3768 SCM_WTA_DISPATCH_2 (g_min, x, y, SCM_ARG1, s_min);
0f2d19dd
JB		 3769}
3770
3771
9de33deb 3772SCM_GPROC1 (s_sum, "+", scm_tc7_asubr, scm_sum, g_sum);
942e5b91
MG		 3773/* "Return the sum of all parameter values. Return 0 if called without\n"
3774 * "any parameters."
3775 */
0f2d19dd 3776SCM
6e8d25a6 3777scm_sum (SCM x, SCM y)
0f2d19dd 3778{
ca46fb90
RB		 3779 if (SCM_UNBNDP (y))
3780 {
3781 if (SCM_NUMBERP (x)) return x;
3782 if (SCM_UNBNDP (x)) return SCM_INUM0;
98cb6e75 3783 SCM_WTA_DISPATCH_1 (g_sum, x, SCM_ARG1, s_sum);
f872b822 3784 }
c209c88e 3785
ca46fb90
RB		 3786 if (SCM_INUMP (x))
3787 {
3788 if (SCM_INUMP (y))
3789 {
3790 long xx = SCM_INUM (x);
3791 long yy = SCM_INUM (y);
3792 long int z = xx + yy;
3793 return SCM_FIXABLE (z) ? SCM_MAKINUM (z) : scm_i_long2big (z);
3794 }
3795 else if (SCM_BIGP (y))
3796 {
3797 SCM_SWAP (x, y);
3798 goto add_big_inum;
3799 }
3800 else if (SCM_REALP (y))
3801 {
3802 long int xx = SCM_INUM (x);
3803 return scm_make_real (xx + SCM_REAL_VALUE (y));
3804 }
3805 else if (SCM_COMPLEXP (y))
3806 {
3807 long int xx = SCM_INUM (x);
3808 return scm_make_complex (xx + SCM_COMPLEX_REAL (y),
3809 SCM_COMPLEX_IMAG (y));
3810 }
f92e85f7
MV		 3811 else if (SCM_FRACTIONP (y))
3812 return scm_make_ratio (scm_sum (SCM_FRACTION_NUMERATOR (y),
3813 scm_product (x, SCM_FRACTION_DENOMINATOR (y))),
3814 SCM_FRACTION_DENOMINATOR (y));
ca46fb90
RB		 3815 else
3816 SCM_WTA_DISPATCH_2 (g_sum, x, y, SCM_ARGn, s_sum);
0aacf84e
MD		 3817 } else if (SCM_BIGP (x))
3818 {
3819 if (SCM_INUMP (y))
3820 {
3821 long int inum;
3822 int bigsgn;
3823 add_big_inum:
3824 inum = SCM_INUM (y);
3825 if (inum == 0)
3826 return x;
3827 bigsgn = mpz_sgn (SCM_I_BIG_MPZ (x));
3828 if (inum < 0)
3829 {
3830 SCM result = scm_i_mkbig ();
3831 mpz_sub_ui (SCM_I_BIG_MPZ (result), SCM_I_BIG_MPZ (x), - inum);
3832 scm_remember_upto_here_1 (x);
3833 /* we know the result will have to be a bignum */
3834 if (bigsgn == -1)
3835 return result;
3836 return scm_i_normbig (result);
3837 }
3838 else
3839 {
3840 SCM result = scm_i_mkbig ();
3841 mpz_add_ui (SCM_I_BIG_MPZ (result), SCM_I_BIG_MPZ (x), inum);
3842 scm_remember_upto_here_1 (x);
3843 /* we know the result will have to be a bignum */
3844 if (bigsgn == 1)
3845 return result;
3846 return scm_i_normbig (result);
3847 }
3848 }
3849 else if (SCM_BIGP (y))
3850 {
3851 SCM result = scm_i_mkbig ();
3852 int sgn_x = mpz_sgn (SCM_I_BIG_MPZ (x));
3853 int sgn_y = mpz_sgn (SCM_I_BIG_MPZ (y));
3854 mpz_add (SCM_I_BIG_MPZ (result),
3855 SCM_I_BIG_MPZ (x),
3856 SCM_I_BIG_MPZ (y));
3857 scm_remember_upto_here_2 (x, y);
3858 /* we know the result will have to be a bignum */
3859 if (sgn_x == sgn_y)
3860 return result;
3861 return scm_i_normbig (result);
3862 }
3863 else if (SCM_REALP (y))
3864 {
3865 double result = mpz_get_d (SCM_I_BIG_MPZ (x)) + SCM_REAL_VALUE (y);
3866 scm_remember_upto_here_1 (x);
3867 return scm_make_real (result);
3868 }
3869 else if (SCM_COMPLEXP (y))
3870 {
3871 double real_part = (mpz_get_d (SCM_I_BIG_MPZ (x))
3872 + SCM_COMPLEX_REAL (y));
3873 scm_remember_upto_here_1 (x);
3874 return scm_make_complex (real_part, SCM_COMPLEX_IMAG (y));
3875 }
f92e85f7
MV		 3876 else if (SCM_FRACTIONP (y))
3877 return scm_make_ratio (scm_sum (SCM_FRACTION_NUMERATOR (y),
3878 scm_product (x, SCM_FRACTION_DENOMINATOR (y))),
3879 SCM_FRACTION_DENOMINATOR (y));
0aacf84e
MD		 3880 else
3881 SCM_WTA_DISPATCH_2 (g_sum, x, y, SCM_ARGn, s_sum);
0f2d19dd 3882 }
0aacf84e
MD		 3883 else if (SCM_REALP (x))
3884 {
3885 if (SCM_INUMP (y))
3886 return scm_make_real (SCM_REAL_VALUE (x) + SCM_INUM (y));
3887 else if (SCM_BIGP (y))
3888 {
3889 double result = mpz_get_d (SCM_I_BIG_MPZ (y)) + SCM_REAL_VALUE (x);
3890 scm_remember_upto_here_1 (y);
3891 return scm_make_real (result);
3892 }
3893 else if (SCM_REALP (y))
3894 return scm_make_real (SCM_REAL_VALUE (x) + SCM_REAL_VALUE (y));
3895 else if (SCM_COMPLEXP (y))
3896 return scm_make_complex (SCM_REAL_VALUE (x) + SCM_COMPLEX_REAL (y),
3897 SCM_COMPLEX_IMAG (y));
f92e85f7
MV		 3898 else if (SCM_FRACTIONP (y))
3899 return scm_make_real (SCM_REAL_VALUE (x) + scm_i_fraction2double (y));
0aacf84e
MD		 3900 else
3901 SCM_WTA_DISPATCH_2 (g_sum, x, y, SCM_ARGn, s_sum);
f872b822 3902 }
0aacf84e
MD		 3903 else if (SCM_COMPLEXP (x))
3904 {
3905 if (SCM_INUMP (y))
3906 return scm_make_complex (SCM_COMPLEX_REAL (x) + SCM_INUM (y),
3907 SCM_COMPLEX_IMAG (x));
3908 else if (SCM_BIGP (y))
3909 {
3910 double real_part = (mpz_get_d (SCM_I_BIG_MPZ (y))
3911 + SCM_COMPLEX_REAL (x));
3912 scm_remember_upto_here_1 (y);
3913 return scm_make_complex (real_part, SCM_COMPLEX_IMAG (x));
3914 }
3915 else if (SCM_REALP (y))
3916 return scm_make_complex (SCM_COMPLEX_REAL (x) + SCM_REAL_VALUE (y),
3917 SCM_COMPLEX_IMAG (x));
3918 else if (SCM_COMPLEXP (y))
3919 return scm_make_complex (SCM_COMPLEX_REAL (x) + SCM_COMPLEX_REAL (y),
3920 SCM_COMPLEX_IMAG (x) + SCM_COMPLEX_IMAG (y));
f92e85f7
MV		 3921 else if (SCM_FRACTIONP (y))
3922 return scm_make_complex (SCM_COMPLEX_REAL (x) + scm_i_fraction2double (y),
3923 SCM_COMPLEX_IMAG (x));
3924 else
3925 SCM_WTA_DISPATCH_2 (g_sum, x, y, SCM_ARGn, s_sum);
3926 }
3927 else if (SCM_FRACTIONP (x))
3928 {
3929 if (SCM_INUMP (y))
3930 return scm_make_ratio (scm_sum (SCM_FRACTION_NUMERATOR (x),
3931 scm_product (y, SCM_FRACTION_DENOMINATOR (x))),
3932 SCM_FRACTION_DENOMINATOR (x));
3933 else if (SCM_BIGP (y))
3934 return scm_make_ratio (scm_sum (SCM_FRACTION_NUMERATOR (x),
3935 scm_product (y, SCM_FRACTION_DENOMINATOR (x))),
3936 SCM_FRACTION_DENOMINATOR (x));
3937 else if (SCM_REALP (y))
3938 return scm_make_real (SCM_REAL_VALUE (y) + scm_i_fraction2double (x));
3939 else if (SCM_COMPLEXP (y))
3940 return scm_make_complex (SCM_COMPLEX_REAL (y) + scm_i_fraction2double (x),
3941 SCM_COMPLEX_IMAG (y));
3942 else if (SCM_FRACTIONP (y))
3943 /* a/b + c/d = (ad + bc) / bd */
3944 return scm_make_ratio (scm_sum (scm_product (SCM_FRACTION_NUMERATOR (x), SCM_FRACTION_DENOMINATOR (y)),
3945 scm_product (SCM_FRACTION_NUMERATOR (y), SCM_FRACTION_DENOMINATOR (x))),
3946 scm_product (SCM_FRACTION_DENOMINATOR (x), SCM_FRACTION_DENOMINATOR (y)));
0aacf84e
MD		 3947 else
3948 SCM_WTA_DISPATCH_2 (g_sum, x, y, SCM_ARGn, s_sum);
98cb6e75 3949 }
0aacf84e 3950 else
98cb6e75 3951 SCM_WTA_DISPATCH_2 (g_sum, x, y, SCM_ARG1, s_sum);
0f2d19dd
JB		 3952}
3953
3954
9de33deb 3955SCM_GPROC1 (s_difference, "-", scm_tc7_asubr, scm_difference, g_difference);
609c3d30
MG		 3956/* If called with one argument @var{z1}, -@var{z1} returned. Otherwise
3957 * the sum of all but the first argument are subtracted from the first
3958 * argument. */
c05e97b7 3959#define FUNC_NAME s_difference
0f2d19dd 3960SCM
6e8d25a6 3961scm_difference (SCM x, SCM y)
0f2d19dd 3962{
ca46fb90
RB		 3963 if (SCM_UNBNDP (y))
3964 {
3965 if (SCM_UNBNDP (x))
3966 SCM_WTA_DISPATCH_0 (g_difference, s_difference);
3967 else
3968 if (SCM_INUMP (x))
3969 {
3970 long xx = -SCM_INUM (x);
3971 if (SCM_FIXABLE (xx))
3972 return SCM_MAKINUM (xx);
3973 else
3974 return scm_i_long2big (xx);
3975 }
3976 else if (SCM_BIGP (x))
3977 /* FIXME: do we really need to normalize here? */
3978 return scm_i_normbig (scm_i_clonebig (x, 0));
3979 else if (SCM_REALP (x))
3980 return scm_make_real (-SCM_REAL_VALUE (x));
3981 else if (SCM_COMPLEXP (x))
3982 return scm_make_complex (-SCM_COMPLEX_REAL (x),
3983 -SCM_COMPLEX_IMAG (x));
f92e85f7
MV		 3984 else if (SCM_FRACTIONP (x))
3985 return scm_make_ratio (scm_difference (SCM_FRACTION_NUMERATOR (x), SCM_UNDEFINED),
3986 SCM_FRACTION_DENOMINATOR (x));
ca46fb90
RB		 3987 else
3988 SCM_WTA_DISPATCH_1 (g_difference, x, SCM_ARG1, s_difference);
f872b822 3989 }
ca46fb90 3990
0aacf84e
MD		 3991 if (SCM_INUMP (x))
3992 {
3993 if (SCM_INUMP (y))
3994 {
3995 long int xx = SCM_INUM (x);
3996 long int yy = SCM_INUM (y);
3997 long int z = xx - yy;
3998 if (SCM_FIXABLE (z))
3999 return SCM_MAKINUM (z);
4000 else
4001 return scm_i_long2big (z);
4002 }
4003 else if (SCM_BIGP (y))
4004 {
4005 /* inum-x - big-y */
4006 long xx = SCM_INUM (x);
ca46fb90 4007
0aacf84e
MD		 4008 if (xx == 0)
4009 return scm_i_clonebig (y, 0);
4010 else
4011 {
4012 int sgn_y = mpz_sgn (SCM_I_BIG_MPZ (y));
4013 SCM result = scm_i_mkbig ();
ca46fb90 4014
0aacf84e
MD		 4015 if (xx >= 0)
4016 mpz_ui_sub (SCM_I_BIG_MPZ (result), xx, SCM_I_BIG_MPZ (y));
4017 else
4018 {
4019 /* x - y == -(y + -x) */
4020 mpz_add_ui (SCM_I_BIG_MPZ (result), SCM_I_BIG_MPZ (y), -xx);
4021 mpz_neg (SCM_I_BIG_MPZ (result), SCM_I_BIG_MPZ (result));
4022 }
4023 scm_remember_upto_here_1 (y);
ca46fb90 4024
0aacf84e
MD		 4025 if ((xx < 0 && (sgn_y > 0)) || ((xx > 0) && sgn_y < 0))
4026 /* we know the result will have to be a bignum */
4027 return result;
4028 else
4029 return scm_i_normbig (result);
4030 }
4031 }
4032 else if (SCM_REALP (y))
4033 {
4034 long int xx = SCM_INUM (x);
4035 return scm_make_real (xx - SCM_REAL_VALUE (y));
4036 }
4037 else if (SCM_COMPLEXP (y))
4038 {
4039 long int xx = SCM_INUM (x);
4040 return scm_make_complex (xx - SCM_COMPLEX_REAL (y),
4041 - SCM_COMPLEX_IMAG (y));
4042 }
f92e85f7
MV		 4043 else if (SCM_FRACTIONP (y))
4044 /* a - b/c = (ac - b) / c */
4045 return scm_make_ratio (scm_difference (scm_product (x, SCM_FRACTION_DENOMINATOR (y)),
4046 SCM_FRACTION_NUMERATOR (y)),
4047 SCM_FRACTION_DENOMINATOR (y));
0aacf84e
MD		 4048 else
4049 SCM_WTA_DISPATCH_2 (g_difference, x, y, SCM_ARGn, s_difference);
f872b822 4050 }
0aacf84e
MD		 4051 else if (SCM_BIGP (x))
4052 {
4053 if (SCM_INUMP (y))
4054 {
4055 /* big-x - inum-y */
4056 long yy = SCM_INUM (y);
4057 int sgn_x = mpz_sgn (SCM_I_BIG_MPZ (x));
ca46fb90 4058
0aacf84e
MD		 4059 scm_remember_upto_here_1 (x);
4060 if (sgn_x == 0)
4061 return SCM_FIXABLE (-yy) ? SCM_MAKINUM (-yy) : scm_long2num (-yy);
4062 else
4063 {
4064 SCM result = scm_i_mkbig ();
ca46fb90 4065
708f22c6
KR		 4066 if (yy >= 0)
4067 mpz_sub_ui (SCM_I_BIG_MPZ (result), SCM_I_BIG_MPZ (x), yy);
4068 else
4069 mpz_add_ui (SCM_I_BIG_MPZ (result), SCM_I_BIG_MPZ (x), -yy);
0aacf84e 4070 scm_remember_upto_here_1 (x);
ca46fb90 4071
0aacf84e
MD		 4072 if ((sgn_x < 0 && (yy > 0)) || ((sgn_x > 0) && yy < 0))
4073 /* we know the result will have to be a bignum */
4074 return result;
4075 else
4076 return scm_i_normbig (result);
4077 }
4078 }
4079 else if (SCM_BIGP (y))
4080 {
4081 int sgn_x = mpz_sgn (SCM_I_BIG_MPZ (x));
4082 int sgn_y = mpz_sgn (SCM_I_BIG_MPZ (y));
4083 SCM result = scm_i_mkbig ();
4084 mpz_sub (SCM_I_BIG_MPZ (result),
4085 SCM_I_BIG_MPZ (x),
4086 SCM_I_BIG_MPZ (y));
4087 scm_remember_upto_here_2 (x, y);
4088 /* we know the result will have to be a bignum */
4089 if ((sgn_x == 1) && (sgn_y == -1))
4090 return result;
4091 if ((sgn_x == -1) && (sgn_y == 1))
4092 return result;
4093 return scm_i_normbig (result);
4094 }
4095 else if (SCM_REALP (y))
4096 {
4097 double result = mpz_get_d (SCM_I_BIG_MPZ (x)) - SCM_REAL_VALUE (y);
4098 scm_remember_upto_here_1 (x);
4099 return scm_make_real (result);
4100 }
4101 else if (SCM_COMPLEXP (y))
4102 {
4103 double real_part = (mpz_get_d (SCM_I_BIG_MPZ (x))
4104 - SCM_COMPLEX_REAL (y));
4105 scm_remember_upto_here_1 (x);
4106 return scm_make_complex (real_part, - SCM_COMPLEX_IMAG (y));
4107 }
f92e85f7
MV		 4108 else if (SCM_FRACTIONP (y))
4109 return scm_make_ratio (scm_difference (scm_product (x, SCM_FRACTION_DENOMINATOR (y)),
4110 SCM_FRACTION_NUMERATOR (y)),
4111 SCM_FRACTION_DENOMINATOR (y));
0aacf84e 4112 else SCM_WTA_DISPATCH_2 (g_difference, x, y, SCM_ARGn, s_difference);
ca46fb90 4113 }
0aacf84e
MD		 4114 else if (SCM_REALP (x))
4115 {
4116 if (SCM_INUMP (y))
4117 return scm_make_real (SCM_REAL_VALUE (x) - SCM_INUM (y));
4118 else if (SCM_BIGP (y))
4119 {
4120 double result = SCM_REAL_VALUE (x) - mpz_get_d (SCM_I_BIG_MPZ (y));
4121 scm_remember_upto_here_1 (x);
4122 return scm_make_real (result);
4123 }
4124 else if (SCM_REALP (y))
4125 return scm_make_real (SCM_REAL_VALUE (x) - SCM_REAL_VALUE (y));
4126 else if (SCM_COMPLEXP (y))
4127 return scm_make_complex (SCM_REAL_VALUE (x) - SCM_COMPLEX_REAL (y),
4128 -SCM_COMPLEX_IMAG (y));
f92e85f7
MV		 4129 else if (SCM_FRACTIONP (y))
4130 return scm_make_real (SCM_REAL_VALUE (x) - scm_i_fraction2double (y));
0aacf84e
MD		 4131 else
4132 SCM_WTA_DISPATCH_2 (g_difference, x, y, SCM_ARGn, s_difference);
98cb6e75 4133 }
0aacf84e
MD		 4134 else if (SCM_COMPLEXP (x))
4135 {
4136 if (SCM_INUMP (y))
4137 return scm_make_complex (SCM_COMPLEX_REAL (x) - SCM_INUM (y),
4138 SCM_COMPLEX_IMAG (x));
4139 else if (SCM_BIGP (y))
4140 {
4141 double real_part = (SCM_COMPLEX_REAL (x)
4142 - mpz_get_d (SCM_I_BIG_MPZ (y)));
4143 scm_remember_upto_here_1 (x);
4144 return scm_make_complex (real_part, SCM_COMPLEX_IMAG (y));
4145 }
4146 else if (SCM_REALP (y))
4147 return scm_make_complex (SCM_COMPLEX_REAL (x) - SCM_REAL_VALUE (y),
4148 SCM_COMPLEX_IMAG (x));
4149 else if (SCM_COMPLEXP (y))
4150 return scm_make_complex (SCM_COMPLEX_REAL (x) - SCM_COMPLEX_REAL (y),
4151 SCM_COMPLEX_IMAG (x) - SCM_COMPLEX_IMAG (y));
f92e85f7
MV		 4152 else if (SCM_FRACTIONP (y))
4153 return scm_make_complex (SCM_COMPLEX_REAL (x) - scm_i_fraction2double (y),
4154 SCM_COMPLEX_IMAG (x));
4155 else
4156 SCM_WTA_DISPATCH_2 (g_difference, x, y, SCM_ARGn, s_difference);
4157 }
4158 else if (SCM_FRACTIONP (x))
4159 {
4160 if (SCM_INUMP (y))
4161 /* a/b - c = (a - cb) / b */
4162 return scm_make_ratio (scm_difference (SCM_FRACTION_NUMERATOR (x),
4163 scm_product(y, SCM_FRACTION_DENOMINATOR (x))),
4164 SCM_FRACTION_DENOMINATOR (x));
4165 else if (SCM_BIGP (y))
4166 return scm_make_ratio (scm_difference (SCM_FRACTION_NUMERATOR (x),
4167 scm_product(y, SCM_FRACTION_DENOMINATOR (x))),
4168 SCM_FRACTION_DENOMINATOR (x));
4169 else if (SCM_REALP (y))
4170 return scm_make_real (scm_i_fraction2double (x) - SCM_REAL_VALUE (y));
4171 else if (SCM_COMPLEXP (y))
4172 return scm_make_complex (scm_i_fraction2double (x) - SCM_COMPLEX_REAL (y),
4173 -SCM_COMPLEX_IMAG (y));
4174 else if (SCM_FRACTIONP (y))
4175 /* a/b - c/d = (ad - bc) / bd */
4176 return scm_make_ratio (scm_difference (scm_product (SCM_FRACTION_NUMERATOR (x), SCM_FRACTION_DENOMINATOR (y)),
4177 scm_product (SCM_FRACTION_NUMERATOR (y), SCM_FRACTION_DENOMINATOR (x))),
4178 scm_product (SCM_FRACTION_DENOMINATOR (x), SCM_FRACTION_DENOMINATOR (y)));
0aacf84e
MD		 4179 else
4180 SCM_WTA_DISPATCH_2 (g_difference, x, y, SCM_ARGn, s_difference);
98cb6e75 4181 }
0aacf84e 4182 else
98cb6e75 4183 SCM_WTA_DISPATCH_2 (g_difference, x, y, SCM_ARG1, s_difference);
0f2d19dd 4184}
c05e97b7 4185#undef FUNC_NAME
0f2d19dd 4186
ca46fb90 4187
9de33deb 4188SCM_GPROC1 (s_product, "*", scm_tc7_asubr, scm_product, g_product);
942e5b91
MG		 4189/* "Return the product of all arguments. If called without arguments,\n"
4190 * "1 is returned."
4191 */
0f2d19dd 4192SCM
6e8d25a6 4193scm_product (SCM x, SCM y)
0f2d19dd 4194{
0aacf84e
MD		 4195 if (SCM_UNBNDP (y))
4196 {
4197 if (SCM_UNBNDP (x))
4198 return SCM_MAKINUM (1L);
4199 else if (SCM_NUMBERP (x))
4200 return x;
4201 else
4202 SCM_WTA_DISPATCH_1 (g_product, x, SCM_ARG1, s_product);
f872b822 4203 }
ca46fb90 4204
0aacf84e
MD		 4205 if (SCM_INUMP (x))
4206 {
4207 long xx;
f4c627b3 4208
0aacf84e
MD		 4209 intbig:
4210 xx = SCM_INUM (x);
f4c627b3 4211
0aacf84e
MD		 4212 switch (xx)
4213 {
ca46fb90
RB		 4214 case 0: return x; break;
4215 case 1: return y; break;
0aacf84e 4216 }
f4c627b3 4217
0aacf84e
MD		 4218 if (SCM_INUMP (y))
4219 {
4220 long yy = SCM_INUM (y);
4221 long kk = xx * yy;
4222 SCM k = SCM_MAKINUM (kk);
4223 if ((kk == SCM_INUM (k)) && (kk / xx == yy))
4224 return k;
4225 else
4226 {
4227 SCM result = scm_i_long2big (xx);
4228 mpz_mul_si (SCM_I_BIG_MPZ (result), SCM_I_BIG_MPZ (result), yy);
4229 return scm_i_normbig (result);
4230 }
4231 }
4232 else if (SCM_BIGP (y))
4233 {
4234 SCM result = scm_i_mkbig ();
4235 mpz_mul_si (SCM_I_BIG_MPZ (result), SCM_I_BIG_MPZ (y), xx);
4236 scm_remember_upto_here_1 (y);
4237 return result;
4238 }
4239 else if (SCM_REALP (y))
4240 return scm_make_real (xx * SCM_REAL_VALUE (y));
4241 else if (SCM_COMPLEXP (y))
4242 return scm_make_complex (xx * SCM_COMPLEX_REAL (y),
4243 xx * SCM_COMPLEX_IMAG (y));
f92e85f7
MV		 4244 else if (SCM_FRACTIONP (y))
4245 return scm_make_ratio (scm_product (x, SCM_FRACTION_NUMERATOR (y)),
4246 SCM_FRACTION_DENOMINATOR (y));
0aacf84e
MD		 4247 else
4248 SCM_WTA_DISPATCH_2 (g_product, x, y, SCM_ARGn, s_product);
f4c627b3 4249 }
0aacf84e
MD		 4250 else if (SCM_BIGP (x))
4251 {
4252 if (SCM_INUMP (y))
4253 {
4254 SCM_SWAP (x, y);
4255 goto intbig;
4256 }
4257 else if (SCM_BIGP (y))
4258 {
4259 SCM result = scm_i_mkbig ();
4260 mpz_mul (SCM_I_BIG_MPZ (result),
4261 SCM_I_BIG_MPZ (x),
4262 SCM_I_BIG_MPZ (y));
4263 scm_remember_upto_here_2 (x, y);
4264 return result;
4265 }
4266 else if (SCM_REALP (y))
4267 {
4268 double result = mpz_get_d (SCM_I_BIG_MPZ (x)) * SCM_REAL_VALUE (y);
4269 scm_remember_upto_here_1 (x);
4270 return scm_make_real (result);
4271 }
4272 else if (SCM_COMPLEXP (y))
4273 {
4274 double z = mpz_get_d (SCM_I_BIG_MPZ (x));
4275 scm_remember_upto_here_1 (x);
4276 return scm_make_complex (z * SCM_COMPLEX_REAL (y),
4277 z * SCM_COMPLEX_IMAG (y));
4278 }
f92e85f7
MV		 4279 else if (SCM_FRACTIONP (y))
4280 return scm_make_ratio (scm_product (x, SCM_FRACTION_NUMERATOR (y)),
4281 SCM_FRACTION_DENOMINATOR (y));
0aacf84e
MD		 4282 else
4283 SCM_WTA_DISPATCH_2 (g_product, x, y, SCM_ARGn, s_product);
f4c627b3 4284 }
0aacf84e
MD		 4285 else if (SCM_REALP (x))
4286 {
4287 if (SCM_INUMP (y))
4288 return scm_make_real (SCM_INUM (y) * SCM_REAL_VALUE (x));
4289 else if (SCM_BIGP (y))
4290 {
4291 double result = mpz_get_d (SCM_I_BIG_MPZ (y)) * SCM_REAL_VALUE (x);
4292 scm_remember_upto_here_1 (y);
4293 return scm_make_real (result);
4294 }
4295 else if (SCM_REALP (y))
4296 return scm_make_real (SCM_REAL_VALUE (x) * SCM_REAL_VALUE (y));
4297 else if (SCM_COMPLEXP (y))
4298 return scm_make_complex (SCM_REAL_VALUE (x) * SCM_COMPLEX_REAL (y),
4299 SCM_REAL_VALUE (x) * SCM_COMPLEX_IMAG (y));
f92e85f7
MV		 4300 else if (SCM_FRACTIONP (y))
4301 return scm_make_real (SCM_REAL_VALUE (x) * scm_i_fraction2double (y));
0aacf84e
MD		 4302 else
4303 SCM_WTA_DISPATCH_2 (g_product, x, y, SCM_ARGn, s_product);
f4c627b3 4304 }
0aacf84e
MD		 4305 else if (SCM_COMPLEXP (x))
4306 {
4307 if (SCM_INUMP (y))
4308 return scm_make_complex (SCM_INUM (y) * SCM_COMPLEX_REAL (x),
4309 SCM_INUM (y) * SCM_COMPLEX_IMAG (x));
4310 else if (SCM_BIGP (y))
4311 {
4312 double z = mpz_get_d (SCM_I_BIG_MPZ (y));
4313 scm_remember_upto_here_1 (y);
76506335
KR		 4314 return scm_make_complex (z * SCM_COMPLEX_REAL (x),
4315 z * SCM_COMPLEX_IMAG (x));
0aacf84e
MD		 4316 }
4317 else if (SCM_REALP (y))
4318 return scm_make_complex (SCM_REAL_VALUE (y) * SCM_COMPLEX_REAL (x),
4319 SCM_REAL_VALUE (y) * SCM_COMPLEX_IMAG (x));
4320 else if (SCM_COMPLEXP (y))
4321 {
4322 return scm_make_complex (SCM_COMPLEX_REAL (x) * SCM_COMPLEX_REAL (y)
4323 - SCM_COMPLEX_IMAG (x) * SCM_COMPLEX_IMAG (y),
4324 SCM_COMPLEX_REAL (x) * SCM_COMPLEX_IMAG (y)
4325 + SCM_COMPLEX_IMAG (x) * SCM_COMPLEX_REAL (y));
4326 }
f92e85f7
MV		 4327 else if (SCM_FRACTIONP (y))
4328 {
4329 double yy = scm_i_fraction2double (y);
4330 return scm_make_complex (yy * SCM_COMPLEX_REAL (x),
4331 yy * SCM_COMPLEX_IMAG (x));
4332 }
4333 else
4334 SCM_WTA_DISPATCH_2 (g_product, x, y, SCM_ARGn, s_product);
4335 }
4336 else if (SCM_FRACTIONP (x))
4337 {
4338 if (SCM_INUMP (y))
4339 return scm_make_ratio (scm_product (y, SCM_FRACTION_NUMERATOR (x)),
4340 SCM_FRACTION_DENOMINATOR (x));
4341 else if (SCM_BIGP (y))
4342 return scm_make_ratio (scm_product (y, SCM_FRACTION_NUMERATOR (x)),
4343 SCM_FRACTION_DENOMINATOR (x));
4344 else if (SCM_REALP (y))
4345 return scm_make_real (scm_i_fraction2double (x) * SCM_REAL_VALUE (y));
4346 else if (SCM_COMPLEXP (y))
4347 {
4348 double xx = scm_i_fraction2double (x);
4349 return scm_make_complex (xx * SCM_COMPLEX_REAL (y),
4350 xx * SCM_COMPLEX_IMAG (y));
4351 }
4352 else if (SCM_FRACTIONP (y))
4353 /* a/b * c/d = ac / bd */
4354 return scm_make_ratio (scm_product (SCM_FRACTION_NUMERATOR (x),
4355 SCM_FRACTION_NUMERATOR (y)),
4356 scm_product (SCM_FRACTION_DENOMINATOR (x),
4357 SCM_FRACTION_DENOMINATOR (y)));
0aacf84e
MD		 4358 else
4359 SCM_WTA_DISPATCH_2 (g_product, x, y, SCM_ARGn, s_product);
f4c627b3 4360 }
0aacf84e 4361 else
f4c627b3 4362 SCM_WTA_DISPATCH_2 (g_product, x, y, SCM_ARG1, s_product);
0f2d19dd
JB		 4363}
4364
0f2d19dd 4365double
6e8d25a6 4366scm_num2dbl (SCM a, const char *why)
f4c627b3 4367#define FUNC_NAME why
0f2d19dd 4368{
0aacf84e 4369 if (SCM_INUMP (a))
0f2d19dd 4370 return (double) SCM_INUM (a);
0aacf84e
MD		 4371 else if (SCM_BIGP (a))
4372 {
4373 double result = mpz_get_d (SCM_I_BIG_MPZ (a));
4374 scm_remember_upto_here_1 (a);
4375 return result;
4376 }
4377 else if (SCM_REALP (a))
f4c627b3 4378 return (SCM_REAL_VALUE (a));
f92e85f7
MV		 4379 else if (SCM_FRACTIONP (a))
4380 return scm_i_fraction2double (a);
0aacf84e 4381 else
f4c627b3 4382 SCM_WRONG_TYPE_ARG (SCM_ARGn, a);
0f2d19dd 4383}
f4c627b3 4384#undef FUNC_NAME
0f2d19dd 4385
7351e207
MV		 4386#if ((defined (HAVE_ISINF) && defined (HAVE_ISNAN)) \
4387 || (defined (HAVE_FINITE) && defined (HAVE_ISNAN)))
4388#define ALLOW_DIVIDE_BY_ZERO
4389/* #define ALLOW_DIVIDE_BY_EXACT_ZERO */
4390#endif
0f2d19dd 4391
ba74ef4e
MV		 4392/* The code below for complex division is adapted from the GNU
4393 libstdc++, which adapted it from f2c's libF77, and is subject to
4394 this copyright: */
4395
4396/****************************************************************
4397Copyright 1990, 1991, 1992, 1993 by AT&T Bell Laboratories and Bellcore.
4398
4399Permission to use, copy, modify, and distribute this software
4400and its documentation for any purpose and without fee is hereby
4401granted, provided that the above copyright notice appear in all
4402copies and that both that the copyright notice and this
4403permission notice and warranty disclaimer appear in supporting
4404documentation, and that the names of AT&T Bell Laboratories or
4405Bellcore or any of their entities not be used in advertising or
4406publicity pertaining to distribution of the software without
4407specific, written prior permission.
4408
4409AT&T and Bellcore disclaim all warranties with regard to this
4410software, including all implied warranties of merchantability
4411and fitness. In no event shall AT&T or Bellcore be liable for
4412any special, indirect or consequential damages or any damages
4413whatsoever resulting from loss of use, data or profits, whether
4414in an action of contract, negligence or other tortious action,
4415arising out of or in connection with the use or performance of
4416this software.
4417****************************************************************/
4418
9de33deb 4419SCM_GPROC1 (s_divide, "/", scm_tc7_asubr, scm_divide, g_divide);
609c3d30
MG		 4420/* Divide the first argument by the product of the remaining
4421 arguments. If called with one argument @var{z1}, 1/@var{z1} is
4422 returned. */
c05e97b7 4423#define FUNC_NAME s_divide
f92e85f7
MV		 4424static SCM
4425scm_i_divide (SCM x, SCM y, int inexact)
0f2d19dd 4426{
f8de44c1
DH		 4427 double a;
4428
0aacf84e
MD		 4429 if (SCM_UNBNDP (y))
4430 {
4431 if (SCM_UNBNDP (x))
4432 SCM_WTA_DISPATCH_0 (g_divide, s_divide);
4433 else if (SCM_INUMP (x))
4434 {
4435 long xx = SCM_INUM (x);
4436 if (xx == 1 || xx == -1)
4437 return x;
7351e207 4438#ifndef ALLOW_DIVIDE_BY_EXACT_ZERO
0aacf84e
MD		 4439 else if (xx == 0)
4440 scm_num_overflow (s_divide);
7351e207 4441#endif
0aacf84e 4442 else
f92e85f7
MV		 4443 {
4444 if (inexact)
4445 return scm_make_real (1.0 / (double) xx);
4446 else return scm_make_ratio (SCM_MAKINUM(1), x);
4447 }
0aacf84e
MD		 4448 }
4449 else if (SCM_BIGP (x))
f92e85f7
MV		 4450 {
4451 if (inexact)
4452 return scm_make_real (1.0 / scm_i_big2dbl (x));
4453 else return scm_make_ratio (SCM_MAKINUM(1), x);
4454 }
0aacf84e
MD		 4455 else if (SCM_REALP (x))
4456 {
4457 double xx = SCM_REAL_VALUE (x);
7351e207 4458#ifndef ALLOW_DIVIDE_BY_ZERO
0aacf84e
MD		 4459 if (xx == 0.0)
4460 scm_num_overflow (s_divide);
4461 else
7351e207 4462#endif
0aacf84e
MD		 4463 return scm_make_real (1.0 / xx);
4464 }
4465 else if (SCM_COMPLEXP (x))
4466 {
4467 double r = SCM_COMPLEX_REAL (x);
4468 double i = SCM_COMPLEX_IMAG (x);
4469 if (r <= i)
4470 {
4471 double t = r / i;
4472 double d = i * (1.0 + t * t);
4473 return scm_make_complex (t / d, -1.0 / d);
4474 }
4475 else
4476 {
4477 double t = i / r;
4478 double d = r * (1.0 + t * t);
4479 return scm_make_complex (1.0 / d, -t / d);
4480 }
4481 }
f92e85f7
MV		 4482 else if (SCM_FRACTIONP (x))
4483 return scm_make_ratio (SCM_FRACTION_DENOMINATOR (x),
4484 SCM_FRACTION_NUMERATOR (x));
0aacf84e
MD		 4485 else
4486 SCM_WTA_DISPATCH_1 (g_divide, x, SCM_ARG1, s_divide);
f8de44c1 4487 }
f8de44c1 4488
0aacf84e
MD		 4489 if (SCM_INUMP (x))
4490 {
4491 long xx = SCM_INUM (x);
4492 if (SCM_INUMP (y))
4493 {
4494 long yy = SCM_INUM (y);
4495 if (yy == 0)
4496 {
7351e207 4497#ifndef ALLOW_DIVIDE_BY_EXACT_ZERO
0aacf84e 4498 scm_num_overflow (s_divide);
7351e207 4499#else
0aacf84e 4500 return scm_make_real ((double) xx / (double) yy);
7351e207 4501#endif
0aacf84e
MD		 4502 }
4503 else if (xx % yy != 0)
f92e85f7
MV		 4504 {
4505 if (inexact)
4506 return scm_make_real ((double) xx / (double) yy);
4507 else return scm_make_ratio (x, y);
4508 }
0aacf84e
MD		 4509 else
4510 {
4511 long z = xx / yy;
4512 if (SCM_FIXABLE (z))
4513 return SCM_MAKINUM (z);
4514 else
4515 return scm_i_long2big (z);
4516 }
f872b822 4517 }
0aacf84e 4518 else if (SCM_BIGP (y))
f92e85f7
MV		 4519 {
4520 if (inexact)
4521 return scm_make_real ((double) xx / scm_i_big2dbl (y));
4522 else return scm_make_ratio (x, y);
4523 }
0aacf84e
MD		 4524 else if (SCM_REALP (y))
4525 {
4526 double yy = SCM_REAL_VALUE (y);
7351e207 4527#ifndef ALLOW_DIVIDE_BY_ZERO
0aacf84e
MD		 4528 if (yy == 0.0)
4529 scm_num_overflow (s_divide);
4530 else
7351e207 4531#endif
0aacf84e 4532 return scm_make_real ((double) xx / yy);
ba74ef4e 4533 }
0aacf84e
MD		 4534 else if (SCM_COMPLEXP (y))
4535 {
4536 a = xx;
4537 complex_div: /* y _must_ be a complex number */
4538 {
4539 double r = SCM_COMPLEX_REAL (y);
4540 double i = SCM_COMPLEX_IMAG (y);
4541 if (r <= i)
4542 {
4543 double t = r / i;
4544 double d = i * (1.0 + t * t);
4545 return scm_make_complex ((a * t) / d, -a / d);
4546 }
4547 else
4548 {
4549 double t = i / r;
4550 double d = r * (1.0 + t * t);
4551 return scm_make_complex (a / d, -(a * t) / d);
4552 }
4553 }
4554 }
f92e85f7
MV		 4555 else if (SCM_FRACTIONP (y))
4556 /* a / b/c = ac / b */
4557 return scm_make_ratio (scm_product (x, SCM_FRACTION_DENOMINATOR (y)),
4558 SCM_FRACTION_NUMERATOR (y));
0aacf84e
MD		 4559 else
4560 SCM_WTA_DISPATCH_2 (g_divide, x, y, SCM_ARGn, s_divide);
f8de44c1 4561 }
0aacf84e
MD		 4562 else if (SCM_BIGP (x))
4563 {
4564 if (SCM_INUMP (y))
4565 {
4566 long int yy = SCM_INUM (y);
4567 if (yy == 0)
4568 {
7351e207 4569#ifndef ALLOW_DIVIDE_BY_EXACT_ZERO
0aacf84e 4570 scm_num_overflow (s_divide);
7351e207 4571#else
0aacf84e
MD		 4572 int sgn = mpz_sgn (SCM_I_BIG_MPZ (x));
4573 scm_remember_upto_here_1 (x);
4574 return (sgn == 0) ? scm_nan () : scm_inf ();
7351e207 4575#endif
0aacf84e
MD		 4576 }
4577 else if (yy == 1)
4578 return x;
4579 else
4580 {
4581 /* FIXME: HMM, what are the relative performance issues here?
4582 We need to test. Is it faster on average to test
4583 divisible_p, then perform whichever operation, or is it
4584 faster to perform the integer div opportunistically and
4585 switch to real if there's a remainder? For now we take the
4586 middle ground: test, then if divisible, use the faster div
4587 func. */
4588
4589 long abs_yy = yy < 0 ? -yy : yy;
4590 int divisible_p = mpz_divisible_ui_p (SCM_I_BIG_MPZ (x), abs_yy);
4591
4592 if (divisible_p)
4593 {
4594 SCM result = scm_i_mkbig ();
4595 mpz_divexact_ui (SCM_I_BIG_MPZ (result), SCM_I_BIG_MPZ (x), abs_yy);
4596 scm_remember_upto_here_1 (x);
4597 if (yy < 0)
4598 mpz_neg (SCM_I_BIG_MPZ (result), SCM_I_BIG_MPZ (result));
4599 return scm_i_normbig (result);
4600 }
4601 else
f92e85f7
MV		 4602 {
4603 if (inexact)
4604 return scm_make_real (scm_i_big2dbl (x) / (double) yy);
4605 else return scm_make_ratio (x, y);
4606 }
0aacf84e
MD		 4607 }
4608 }
4609 else if (SCM_BIGP (y))
4610 {
4611 int y_is_zero = (mpz_sgn (SCM_I_BIG_MPZ (y)) == 0);
4612 if (y_is_zero)
4613 {
ca46fb90 4614#ifndef ALLOW_DIVIDE_BY_EXACT_ZERO
0aacf84e 4615 scm_num_overflow (s_divide);
f872b822 4616#else
0aacf84e
MD		 4617 int sgn = mpz_sgn (SCM_I_BIG_MPZ (x));
4618 scm_remember_upto_here_1 (x);
4619 return (sgn == 0) ? scm_nan () : scm_inf ();
f872b822 4620#endif
0aacf84e
MD		 4621 }
4622 else
4623 {
4624 /* big_x / big_y */
4625 int divisible_p = mpz_divisible_p (SCM_I_BIG_MPZ (x),
4626 SCM_I_BIG_MPZ (y));
4627 if (divisible_p)
4628 {
4629 SCM result = scm_i_mkbig ();
4630 mpz_divexact (SCM_I_BIG_MPZ (result),
4631 SCM_I_BIG_MPZ (x),
4632 SCM_I_BIG_MPZ (y));
4633 scm_remember_upto_here_2 (x, y);
4634 return scm_i_normbig (result);
4635 }
4636 else
4637 {
f92e85f7
MV		 4638 if (inexact)
4639 {
4640 double dbx = mpz_get_d (SCM_I_BIG_MPZ (x));
4641 double dby = mpz_get_d (SCM_I_BIG_MPZ (y));
4642 scm_remember_upto_here_2 (x, y);
4643 return scm_make_real (dbx / dby);
4644 }
4645 else return scm_make_ratio (x, y);
0aacf84e
MD		 4646 }
4647 }
4648 }
4649 else if (SCM_REALP (y))
4650 {
4651 double yy = SCM_REAL_VALUE (y);
7351e207 4652#ifndef ALLOW_DIVIDE_BY_ZERO
0aacf84e
MD		 4653 if (yy == 0.0)
4654 scm_num_overflow (s_divide);
4655 else
7351e207 4656#endif
0aacf84e
MD		 4657 return scm_make_real (scm_i_big2dbl (x) / yy);
4658 }
4659 else if (SCM_COMPLEXP (y))
4660 {
4661 a = scm_i_big2dbl (x);
4662 goto complex_div;
4663 }
f92e85f7
MV		 4664 else if (SCM_FRACTIONP (y))
4665 return scm_make_ratio (scm_product (x, SCM_FRACTION_DENOMINATOR (y)),
4666 SCM_FRACTION_NUMERATOR (y));
0aacf84e
MD		 4667 else
4668 SCM_WTA_DISPATCH_2 (g_divide, x, y, SCM_ARGn, s_divide);
f872b822 4669 }
0aacf84e
MD		 4670 else if (SCM_REALP (x))
4671 {
4672 double rx = SCM_REAL_VALUE (x);
4673 if (SCM_INUMP (y))
4674 {
4675 long int yy = SCM_INUM (y);
7351e207 4676#ifndef ALLOW_DIVIDE_BY_EXACT_ZERO
0aacf84e
MD		 4677 if (yy == 0)
4678 scm_num_overflow (s_divide);
4679 else
7351e207 4680#endif
0aacf84e
MD		 4681 return scm_make_real (rx / (double) yy);
4682 }
4683 else if (SCM_BIGP (y))
4684 {
4685 double dby = mpz_get_d (SCM_I_BIG_MPZ (y));
4686 scm_remember_upto_here_1 (y);
4687 return scm_make_real (rx / dby);
4688 }
4689 else if (SCM_REALP (y))
4690 {
4691 double yy = SCM_REAL_VALUE (y);
7351e207 4692#ifndef ALLOW_DIVIDE_BY_ZERO
0aacf84e
MD		 4693 if (yy == 0.0)
4694 scm_num_overflow (s_divide);
4695 else
7351e207 4696#endif
0aacf84e
MD		 4697 return scm_make_real (rx / yy);
4698 }
4699 else if (SCM_COMPLEXP (y))
4700 {
4701 a = rx;
4702 goto complex_div;
4703 }
f92e85f7
MV		 4704 else if (SCM_FRACTIONP (y))
4705 return scm_make_real (rx / scm_i_fraction2double (y));
0aacf84e
MD		 4706 else
4707 SCM_WTA_DISPATCH_2 (g_divide, x, y, SCM_ARGn, s_divide);
f872b822 4708 }
0aacf84e
MD		 4709 else if (SCM_COMPLEXP (x))
4710 {
4711 double rx = SCM_COMPLEX_REAL (x);
4712 double ix = SCM_COMPLEX_IMAG (x);
4713 if (SCM_INUMP (y))
4714 {
4715 long int yy = SCM_INUM (y);
7351e207 4716#ifndef ALLOW_DIVIDE_BY_EXACT_ZERO
0aacf84e
MD		 4717 if (yy == 0)
4718 scm_num_overflow (s_divide);
4719 else
7351e207 4720#endif
0aacf84e
MD		 4721 {
4722 double d = yy;
4723 return scm_make_complex (rx / d, ix / d);
4724 }
4725 }
4726 else if (SCM_BIGP (y))
4727 {
4728 double dby = mpz_get_d (SCM_I_BIG_MPZ (y));
4729 scm_remember_upto_here_1 (y);
4730 return scm_make_complex (rx / dby, ix / dby);
4731 }
4732 else if (SCM_REALP (y))
4733 {
4734 double yy = SCM_REAL_VALUE (y);
7351e207 4735#ifndef ALLOW_DIVIDE_BY_ZERO
0aacf84e
MD		 4736 if (yy == 0.0)
4737 scm_num_overflow (s_divide);
4738 else
7351e207 4739#endif
0aacf84e
MD		 4740 return scm_make_complex (rx / yy, ix / yy);
4741 }
4742 else if (SCM_COMPLEXP (y))
4743 {
4744 double ry = SCM_COMPLEX_REAL (y);
4745 double iy = SCM_COMPLEX_IMAG (y);
4746 if (ry <= iy)
4747 {
4748 double t = ry / iy;
4749 double d = iy * (1.0 + t * t);
4750 return scm_make_complex ((rx * t + ix) / d, (ix * t - rx) / d);
4751 }
4752 else
4753 {
4754 double t = iy / ry;
4755 double d = ry * (1.0 + t * t);
4756 return scm_make_complex ((rx + ix * t) / d, (ix - rx * t) / d);
4757 }
4758 }
f92e85f7
MV		 4759 else if (SCM_FRACTIONP (y))
4760 {
4761 double yy = scm_i_fraction2double (y);
4762 return scm_make_complex (rx / yy, ix / yy);
4763 }
0aacf84e
MD		 4764 else
4765 SCM_WTA_DISPATCH_2 (g_divide, x, y, SCM_ARGn, s_divide);
f8de44c1 4766 }
f92e85f7
MV		 4767 else if (SCM_FRACTIONP (x))
4768 {
4769 if (SCM_INUMP (y))
4770 {
4771 long int yy = SCM_INUM (y);
4772#ifndef ALLOW_DIVIDE_BY_EXACT_ZERO
4773 if (yy == 0)
4774 scm_num_overflow (s_divide);
4775 else
4776#endif
4777 return scm_make_ratio (SCM_FRACTION_NUMERATOR (x),
4778 scm_product (SCM_FRACTION_DENOMINATOR (x), y));
4779 }
4780 else if (SCM_BIGP (y))
4781 {
4782 return scm_make_ratio (SCM_FRACTION_NUMERATOR (x),
4783 scm_product (SCM_FRACTION_DENOMINATOR (x), y));
4784 }
4785 else if (SCM_REALP (y))
4786 {
4787 double yy = SCM_REAL_VALUE (y);
4788#ifndef ALLOW_DIVIDE_BY_ZERO
4789 if (yy == 0.0)
4790 scm_num_overflow (s_divide);
4791 else
4792#endif
4793 return scm_make_real (scm_i_fraction2double (x) / yy);
4794 }
4795 else if (SCM_COMPLEXP (y))
4796 {
4797 a = scm_i_fraction2double (x);
4798 goto complex_div;
4799 }
4800 else if (SCM_FRACTIONP (y))
4801 return scm_make_ratio (scm_product (SCM_FRACTION_NUMERATOR (x), SCM_FRACTION_DENOMINATOR (y)),
4802 scm_product (SCM_FRACTION_NUMERATOR (y), SCM_FRACTION_DENOMINATOR (x)));
4803 else
4804 SCM_WTA_DISPATCH_2 (g_divide, x, y, SCM_ARGn, s_divide);
4805 }
0aacf84e 4806 else
f8de44c1 4807 SCM_WTA_DISPATCH_2 (g_divide, x, y, SCM_ARG1, s_divide);
0f2d19dd 4808}
f92e85f7
MV		 4809
4810SCM
4811scm_divide (SCM x, SCM y)
4812{
4813 return scm_i_divide (x, y, 0);
4814}
4815
4816static SCM scm_divide2real (SCM x, SCM y)
4817{
4818 return scm_i_divide (x, y, 1);
4819}
c05e97b7 4820#undef FUNC_NAME
0f2d19dd 4821
fa605590 4822
0f2d19dd 4823double
6e8d25a6 4824scm_asinh (double x)
0f2d19dd 4825{
fa605590
KR		 4826#if HAVE_ASINH
4827 return asinh (x);
4828#else
4829#define asinh scm_asinh
f872b822 4830 return log (x + sqrt (x * x + 1));
fa605590 4831#endif
0f2d19dd 4832}
fa605590
KR		 4833SCM_GPROC1 (s_asinh, "$asinh", scm_tc7_dsubr, (SCM (*)()) asinh, g_asinh);
4834/* "Return the inverse hyperbolic sine of @var{x}."
4835 */
0f2d19dd
JB		 4836
4837
0f2d19dd 4838double
6e8d25a6 4839scm_acosh (double x)
0f2d19dd 4840{
fa605590
KR		 4841#if HAVE_ACOSH
4842 return acosh (x);
4843#else
4844#define acosh scm_acosh
f872b822 4845 return log (x + sqrt (x * x - 1));
fa605590 4846#endif
0f2d19dd 4847}
fa605590
KR		 4848SCM_GPROC1 (s_acosh, "$acosh", scm_tc7_dsubr, (SCM (*)()) acosh, g_acosh);
4849/* "Return the inverse hyperbolic cosine of @var{x}."
4850 */
0f2d19dd
JB		 4851
4852
0f2d19dd 4853double
6e8d25a6 4854scm_atanh (double x)
0f2d19dd 4855{
fa605590
KR		 4856#if HAVE_ATANH
4857 return atanh (x);
4858#else
4859#define atanh scm_atanh
f872b822 4860 return 0.5 * log ((1 + x) / (1 - x));
fa605590 4861#endif
0f2d19dd 4862}
fa605590
KR		 4863SCM_GPROC1 (s_atanh, "$atanh", scm_tc7_dsubr, (SCM (*)()) atanh, g_atanh);
4864/* "Return the inverse hyperbolic tangent of @var{x}."
4865 */
0f2d19dd
JB		 4866
4867
f92e85f7
MV		 4868/* XXX - eventually, we should remove this definition of scm_round and
4869 rename scm_round_number to scm_round. Likewise for scm_truncate
4870 and scm_truncate_number.
4871 */
4872
0f2d19dd 4873double
6e8d25a6 4874scm_truncate (double x)
0f2d19dd 4875{
fa605590
KR		 4876#if HAVE_TRUNC
4877 return trunc (x);
4878#else
4879#define trunc scm_truncate
f872b822
MD		 4880 if (x < 0.0)
4881 return -floor (-x);
4882 return floor (x);
fa605590 4883#endif
0f2d19dd 4884}
0f2d19dd 4885
6187f48b
KR		 4886/* scm_round is done using floor(x+0.5) to round to nearest and with
4887 half-way case (ie. when x is an integer plus 0.5) going upwards. Then
4888 half-way cases are identified and adjusted down if the round-upwards
4889 didn't give the desired even integer.
4890
4891 "plus_half == result" identifies a half-way case. If plus_half, which is
4892 x + 0.5, is an integer then x must be an integer plus 0.5.
4893
4894 An odd "result" value is identified with result/2 != floor(result/2).
4895 This is done with plus_half, since that value is ready for use sooner in
4896 a pipelined cpu, and we're already requiring plus_half == result.
4897
4898 Note however that we need to be careful when x is big and already an
4899 integer. In that case "x+0.5" may round to an adjacent integer, causing
4900 us to return such a value, incorrectly. For instance if the hardware is
4901 in the usual default nearest-even rounding, then for x = 0x1FFFFFFFFFFFFF
4902 (ie. 53 one bits) we will have x+0.5 = 0x20000000000000 and that value
4903 returned. Or if the hardware is in round-upwards mode, then other bigger
4904 values like say x == 2^128 will see x+0.5 rounding up to the next higher
4905 representable value, 2^128+2^76 (or whatever), again incorrect.
4906
4907 These bad roundings of x+0.5 are avoided by testing at the start whether
4908 x is already an integer. If it is then clearly that's the desired result
4909 already. And if it's not then the exponent must be small enough to allow
4910 an 0.5 to be represented, and hence added without a bad rounding. */
4911
0f2d19dd 4912double
6e8d25a6 4913scm_round (double x)
0f2d19dd 4914{
6187f48b
KR		 4915 double plus_half, result;
4916
4917 if (x == floor (x))
4918 return x;
4919
4920 plus_half = x + 0.5;
4921 result = floor (plus_half);
0f2d19dd 4922 /* Adjust so that the scm_round is towards even. */
0aacf84e
MD		 4923 return ((plus_half == result && plus_half / 2 != floor (plus_half / 2))
4924 ? result - 1
4925 : result);
0f2d19dd
JB		 4926}
4927
f92e85f7
MV		 4928SCM_DEFINE (scm_truncate_number, "truncate", 1, 0, 0,
4929 (SCM x),
4930 "Round the number @var{x} towards zero.")
4931#define FUNC_NAME s_scm_truncate_number
4932{
4933 if (SCM_FALSEP (scm_negative_p (x)))
4934 return scm_floor (x);
4935 else
4936 return scm_ceiling (x);
4937}
4938#undef FUNC_NAME
4939
4940static SCM exactly_one_half;
4941
4942SCM_DEFINE (scm_round_number, "round", 1, 0, 0,
4943 (SCM x),
4944 "Round the number @var{x} towards the nearest integer. "
4945 "When it is exactly halfway between two integers, "
4946 "round towards the even one.")
4947#define FUNC_NAME s_scm_round_number
4948{
4949 SCM plus_half = scm_sum (x, exactly_one_half);
4950 SCM result = scm_floor (plus_half);
4951 /* Adjust so that the scm_round is towards even. */
4952 if (!SCM_FALSEP (scm_num_eq_p (plus_half, result))
4953 && !SCM_FALSEP (scm_odd_p (result)))
4954 return scm_difference (result, SCM_MAKINUM (1));
4955 else
4956 return result;
4957}
4958#undef FUNC_NAME
4959
4960SCM_PRIMITIVE_GENERIC (scm_floor, "floor", 1, 0, 0,
4961 (SCM x),
4962 "Round the number @var{x} towards minus infinity.")
4963#define FUNC_NAME s_scm_floor
4964{
4965 if (SCM_INUMP (x) || SCM_BIGP (x))
4966 return x;
4967 else if (SCM_REALP (x))
4968 return scm_make_real (floor (SCM_REAL_VALUE (x)));
4969 else if (SCM_FRACTIONP (x))
4970 {
4971 SCM q = scm_quotient (SCM_FRACTION_NUMERATOR (x),
4972 SCM_FRACTION_DENOMINATOR (x));
4973 if (SCM_FALSEP (scm_negative_p (x)))
4974 {
4975 /* For positive x, rounding towards zero is correct. */
4976 return q;
4977 }
4978 else
4979 {
4980 /* For negative x, we need to return q-1 unless x is an
4981 integer. But fractions are never integer, per our
4982 assumptions. */
4983 return scm_difference (q, SCM_MAKINUM (1));
4984 }
4985 }
4986 else
4987 SCM_WTA_DISPATCH_1 (g_scm_floor, x, 1, s_scm_floor);
4988}
4989#undef FUNC_NAME
4990
4991SCM_PRIMITIVE_GENERIC (scm_ceiling, "ceiling", 1, 0, 0,
4992 (SCM x),
4993 "Round the number @var{x} towards infinity.")
4994#define FUNC_NAME s_scm_ceiling
4995{
4996 if (SCM_INUMP (x) || SCM_BIGP (x))
4997 return x;
4998 else if (SCM_REALP (x))
4999 return scm_make_real (ceil (SCM_REAL_VALUE (x)));
5000 else if (SCM_FRACTIONP (x))
5001 {
5002 SCM q = scm_quotient (SCM_FRACTION_NUMERATOR (x),
5003 SCM_FRACTION_DENOMINATOR (x));
5004 if (SCM_FALSEP (scm_positive_p (x)))
5005 {
5006 /* For negative x, rounding towards zero is correct. */
5007 return q;
5008 }
5009 else
5010 {
5011 /* For positive x, we need to return q+1 unless x is an
5012 integer. But fractions are never integer, per our
5013 assumptions. */
5014 return scm_sum (q, SCM_MAKINUM (1));
5015 }
5016 }
5017 else
5018 SCM_WTA_DISPATCH_1 (g_scm_ceiling, x, 1, s_scm_ceiling);
5019}
5020#undef FUNC_NAME
0f2d19dd 5021
14b18ed6 5022SCM_GPROC1 (s_i_sqrt, "$sqrt", scm_tc7_dsubr, (SCM (*)()) sqrt, g_i_sqrt);
942e5b91
MG		 5023/* "Return the square root of the real number @var{x}."
5024 */
14b18ed6 5025SCM_GPROC1 (s_i_abs, "$abs", scm_tc7_dsubr, (SCM (*)()) fabs, g_i_abs);
942e5b91
MG		 5026/* "Return the absolute value of the real number @var{x}."
5027 */
14b18ed6 5028SCM_GPROC1 (s_i_exp, "$exp", scm_tc7_dsubr, (SCM (*)()) exp, g_i_exp);
942e5b91
MG		 5029/* "Return the @var{x}th power of e."
5030 */
14b18ed6 5031SCM_GPROC1 (s_i_log, "$log", scm_tc7_dsubr, (SCM (*)()) log, g_i_log);
b3fcac34 5032/* "Return the natural logarithm of the real number @var{x}."
942e5b91 5033 */
14b18ed6 5034SCM_GPROC1 (s_i_sin, "$sin", scm_tc7_dsubr, (SCM (*)()) sin, g_i_sin);
942e5b91
MG		 5035/* "Return the sine of the real number @var{x}."
5036 */
14b18ed6 5037SCM_GPROC1 (s_i_cos, "$cos", scm_tc7_dsubr, (SCM (*)()) cos, g_i_cos);
942e5b91
MG		 5038/* "Return the cosine of the real number @var{x}."
5039 */
14b18ed6 5040SCM_GPROC1 (s_i_tan, "$tan", scm_tc7_dsubr, (SCM (*)()) tan, g_i_tan);
942e5b91
MG		 5041/* "Return the tangent of the real number @var{x}."
5042 */
14b18ed6 5043SCM_GPROC1 (s_i_asin, "$asin", scm_tc7_dsubr, (SCM (*)()) asin, g_i_asin);
942e5b91
MG		 5044/* "Return the arc sine of the real number @var{x}."
5045 */
14b18ed6 5046SCM_GPROC1 (s_i_acos, "$acos", scm_tc7_dsubr, (SCM (*)()) acos, g_i_acos);
942e5b91
MG		 5047/* "Return the arc cosine of the real number @var{x}."
5048 */
14b18ed6 5049SCM_GPROC1 (s_i_atan, "$atan", scm_tc7_dsubr, (SCM (*)()) atan, g_i_atan);
942e5b91
MG		 5050/* "Return the arc tangent of the real number @var{x}."
5051 */
14b18ed6 5052SCM_GPROC1 (s_i_sinh, "$sinh", scm_tc7_dsubr, (SCM (*)()) sinh, g_i_sinh);
942e5b91
MG		 5053/* "Return the hyperbolic sine of the real number @var{x}."
5054 */
14b18ed6 5055SCM_GPROC1 (s_i_cosh, "$cosh", scm_tc7_dsubr, (SCM (*)()) cosh, g_i_cosh);
942e5b91
MG		 5056/* "Return the hyperbolic cosine of the real number @var{x}."
5057 */
14b18ed6 5058SCM_GPROC1 (s_i_tanh, "$tanh", scm_tc7_dsubr, (SCM (*)()) tanh, g_i_tanh);
942e5b91
MG		 5059/* "Return the hyperbolic tangent of the real number @var{x}."
5060 */
f872b822
MD		 5061
5062struct dpair
5063{
5064 double x, y;
5065};
5066
27c37006
NJ		 5067static void scm_two_doubles (SCM x,
5068 SCM y,
3eeba8d4
JB		 5069 const char *sstring,
5070 struct dpair * xy);
f872b822
MD		 5071
5072static void
27c37006
NJ		 5073scm_two_doubles (SCM x, SCM y, const char *sstring, struct dpair *xy)
5074{
0aacf84e 5075 if (SCM_INUMP (x))
27c37006 5076 xy->x = SCM_INUM (x);
0aacf84e 5077 else if (SCM_BIGP (x))
1be6b49c 5078 xy->x = scm_i_big2dbl (x);
0aacf84e 5079 else if (SCM_REALP (x))
27c37006 5080 xy->x = SCM_REAL_VALUE (x);
f92e85f7
MV		 5081 else if (SCM_FRACTIONP (x))
5082 xy->x = scm_i_fraction2double (x);
0aacf84e 5083 else
27c37006 5084 scm_wrong_type_arg (sstring, SCM_ARG1, x);
98cb6e75 5085
0aacf84e 5086 if (SCM_INUMP (y))
27c37006 5087 xy->y = SCM_INUM (y);
0aacf84e 5088 else if (SCM_BIGP (y))
1be6b49c 5089 xy->y = scm_i_big2dbl (y);
0aacf84e 5090 else if (SCM_REALP (y))
27c37006 5091 xy->y = SCM_REAL_VALUE (y);
f92e85f7
MV		 5092 else if (SCM_FRACTIONP (y))
5093 xy->y = scm_i_fraction2double (y);
0aacf84e 5094 else
27c37006 5095 scm_wrong_type_arg (sstring, SCM_ARG2, y);
0f2d19dd
JB		 5096}
5097
5098
a1ec6916 5099SCM_DEFINE (scm_sys_expt, "$expt", 2, 0, 0,
27c37006
NJ		 5100 (SCM x, SCM y),
5101 "Return @var{x} raised to the power of @var{y}. This\n"
0137a31b 5102 "procedure does not accept complex arguments.")
1bbd0b84 5103#define FUNC_NAME s_scm_sys_expt
0f2d19dd
JB		 5104{
5105 struct dpair xy;
27c37006 5106 scm_two_doubles (x, y, FUNC_NAME, &xy);
f8de44c1 5107 return scm_make_real (pow (xy.x, xy.y));
0f2d19dd 5108}
1bbd0b84 5109#undef FUNC_NAME
0f2d19dd
JB		 5110
5111
a1ec6916 5112SCM_DEFINE (scm_sys_atan2, "$atan2", 2, 0, 0,
27c37006
NJ		 5113 (SCM x, SCM y),
5114 "Return the arc tangent of the two arguments @var{x} and\n"
5115 "@var{y}. This is similar to calculating the arc tangent of\n"
5116 "@var{x} / @var{y}, except that the signs of both arguments\n"
0137a31b
MG		 5117 "are used to determine the quadrant of the result. This\n"
5118 "procedure does not accept complex arguments.")
1bbd0b84 5119#define FUNC_NAME s_scm_sys_atan2
0f2d19dd
JB		 5120{
5121 struct dpair xy;
27c37006 5122 scm_two_doubles (x, y, FUNC_NAME, &xy);
f8de44c1 5123 return scm_make_real (atan2 (xy.x, xy.y));
0f2d19dd 5124}
1bbd0b84 5125#undef FUNC_NAME
0f2d19dd
JB		 5126
5127
a1ec6916 5128SCM_DEFINE (scm_make_rectangular, "make-rectangular", 2, 0, 0,
bb628794 5129 (SCM real, SCM imaginary),
942e5b91
MG		 5130 "Return a complex number constructed of the given @var{real} and\n"
5131 "@var{imaginary} parts.")
1bbd0b84 5132#define FUNC_NAME s_scm_make_rectangular
0f2d19dd
JB		 5133{
5134 struct dpair xy;
bb628794 5135 scm_two_doubles (real, imaginary, FUNC_NAME, &xy);
f8de44c1 5136 return scm_make_complex (xy.x, xy.y);
0f2d19dd 5137}
1bbd0b84 5138#undef FUNC_NAME
0f2d19dd
JB		 5139
5140
5141
a1ec6916 5142SCM_DEFINE (scm_make_polar, "make-polar", 2, 0, 0,
27c37006 5143 (SCM x, SCM y),
942e5b91 5144 "Return the complex number @var{x} * e^(i * @var{y}).")
1bbd0b84 5145#define FUNC_NAME s_scm_make_polar
0f2d19dd
JB		 5146{
5147 struct dpair xy;
6efadd7c 5148 double s, c;
27c37006 5149 scm_two_doubles (x, y, FUNC_NAME, &xy);
6efadd7c
KR		 5150#if HAVE_SINCOS
5151 sincos (xy.y, &s, &c);
5152#else
5153 s = sin (xy.y);
5154 c = cos (xy.y);
5155#endif
5156 return scm_make_complex (xy.x * c, xy.x * s);
0f2d19dd 5157}
1bbd0b84 5158#undef FUNC_NAME
0f2d19dd
JB		 5159
5160
152f82bf 5161SCM_GPROC (s_real_part, "real-part", 1, 0, 0, scm_real_part, g_real_part);
942e5b91
MG		 5162/* "Return the real part of the number @var{z}."
5163 */
0f2d19dd 5164SCM
6e8d25a6 5165scm_real_part (SCM z)
0f2d19dd 5166{
0aacf84e 5167 if (SCM_INUMP (z))
c2ff8ab0 5168 return z;
0aacf84e 5169 else if (SCM_BIGP (z))
c2ff8ab0 5170 return z;
0aacf84e 5171 else if (SCM_REALP (z))
c2ff8ab0 5172 return z;
0aacf84e 5173 else if (SCM_COMPLEXP (z))
c2ff8ab0 5174 return scm_make_real (SCM_COMPLEX_REAL (z));
f92e85f7 5175 else if (SCM_FRACTIONP (z))
2fa2d879 5176 return z;
0aacf84e 5177 else
c2ff8ab0 5178 SCM_WTA_DISPATCH_1 (g_real_part, z, SCM_ARG1, s_real_part);
0f2d19dd
JB		 5179}
5180
5181
152f82bf 5182SCM_GPROC (s_imag_part, "imag-part", 1, 0, 0, scm_imag_part, g_imag_part);
942e5b91
MG		 5183/* "Return the imaginary part of the number @var{z}."
5184 */
0f2d19dd 5185SCM
6e8d25a6 5186scm_imag_part (SCM z)
0f2d19dd 5187{
0aacf84e 5188 if (SCM_INUMP (z))
f872b822 5189 return SCM_INUM0;
0aacf84e 5190 else if (SCM_BIGP (z))
f872b822 5191 return SCM_INUM0;
0aacf84e 5192 else if (SCM_REALP (z))
c2ff8ab0 5193 return scm_flo0;
0aacf84e 5194 else if (SCM_COMPLEXP (z))
c2ff8ab0 5195 return scm_make_real (SCM_COMPLEX_IMAG (z));
f92e85f7
MV		 5196 else if (SCM_FRACTIONP (z))
5197 return SCM_INUM0;
0aacf84e 5198 else
c2ff8ab0 5199 SCM_WTA_DISPATCH_1 (g_imag_part, z, SCM_ARG1, s_imag_part);
0f2d19dd
JB		 5200}
5201
f92e85f7
MV		 5202SCM_GPROC (s_numerator, "numerator", 1, 0, 0, scm_numerator, g_numerator);
5203/* "Return the numerator of the number @var{z}."
5204 */
5205SCM
5206scm_numerator (SCM z)
5207{
5208 if (SCM_INUMP (z))
5209 return z;
5210 else if (SCM_BIGP (z))
5211 return z;
5212 else if (SCM_FRACTIONP (z))
5213 {
5214 scm_i_fraction_reduce (z);
5215 return SCM_FRACTION_NUMERATOR (z);
5216 }
5217 else if (SCM_REALP (z))
5218 return scm_exact_to_inexact (scm_numerator (scm_inexact_to_exact (z)));
5219 else
5220 SCM_WTA_DISPATCH_1 (g_numerator, z, SCM_ARG1, s_numerator);
5221}
5222
5223
5224SCM_GPROC (s_denominator, "denominator", 1, 0, 0, scm_denominator, g_denominator);
5225/* "Return the denominator of the number @var{z}."
5226 */
5227SCM
5228scm_denominator (SCM z)
5229{
5230 if (SCM_INUMP (z))
5231 return SCM_MAKINUM (1);
5232 else if (SCM_BIGP (z))
5233 return SCM_MAKINUM (1);
5234 else if (SCM_FRACTIONP (z))
5235 {
5236 scm_i_fraction_reduce (z);
5237 return SCM_FRACTION_DENOMINATOR (z);
5238 }
5239 else if (SCM_REALP (z))
5240 return scm_exact_to_inexact (scm_denominator (scm_inexact_to_exact (z)));
5241 else
5242 SCM_WTA_DISPATCH_1 (g_denominator, z, SCM_ARG1, s_denominator);
5243}
0f2d19dd 5244
9de33deb 5245SCM_GPROC (s_magnitude, "magnitude", 1, 0, 0, scm_magnitude, g_magnitude);
942e5b91
MG		 5246/* "Return the magnitude of the number @var{z}. This is the same as\n"
5247 * "@code{abs} for real arguments, but also allows complex numbers."
5248 */
0f2d19dd 5249SCM
6e8d25a6 5250scm_magnitude (SCM z)
0f2d19dd 5251{
0aacf84e
MD		 5252 if (SCM_INUMP (z))
5253 {
5254 long int zz = SCM_INUM (z);
5255 if (zz >= 0)
5256 return z;
5257 else if (SCM_POSFIXABLE (-zz))
5258 return SCM_MAKINUM (-zz);
5259 else
5260 return scm_i_long2big (-zz);
5986c47d 5261 }
0aacf84e
MD		 5262 else if (SCM_BIGP (z))
5263 {
5264 int sgn = mpz_sgn (SCM_I_BIG_MPZ (z));
5265 scm_remember_upto_here_1 (z);
5266 if (sgn < 0)
5267 return scm_i_clonebig (z, 0);
5268 else
5269 return z;
5986c47d 5270 }
0aacf84e 5271 else if (SCM_REALP (z))
c2ff8ab0 5272 return scm_make_real (fabs (SCM_REAL_VALUE (z)));
0aacf84e 5273 else if (SCM_COMPLEXP (z))
6efadd7c 5274 return scm_make_real (hypot (SCM_COMPLEX_REAL (z), SCM_COMPLEX_IMAG (z)));
f92e85f7
MV		 5275 else if (SCM_FRACTIONP (z))
5276 {
5277 if (SCM_FALSEP (scm_negative_p (SCM_FRACTION_NUMERATOR (z))))
5278 return z;
5279 return scm_make_ratio (scm_difference (SCM_FRACTION_NUMERATOR (z), SCM_UNDEFINED),
5280 SCM_FRACTION_DENOMINATOR (z));
5281 }
0aacf84e 5282 else
c2ff8ab0 5283 SCM_WTA_DISPATCH_1 (g_magnitude, z, SCM_ARG1, s_magnitude);
0f2d19dd
JB		 5284}
5285
5286
9de33deb 5287SCM_GPROC (s_angle, "angle", 1, 0, 0, scm_angle, g_angle);
942e5b91
MG		 5288/* "Return the angle of the complex number @var{z}."
5289 */
0f2d19dd 5290SCM
6e8d25a6 5291scm_angle (SCM z)
0f2d19dd 5292{
c8ae173e
KR		 5293 /* atan(0,-1) is pi and it'd be possible to have that as a constant like
5294 scm_flo0 to save allocating a new flonum with scm_make_real each time.
5295 But if atan2 follows the floating point rounding mode, then the value
5296 is not a constant. Maybe it'd be close enough though. */
0aacf84e
MD		 5297 if (SCM_INUMP (z))
5298 {
5299 if (SCM_INUM (z) >= 0)
c8ae173e 5300 return scm_flo0;
0aacf84e
MD		 5301 else
5302 return scm_make_real (atan2 (0.0, -1.0));
f872b822 5303 }
0aacf84e
MD		 5304 else if (SCM_BIGP (z))
5305 {
5306 int sgn = mpz_sgn (SCM_I_BIG_MPZ (z));
5307 scm_remember_upto_here_1 (z);
5308 if (sgn < 0)
5309 return scm_make_real (atan2 (0.0, -1.0));
5310 else
c8ae173e 5311 return scm_flo0;
0f2d19dd 5312 }
0aacf84e 5313 else if (SCM_REALP (z))
c8ae173e
KR		 5314 {
5315 if (SCM_REAL_VALUE (z) >= 0)
5316 return scm_flo0;
5317 else
5318 return scm_make_real (atan2 (0.0, -1.0));
5319 }
0aacf84e 5320 else if (SCM_COMPLEXP (z))
f4c627b3 5321 return scm_make_real (atan2 (SCM_COMPLEX_IMAG (z), SCM_COMPLEX_REAL (z)));
f92e85f7
MV		 5322 else if (SCM_FRACTIONP (z))
5323 {
5324 if (SCM_FALSEP (scm_negative_p (SCM_FRACTION_NUMERATOR (z))))
5325 return scm_flo0;
5326 else return scm_make_real (atan2 (0.0, -1.0));
5327 }
0aacf84e 5328 else
f4c627b3 5329 SCM_WTA_DISPATCH_1 (g_angle, z, SCM_ARG1, s_angle);
0f2d19dd
JB		 5330}
5331
5332
3c9a524f
DH		 5333SCM_GPROC (s_exact_to_inexact, "exact->inexact", 1, 0, 0, scm_exact_to_inexact, g_exact_to_inexact);
5334/* Convert the number @var{x} to its inexact representation.\n"
5335 */
5336SCM
5337scm_exact_to_inexact (SCM z)
5338{
5339 if (SCM_INUMP (z))
5340 return scm_make_real ((double) SCM_INUM (z));
5341 else if (SCM_BIGP (z))
5342 return scm_make_real (scm_i_big2dbl (z));
f92e85f7
MV		 5343 else if (SCM_FRACTIONP (z))
5344 return scm_make_real (scm_i_fraction2double (z));
3c9a524f
DH		 5345 else if (SCM_INEXACTP (z))
5346 return z;
5347 else
5348 SCM_WTA_DISPATCH_1 (g_exact_to_inexact, z, 1, s_exact_to_inexact);
5349}
5350
5351
a1ec6916 5352SCM_DEFINE (scm_inexact_to_exact, "inexact->exact", 1, 0, 0,
1bbd0b84 5353 (SCM z),
1e6808ea 5354 "Return an exact number that is numerically closest to @var{z}.")
1bbd0b84 5355#define FUNC_NAME s_scm_inexact_to_exact
0f2d19dd 5356{
0aacf84e 5357 if (SCM_INUMP (z))
f872b822 5358 return z;
0aacf84e 5359 else if (SCM_BIGP (z))
f872b822 5360 return z;
0aacf84e
MD		 5361 else if (SCM_REALP (z))
5362 {
f92e85f7
MV		 5363 if (xisinf (SCM_REAL_VALUE (z)) || xisnan (SCM_REAL_VALUE (z)))
5364 SCM_OUT_OF_RANGE (1, z);
2be24db4 5365 else
f92e85f7
MV		 5366 {
5367 mpq_t frac;
5368 SCM q;
5369
5370 mpq_init (frac);
5371 mpq_set_d (frac, SCM_REAL_VALUE (z));
5372 q = scm_make_ratio (scm_i_mpz2num (mpq_numref (frac)),
5373 scm_i_mpz2num (mpq_denref (frac)));
5374
5375 /* When scm_make_ratio throws, we leak the memory allocated
5376 for frac...
5377 */
5378 mpq_clear (frac);
5379 return q;
5380 }
c2ff8ab0 5381 }
f92e85f7
MV		 5382 else if (SCM_FRACTIONP (z))
5383 return z;
0aacf84e 5384 else
c2ff8ab0 5385 SCM_WRONG_TYPE_ARG (1, z);
0f2d19dd 5386}
1bbd0b84 5387#undef FUNC_NAME
0f2d19dd 5388
f92e85f7
MV		 5389SCM_DEFINE (scm_rationalize, "rationalize", 2, 0, 0,
5390 (SCM x, SCM err),
5391 "Return an exact number that is within @var{err} of @var{x}.")
5392#define FUNC_NAME s_scm_rationalize
5393{
5394 if (SCM_INUMP (x))
5395 return x;
5396 else if (SCM_BIGP (x))
5397 return x;
5398 else if ((SCM_REALP (x)) || SCM_FRACTIONP (x))
5399 {
5400 /* Use continued fractions to find closest ratio. All
5401 arithmetic is done with exact numbers.
5402 */
5403
5404 SCM ex = scm_inexact_to_exact (x);
5405 SCM int_part = scm_floor (ex);
5406 SCM tt = SCM_MAKINUM (1);
5407 SCM a1 = SCM_MAKINUM (0), a2 = SCM_MAKINUM (1), a = SCM_MAKINUM (0);
5408 SCM b1 = SCM_MAKINUM (1), b2 = SCM_MAKINUM (0), b = SCM_MAKINUM (0);
5409 SCM rx;
5410 int i = 0;
5411
5412 if (!SCM_FALSEP (scm_num_eq_p (ex, int_part)))
5413 return ex;
5414
5415 ex = scm_difference (ex, int_part); /* x = x-int_part */
5416 rx = scm_divide (ex, SCM_UNDEFINED); /* rx = 1/x */
5417
5418 /* We stop after a million iterations just to be absolutely sure
5419 that we don't go into an infinite loop. The process normally
5420 converges after less than a dozen iterations.
5421 */
5422
5423 err = scm_abs (err);
5424 while (++i < 1000000)
5425 {
5426 a = scm_sum (scm_product (a1, tt), a2); /* a = a1*tt + a2 */
5427 b = scm_sum (scm_product (b1, tt), b2); /* b = b1*tt + b2 */
5428 if (SCM_FALSEP (scm_zero_p (b)) && /* b != 0 */
5429 SCM_FALSEP
5430 (scm_gr_p (scm_abs (scm_difference (ex, scm_divide (a, b))),
5431 err))) /* abs(x-a/b) <= err */
02164269
MV		 5432 {
5433 SCM res = scm_sum (int_part, scm_divide (a, b));
5434 if (SCM_FALSEP (scm_exact_p (x))
5435 || SCM_FALSEP (scm_exact_p (err)))
5436 return scm_exact_to_inexact (res);
5437 else
5438 return res;
5439 }
f92e85f7
MV		 5440 rx = scm_divide (scm_difference (rx, tt), /* rx = 1/(rx - tt) */
5441 SCM_UNDEFINED);
5442 tt = scm_floor (rx); /* tt = floor (rx) */
5443 a2 = a1;
5444 b2 = b1;
5445 a1 = a;
5446 b1 = b;
5447 }
5448 scm_num_overflow (s_scm_rationalize);
5449 }
5450 else
5451 SCM_WRONG_TYPE_ARG (1, x);
5452}
5453#undef FUNC_NAME
5454
87617347 5455/* if you need to change this, change test-num2integral.c as well */
ee33d62a 5456#if SCM_SIZEOF_LONG_LONG != 0
1be6b49c
ML		 5457# ifndef LLONG_MAX
5458# define ULLONG_MAX ((unsigned long long) (-1))
5459# define LLONG_MAX ((long long) (ULLONG_MAX >> 1))
5460# define LLONG_MIN (~LLONG_MAX)
5461# endif
f872b822 5462#endif
0f2d19dd 5463
3d2e8ceb
MV		 5464/* Parameters for creating integer conversion routines.
5465
5466 Define the following preprocessor macros before including
5467 "libguile/num2integral.i.c":
5468
5469 NUM2INTEGRAL - the name of the function for converting from a
ca46fb90
RB		 5470 Scheme object to the integral type. This function will be
5471 defined when including "num2integral.i.c".
3d2e8ceb
MV		 5472
5473 INTEGRAL2NUM - the name of the function for converting from the
ca46fb90 5474 integral type to a Scheme object. This function will be defined.
3d2e8ceb
MV		 5475
5476 INTEGRAL2BIG - the name of an internal function that createas a
ca46fb90
RB		 5477 bignum from the integral type. This function will be defined.
5478 The name should start with "scm_i_".
5479
5480 ITYPE - the name of the integral type.
5481
9dd023e1
MV		 5482 UNSIGNED - Define this to 1 when ITYPE is an unsigned type. Define
5483 it to 0 otherwise.
ca46fb90
RB		 5484
5485 UNSIGNED_ITYPE - the name of the the unsigned variant of the
5486 integral type. If you don't define this, it defaults to
5487 "unsigned ITYPE" for signed types and simply "ITYPE" for unsigned
5488 ones.
5489
5490 SIZEOF_ITYPE - an expression giving the size of the integral type
5491 in bytes. This expression must be computable by the
5492 preprocessor. (SIZEOF_FOO values are calculated by configure.in
5493 for common types).
5494
3d2e8ceb
MV		 5495*/
5496
1be6b49c
ML		 5497#define NUM2INTEGRAL scm_num2short
5498#define INTEGRAL2NUM scm_short2num
5499#define INTEGRAL2BIG scm_i_short2big
ca46fb90 5500#define UNSIGNED 0
1be6b49c 5501#define ITYPE short
3d2e8ceb 5502#define SIZEOF_ITYPE SIZEOF_SHORT
1be6b49c
ML		 5503#include "libguile/num2integral.i.c"
5504
5505#define NUM2INTEGRAL scm_num2ushort
5506#define INTEGRAL2NUM scm_ushort2num
5507#define INTEGRAL2BIG scm_i_ushort2big
ca46fb90 5508#define UNSIGNED 1
1be6b49c 5509#define ITYPE unsigned short
ca46fb90 5510#define SIZEOF_ITYPE SIZEOF_UNSIGNED_SHORT
1be6b49c
ML		 5511#include "libguile/num2integral.i.c"
5512
5513#define NUM2INTEGRAL scm_num2int
5514#define INTEGRAL2NUM scm_int2num
5515#define INTEGRAL2BIG scm_i_int2big
ca46fb90 5516#define UNSIGNED 0
1be6b49c 5517#define ITYPE int
3d2e8ceb 5518#define SIZEOF_ITYPE SIZEOF_INT
1be6b49c
ML		 5519#include "libguile/num2integral.i.c"
5520
5521#define NUM2INTEGRAL scm_num2uint
5522#define INTEGRAL2NUM scm_uint2num
5523#define INTEGRAL2BIG scm_i_uint2big
ca46fb90 5524#define UNSIGNED 1
1be6b49c 5525#define ITYPE unsigned int
ca46fb90 5526#define SIZEOF_ITYPE SIZEOF_UNSIGNED_INT
1be6b49c
ML		 5527#include "libguile/num2integral.i.c"
5528
5529#define NUM2INTEGRAL scm_num2long
5530#define INTEGRAL2NUM scm_long2num
5531#define INTEGRAL2BIG scm_i_long2big
ca46fb90 5532#define UNSIGNED 0
1be6b49c 5533#define ITYPE long
3d2e8ceb 5534#define SIZEOF_ITYPE SIZEOF_LONG
1be6b49c
ML		 5535#include "libguile/num2integral.i.c"
5536
5537#define NUM2INTEGRAL scm_num2ulong
5538#define INTEGRAL2NUM scm_ulong2num
5539#define INTEGRAL2BIG scm_i_ulong2big
ca46fb90 5540#define UNSIGNED 1
1be6b49c 5541#define ITYPE unsigned long
ca46fb90 5542#define SIZEOF_ITYPE SIZEOF_UNSIGNED_LONG
1be6b49c
ML		 5543#include "libguile/num2integral.i.c"
5544
1be6b49c
ML		 5545#define NUM2INTEGRAL scm_num2ptrdiff
5546#define INTEGRAL2NUM scm_ptrdiff2num
5547#define INTEGRAL2BIG scm_i_ptrdiff2big
ca46fb90 5548#define UNSIGNED 0
ee33d62a 5549#define ITYPE scm_t_ptrdiff
3d2e8ceb 5550#define UNSIGNED_ITYPE size_t
ee33d62a 5551#define SIZEOF_ITYPE SCM_SIZEOF_SCM_T_PTRDIFF
1be6b49c
ML		 5552#include "libguile/num2integral.i.c"
5553
5554#define NUM2INTEGRAL scm_num2size
5555#define INTEGRAL2NUM scm_size2num
5556#define INTEGRAL2BIG scm_i_size2big
ca46fb90 5557#define UNSIGNED 1
1be6b49c 5558#define ITYPE size_t
3d2e8ceb 5559#define SIZEOF_ITYPE SIZEOF_SIZE_T
1be6b49c 5560#include "libguile/num2integral.i.c"
0f2d19dd 5561
ee33d62a 5562#if SCM_SIZEOF_LONG_LONG != 0
1cc91f1b 5563
caf08e65
MV		 5564#ifndef ULONG_LONG_MAX
5565#define ULONG_LONG_MAX (~0ULL)
5566#endif
5567
1be6b49c
ML		 5568#define NUM2INTEGRAL scm_num2long_long
5569#define INTEGRAL2NUM scm_long_long2num
5570#define INTEGRAL2BIG scm_i_long_long2big
ca46fb90 5571#define UNSIGNED 0
1be6b49c 5572#define ITYPE long long
3d2e8ceb 5573#define SIZEOF_ITYPE SIZEOF_LONG_LONG
1be6b49c
ML		 5574#include "libguile/num2integral.i.c"
5575
5576#define NUM2INTEGRAL scm_num2ulong_long
5577#define INTEGRAL2NUM scm_ulong_long2num
5578#define INTEGRAL2BIG scm_i_ulong_long2big
ca46fb90 5579#define UNSIGNED 1
1be6b49c 5580#define ITYPE unsigned long long
ca46fb90 5581#define SIZEOF_ITYPE SIZEOF_UNSIGNED_LONG_LONG
1be6b49c 5582#include "libguile/num2integral.i.c"
0f2d19dd 5583
ee33d62a 5584#endif /* SCM_SIZEOF_LONG_LONG != 0 */
caf08e65 5585
5437598b
MD		 5586#define NUM2FLOAT scm_num2float
5587#define FLOAT2NUM scm_float2num
5588#define FTYPE float
5589#include "libguile/num2float.i.c"
5590
5591#define NUM2FLOAT scm_num2double
5592#define FLOAT2NUM scm_double2num
5593#define FTYPE double
5594#include "libguile/num2float.i.c"
5595
1be6b49c 5596#ifdef GUILE_DEBUG
caf08e65 5597
6063dc1d
SJ		 5598#ifndef SIZE_MAX
5599#define SIZE_MAX ((size_t) (-1))
5600#endif
5601#ifndef PTRDIFF_MIN
5602#define PTRDIFF_MIN \
b4fb7de8
RB		 5603 ((scm_t_ptrdiff) ((scm_t_ptrdiff) 1 \
5604 << ((sizeof (scm_t_ptrdiff) * SCM_CHAR_BIT) - 1)))
6063dc1d
SJ		 5605#endif
5606#ifndef PTRDIFF_MAX
5607#define PTRDIFF_MAX (~ PTRDIFF_MIN)
5608#endif
5609
0aacf84e
MD		 5610#define CHECK(type, v) \
5611 do \
5612 { \
5613 if ((v) != scm_num2##type (scm_##type##2num (v), 1, "check_sanity")) \
5614 abort (); \
5615 } \
5616 while (0)
caf08e65 5617
1be6b49c
ML		 5618static void
5619check_sanity ()
5620{
5621 CHECK (short, 0);
5622 CHECK (ushort, 0U);
5623 CHECK (int, 0);
5624 CHECK (uint, 0U);
5625 CHECK (long, 0L);
5626 CHECK (ulong, 0UL);
5627 CHECK (size, 0);
5628 CHECK (ptrdiff, 0);
5629
5630 CHECK (short, -1);
5631 CHECK (int, -1);
5632 CHECK (long, -1L);
5633 CHECK (ptrdiff, -1);
5634
5635 CHECK (short, SHRT_MAX);
5636 CHECK (short, SHRT_MIN);
5637 CHECK (ushort, USHRT_MAX);
5638 CHECK (int, INT_MAX);
5639 CHECK (int, INT_MIN);
5640 CHECK (uint, UINT_MAX);
5641 CHECK (long, LONG_MAX);
5642 CHECK (long, LONG_MIN);
5643 CHECK (ulong, ULONG_MAX);
5644 CHECK (size, SIZE_MAX);
5645 CHECK (ptrdiff, PTRDIFF_MAX);
5646 CHECK (ptrdiff, PTRDIFF_MIN);
0f2d19dd 5647
ee33d62a 5648#if SCM_SIZEOF_LONG_LONG != 0
1be6b49c
ML		 5649 CHECK (long_long, 0LL);
5650 CHECK (ulong_long, 0ULL);
1be6b49c 5651 CHECK (long_long, -1LL);
1be6b49c
ML		 5652 CHECK (long_long, LLONG_MAX);
5653 CHECK (long_long, LLONG_MIN);
5654 CHECK (ulong_long, ULLONG_MAX);
5655#endif
0f2d19dd
JB		 5656}
5657
b10586f0
ML		 5658#undef CHECK
5659
5660#define CHECK \
5661 scm_internal_catch (SCM_BOOL_T, check_body, &data, check_handler, &data); \
5662 if (!SCM_FALSEP (data)) abort();
5663
5664static SCM
5665check_body (void *data)
5666{
5667 SCM num = *(SCM *) data;
5668 scm_num2ulong (num, 1, NULL);
5669
5670 return SCM_UNSPECIFIED;
5671}
5672
5673static SCM
5674check_handler (void *data, SCM tag, SCM throw_args)
5675{
5676 SCM *num = (SCM *) data;
5677 *num = SCM_BOOL_F;
5678
5679 return SCM_UNSPECIFIED;
5680}
5681
5682SCM_DEFINE (scm_sys_check_number_conversions, "%check-number-conversions", 0, 0, 0,
b4e15479 5683 (void),
b10586f0
ML		 5684 "Number conversion sanity checking.")
5685#define FUNC_NAME s_scm_sys_check_number_conversions
5686{
5687 SCM data = SCM_MAKINUM (-1);
5688 CHECK;
5689 data = scm_int2num (INT_MIN);
5690 CHECK;
5691 data = scm_ulong2num (ULONG_MAX);
5692 data = scm_difference (SCM_INUM0, data);
5693 CHECK;
5694 data = scm_ulong2num (ULONG_MAX);
5695 data = scm_sum (SCM_MAKINUM (1), data); data = scm_difference (SCM_INUM0, data);
5696 CHECK;
5697 data = scm_int2num (-10000); data = scm_product (data, data); data = scm_product (data, data);
5698 CHECK;
5699
5700 return SCM_UNSPECIFIED;
5701}
5702#undef FUNC_NAME
5703
1be6b49c 5704#endif
0f2d19dd 5705
0f2d19dd
JB		 5706void
5707scm_init_numbers ()
0f2d19dd 5708{
713a4259
KR		 5709 mpz_init_set_si (z_negative_one, -1);
5710
a261c0e9
DH		 5711 /* It may be possible to tune the performance of some algorithms by using
5712 * the following constants to avoid the creation of bignums. Please, before
5713 * using these values, remember the two rules of program optimization:
5714 * 1st Rule: Don't do it. 2nd Rule (experts only): Don't do it yet. */
86d31dfe
MV		 5715 scm_c_define ("most-positive-fixnum",
5716 SCM_MAKINUM (SCM_MOST_POSITIVE_FIXNUM));
5717 scm_c_define ("most-negative-fixnum",
5718 SCM_MAKINUM (SCM_MOST_NEGATIVE_FIXNUM));
a261c0e9 5719
f3ae5d60
MD		 5720 scm_add_feature ("complex");
5721 scm_add_feature ("inexact");
5986c47d 5722 scm_flo0 = scm_make_real (0.0);
f872b822 5723#ifdef DBL_DIG
0f2d19dd 5724 scm_dblprec = (DBL_DIG > 20) ? 20 : DBL_DIG;
f872b822 5725#else
0f2d19dd
JB		 5726 { /* determine floating point precision */
5727 double f = 0.1;
f872b822 5728 double fsum = 1.0 + f;
0aacf84e
MD		 5729 while (fsum != 1.0)
5730 {
5731 if (++scm_dblprec > 20)
5732 fsum = 1.0;
5733 else
5734 {
5735 f /= 10.0;
5736 fsum = f + 1.0;
5737 }
f872b822
MD		 5738 }
5739 scm_dblprec = scm_dblprec - 1;
0f2d19dd 5740 }
f872b822 5741#endif /* DBL_DIG */
1be6b49c
ML		 5742
5743#ifdef GUILE_DEBUG
5744 check_sanity ();
5745#endif
f92e85f7
MV		 5746
5747 exactly_one_half = scm_permanent_object (scm_divide (SCM_MAKINUM (1),
5748 SCM_MAKINUM (2)));
a0599745 5749#include "libguile/numbers.x"
0f2d19dd 5750}
89e00824
ML		 5751
5752/*
5753 Local Variables:
5754 c-file-style: "gnu"
5755 End:
5756*/
Unnamed repository; edit this file 'description' to name the repository.