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 */
0b799eea 94#define FLOBUFLEN (40+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))
0d75f6d8
KR		 1486 {
1487 /* bits above what's in an inum follow the sign bit */
20fcc8ed 1488 iindex = min (iindex, SCM_LONG_BIT - 1);
0d75f6d8
KR		 1489 return SCM_BOOL ((1L << iindex) & SCM_INUM (j));
1490 }
0aacf84e
MD		 1491 else if (SCM_BIGP (j))
1492 {
1493 int val = mpz_tstbit (SCM_I_BIG_MPZ (j), iindex);
1494 scm_remember_upto_here_1 (j);
1495 return SCM_BOOL (val);
1496 }
1497 else
78166ad5 1498 SCM_WRONG_TYPE_ARG (SCM_ARG2, j);
0f2d19dd 1499}
1bbd0b84 1500#undef FUNC_NAME
0f2d19dd 1501
78166ad5 1502
a1ec6916 1503SCM_DEFINE (scm_lognot, "lognot", 1, 0, 0,
1bbd0b84 1504 (SCM n),
4d814788 1505 "Return the integer which is the ones-complement of the integer\n"
1e6808ea
MG		 1506 "argument.\n"
1507 "\n"
b380b885
MD		 1508 "@lisp\n"
1509 "(number->string (lognot #b10000000) 2)\n"
1510 " @result{} \"-10000001\"\n"
1511 "(number->string (lognot #b0) 2)\n"
1512 " @result{} \"-1\"\n"
1e6808ea 1513 "@end lisp")
1bbd0b84 1514#define FUNC_NAME s_scm_lognot
0f2d19dd 1515{
f9811f9f
KR		 1516 if (SCM_INUMP (n)) {
1517 /* No overflow here, just need to toggle all the bits making up the inum.
1518 Enhancement: No need to strip the tag and add it back, could just xor
1519 a block of 1 bits, if that worked with the various debug versions of
1520 the SCM typedef. */
1521 return SCM_MAKINUM (~ SCM_INUM (n));
1522
1523 } else if (SCM_BIGP (n)) {
1524 SCM result = scm_i_mkbig ();
1525 mpz_com (SCM_I_BIG_MPZ (result), SCM_I_BIG_MPZ (n));
1526 scm_remember_upto_here_1 (n);
1527 return result;
1528
1529 } else {
1530 SCM_WRONG_TYPE_ARG (SCM_ARG1, n);
1531 }
0f2d19dd 1532}
1bbd0b84 1533#undef FUNC_NAME
0f2d19dd 1534
518b7508
KR		 1535/* returns 0 if IN is not an integer. OUT must already be
1536 initialized. */
1537static int
1538coerce_to_big (SCM in, mpz_t out)
1539{
1540 if (SCM_BIGP (in))
1541 mpz_set (out, SCM_I_BIG_MPZ (in));
1542 else if (SCM_INUMP (in))
1543 mpz_set_si (out, SCM_INUM (in));
1544 else
1545 return 0;
1546
1547 return 1;
1548}
1549
d885e204 1550SCM_DEFINE (scm_modulo_expt, "modulo-expt", 3, 0, 0,
518b7508
KR		 1551 (SCM n, SCM k, SCM m),
1552 "Return @var{n} raised to the integer exponent\n"
1553 "@var{k}, modulo @var{m}.\n"
1554 "\n"
1555 "@lisp\n"
1556 "(modulo-expt 2 3 5)\n"
1557 " @result{} 3\n"
1558 "@end lisp")
d885e204 1559#define FUNC_NAME s_scm_modulo_expt
518b7508
KR		 1560{
1561 mpz_t n_tmp;
1562 mpz_t k_tmp;
1563 mpz_t m_tmp;
1564
1565 /* There are two classes of error we might encounter --
1566 1) Math errors, which we'll report by calling scm_num_overflow,
1567 and
1568 2) wrong-type errors, which of course we'll report by calling
1569 SCM_WRONG_TYPE_ARG.
1570 We don't report those errors immediately, however; instead we do
1571 some cleanup first. These variables tell us which error (if
1572 any) we should report after cleaning up.
1573 */
1574 int report_overflow = 0;
1575
1576 int position_of_wrong_type = 0;
1577 SCM value_of_wrong_type = SCM_INUM0;
1578
1579 SCM result = SCM_UNDEFINED;
1580
1581 mpz_init (n_tmp);
1582 mpz_init (k_tmp);
1583 mpz_init (m_tmp);
1584
1585 if (SCM_EQ_P (m, SCM_INUM0))
1586 {
1587 report_overflow = 1;
1588 goto cleanup;
1589 }
1590
1591 if (!coerce_to_big (n, n_tmp))
1592 {
1593 value_of_wrong_type = n;
1594 position_of_wrong_type = 1;
1595 goto cleanup;
1596 }
1597
1598 if (!coerce_to_big (k, k_tmp))
1599 {
1600 value_of_wrong_type = k;
1601 position_of_wrong_type = 2;
1602 goto cleanup;
1603 }
1604
1605 if (!coerce_to_big (m, m_tmp))
1606 {
1607 value_of_wrong_type = m;
1608 position_of_wrong_type = 3;
1609 goto cleanup;
1610 }
1611
1612 /* if the exponent K is negative, and we simply call mpz_powm, we
1613 will get a divide-by-zero exception when an inverse 1/n mod m
1614 doesn't exist (or is not unique). Since exceptions are hard to
1615 handle, we'll attempt the inversion "by hand" -- that way, we get
1616 a simple failure code, which is easy to handle. */
1617
1618 if (-1 == mpz_sgn (k_tmp))
1619 {
1620 if (!mpz_invert (n_tmp, n_tmp, m_tmp))
1621 {
1622 report_overflow = 1;
1623 goto cleanup;
1624 }
1625 mpz_neg (k_tmp, k_tmp);
1626 }
1627
1628 result = scm_i_mkbig ();
1629 mpz_powm (SCM_I_BIG_MPZ (result),
1630 n_tmp,
1631 k_tmp,
1632 m_tmp);
b7b8c575
KR		 1633
1634 if (mpz_sgn (m_tmp) < 0 && mpz_sgn (SCM_I_BIG_MPZ (result)) != 0)
1635 mpz_add (SCM_I_BIG_MPZ (result), SCM_I_BIG_MPZ (result), m_tmp);
1636
518b7508
KR		 1637 cleanup:
1638 mpz_clear (m_tmp);
1639 mpz_clear (k_tmp);
1640 mpz_clear (n_tmp);
1641
1642 if (report_overflow)
1643 scm_num_overflow (FUNC_NAME);
1644
1645 if (position_of_wrong_type)
1646 SCM_WRONG_TYPE_ARG (position_of_wrong_type,
1647 value_of_wrong_type);
1648
1649 return scm_i_normbig (result);
1650}
1651#undef FUNC_NAME
1652
a1ec6916 1653SCM_DEFINE (scm_integer_expt, "integer-expt", 2, 0, 0,
2cd04b42 1654 (SCM n, SCM k),
1e6808ea
MG		 1655 "Return @var{n} raised to the non-negative integer exponent\n"
1656 "@var{k}.\n"
1657 "\n"
b380b885
MD		 1658 "@lisp\n"
1659 "(integer-expt 2 5)\n"
1660 " @result{} 32\n"
1661 "(integer-expt -3 3)\n"
1662 " @result{} -27\n"
1663 "@end lisp")
1bbd0b84 1664#define FUNC_NAME s_scm_integer_expt
0f2d19dd 1665{
1c35cb19
RB		 1666 long i2 = 0;
1667 SCM z_i2 = SCM_BOOL_F;
1668 int i2_is_big = 0;
f872b822 1669 SCM acc = SCM_MAKINUM (1L);
ca46fb90 1670
d57ed702 1671 /* 0^0 == 1 according to R5RS */
4260a7fc 1672 if (SCM_EQ_P (n, SCM_INUM0) || SCM_EQ_P (n, acc))
7b3381f4 1673 return SCM_FALSEP (scm_zero_p(k)) ? n : acc;
4260a7fc
DH		 1674 else if (SCM_EQ_P (n, SCM_MAKINUM (-1L)))
1675 return SCM_FALSEP (scm_even_p (k)) ? n : acc;
ca46fb90 1676
ca46fb90
RB		 1677 if (SCM_INUMP (k))
1678 i2 = SCM_INUM (k);
1679 else if (SCM_BIGP (k))
1680 {
1681 z_i2 = scm_i_clonebig (k, 1);
ca46fb90
RB		 1682 scm_remember_upto_here_1 (k);
1683 i2_is_big = 1;
1684 }
1685 else if (SCM_REALP (k))
2830fd91
MD		 1686 {
1687 double r = SCM_REAL_VALUE (k);
ca46fb90
RB		 1688 if (floor (r) != r)
1689 SCM_WRONG_TYPE_ARG (2, k);
1690 if ((r > SCM_MOST_POSITIVE_FIXNUM) || (r < SCM_MOST_NEGATIVE_FIXNUM))
1691 {
1692 z_i2 = scm_i_mkbig ();
753ac1e7 1693 mpz_set_d (SCM_I_BIG_MPZ (z_i2), r);
ca46fb90
RB		 1694 i2_is_big = 1;
1695 }
1696 else
1697 {
1698 i2 = r;
1699 }
2830fd91
MD		 1700 }
1701 else
ca46fb90
RB		 1702 SCM_WRONG_TYPE_ARG (2, k);
1703
1704 if (i2_is_big)
f872b822 1705 {
ca46fb90
RB		 1706 if (mpz_sgn(SCM_I_BIG_MPZ (z_i2)) == -1)
1707 {
1708 mpz_neg (SCM_I_BIG_MPZ (z_i2), SCM_I_BIG_MPZ (z_i2));
1709 n = scm_divide (n, SCM_UNDEFINED);
1710 }
1711 while (1)
1712 {
1713 if (mpz_sgn(SCM_I_BIG_MPZ (z_i2)) == 0)
1714 {
ca46fb90
RB		 1715 return acc;
1716 }
1717 if (mpz_cmp_ui(SCM_I_BIG_MPZ (z_i2), 1) == 0)
1718 {
ca46fb90
RB		 1719 return scm_product (acc, n);
1720 }
1721 if (mpz_tstbit(SCM_I_BIG_MPZ (z_i2), 0))
1722 acc = scm_product (acc, n);
1723 n = scm_product (n, n);
1724 mpz_fdiv_q_2exp (SCM_I_BIG_MPZ (z_i2), SCM_I_BIG_MPZ (z_i2), 1);
1725 }
f872b822 1726 }
ca46fb90 1727 else
f872b822 1728 {
ca46fb90
RB		 1729 if (i2 < 0)
1730 {
1731 i2 = -i2;
1732 n = scm_divide (n, SCM_UNDEFINED);
1733 }
1734 while (1)
1735 {
1736 if (0 == i2)
1737 return acc;
1738 if (1 == i2)
1739 return scm_product (acc, n);
1740 if (i2 & 1)
1741 acc = scm_product (acc, n);
1742 n = scm_product (n, n);
1743 i2 >>= 1;
1744 }
f872b822 1745 }
0f2d19dd 1746}
1bbd0b84 1747#undef FUNC_NAME
0f2d19dd 1748
a1ec6916 1749SCM_DEFINE (scm_ash, "ash", 2, 0, 0,
1bbd0b84 1750 (SCM n, SCM cnt),
32f19569
KR		 1751 "Return @var{n} shifted left by @var{cnt} bits, or shifted right\n"
1752 "if @var{cnt} is negative. This is an ``arithmetic'' shift.\n"
1e6808ea 1753 "\n"
e7644cb2 1754 "This is effectively a multiplication by 2^@var{cnt}, and when\n"
32f19569
KR		 1755 "@var{cnt} is negative it's a division, rounded towards negative\n"
1756 "infinity. (Note that this is not the same rounding as\n"
1757 "@code{quotient} does.)\n"
1758 "\n"
1759 "With @var{n} viewed as an infinite precision twos complement,\n"
1760 "@code{ash} means a left shift introducing zero bits, or a right\n"
1761 "shift dropping bits.\n"
1e6808ea 1762 "\n"
b380b885 1763 "@lisp\n"
1e6808ea
MG		 1764 "(number->string (ash #b1 3) 2) @result{} \"1000\"\n"
1765 "(number->string (ash #b1010 -1) 2) @result{} \"101\"\n"
32f19569
KR		 1766 "\n"
1767 ";; -23 is bits ...11101001, -6 is bits ...111010\n"
1768 "(ash -23 -2) @result{} -6\n"
a3c8b9fc 1769 "@end lisp")
1bbd0b84 1770#define FUNC_NAME s_scm_ash
0f2d19dd 1771{
3ab9f56e
DH		 1772 long bits_to_shift;
1773
3ab9f56e
DH		 1774 SCM_VALIDATE_INUM (2, cnt);
1775
1776 bits_to_shift = SCM_INUM (cnt);
ca46fb90
RB		 1777
1778 if (bits_to_shift < 0)
1779 {
1780 /* Shift right by abs(cnt) bits. This is realized as a division
1781 by div:=2^abs(cnt). However, to guarantee the floor
1782 rounding, negative values require some special treatment.
1783 */
1784 SCM div = scm_integer_expt (SCM_MAKINUM (2),
1785 SCM_MAKINUM (-bits_to_shift));
f92e85f7
MV		 1786
1787 /* scm_quotient assumes its arguments are integers, but it's legal to (ash 1/2 -1) */
ca46fb90
RB		 1788 if (SCM_FALSEP (scm_negative_p (n)))
1789 return scm_quotient (n, div);
1790 else
1791 return scm_sum (SCM_MAKINUM (-1L),
1792 scm_quotient (scm_sum (SCM_MAKINUM (1L), n), div));
1793 }
1794 else
3ab9f56e 1795 /* Shift left is done by multiplication with 2^CNT */
f872b822 1796 return scm_product (n, scm_integer_expt (SCM_MAKINUM (2), cnt));
0f2d19dd 1797}
1bbd0b84 1798#undef FUNC_NAME
0f2d19dd 1799
3c9f20f8 1800
a1ec6916 1801SCM_DEFINE (scm_bit_extract, "bit-extract", 3, 0, 0,
1bbd0b84 1802 (SCM n, SCM start, SCM end),
1e6808ea
MG		 1803 "Return the integer composed of the @var{start} (inclusive)\n"
1804 "through @var{end} (exclusive) bits of @var{n}. The\n"
1805 "@var{start}th bit becomes the 0-th bit in the result.\n"
1806 "\n"
b380b885
MD		 1807 "@lisp\n"
1808 "(number->string (bit-extract #b1101101010 0 4) 2)\n"
1809 " @result{} \"1010\"\n"
1810 "(number->string (bit-extract #b1101101010 4 9) 2)\n"
1811 " @result{} \"10110\"\n"
1812 "@end lisp")
1bbd0b84 1813#define FUNC_NAME s_scm_bit_extract
0f2d19dd 1814{
7f848242 1815 unsigned long int istart, iend, bits;
34d19ef6 1816 SCM_VALIDATE_INUM_MIN_COPY (2, start,0, istart);
c1bfcf60
GB		 1817 SCM_VALIDATE_INUM_MIN_COPY (3, end, 0, iend);
1818 SCM_ASSERT_RANGE (3, end, (iend >= istart));
78166ad5 1819
7f848242
KR		 1820 /* how many bits to keep */
1821 bits = iend - istart;
1822
0aacf84e
MD		 1823 if (SCM_INUMP (n))
1824 {
1825 long int in = SCM_INUM (n);
7f848242
KR		 1826
1827 /* When istart>=SCM_I_FIXNUM_BIT we can just limit the shift to
d77ad560 1828 SCM_I_FIXNUM_BIT-1 to get either 0 or -1 per the sign of "in". */
857ae6af 1829 in = SCM_SRS (in, min (istart, SCM_I_FIXNUM_BIT-1));
ac0c002c 1830
0aacf84e
MD		 1831 if (in < 0 && bits >= SCM_I_FIXNUM_BIT)
1832 {
1833 /* Since we emulate two's complement encoded numbers, this
1834 * special case requires us to produce a result that has
7f848242 1835 * more bits than can be stored in a fixnum.
0aacf84e 1836 */
7f848242
KR		 1837 SCM result = scm_i_long2big (in);
1838 mpz_fdiv_r_2exp (SCM_I_BIG_MPZ (result), SCM_I_BIG_MPZ (result),
1839 bits);
1840 return result;
0aacf84e 1841 }
ac0c002c 1842
7f848242 1843 /* mask down to requisite bits */
857ae6af 1844 bits = min (bits, SCM_I_FIXNUM_BIT);
7f848242 1845 return SCM_MAKINUM (in & ((1L << bits) - 1));
0aacf84e
MD		 1846 }
1847 else if (SCM_BIGP (n))
ac0c002c 1848 {
7f848242
KR		 1849 SCM result;
1850 if (bits == 1)
1851 {
1852 result = SCM_MAKINUM (mpz_tstbit (SCM_I_BIG_MPZ (n), istart));
1853 }
1854 else
1855 {
1856 /* ENHANCE-ME: It'd be nice not to allocate a new bignum when
1857 bits<SCM_I_FIXNUM_BIT. Would want some help from GMP to get
1858 such bits into a ulong. */
1859 result = scm_i_mkbig ();
1860 mpz_fdiv_q_2exp (SCM_I_BIG_MPZ(result), SCM_I_BIG_MPZ(n), istart);
1861 mpz_fdiv_r_2exp (SCM_I_BIG_MPZ(result), SCM_I_BIG_MPZ(result), bits);
1862 result = scm_i_normbig (result);
1863 }
1864 scm_remember_upto_here_1 (n);
1865 return result;
ac0c002c 1866 }
0aacf84e 1867 else
78166ad5 1868 SCM_WRONG_TYPE_ARG (SCM_ARG1, n);
0f2d19dd 1869}
1bbd0b84 1870#undef FUNC_NAME
0f2d19dd 1871
7f848242 1872
e4755e5c
JB		 1873static const char scm_logtab[] = {
1874 0, 1, 1, 2, 1, 2, 2, 3, 1, 2, 2, 3, 2, 3, 3, 4
1875};
1cc91f1b 1876
a1ec6916 1877SCM_DEFINE (scm_logcount, "logcount", 1, 0, 0,
1bbd0b84 1878 (SCM n),
1e6808ea
MG		 1879 "Return the number of bits in integer @var{n}. If integer is\n"
1880 "positive, the 1-bits in its binary representation are counted.\n"
1881 "If negative, the 0-bits in its two's-complement binary\n"
1882 "representation are counted. If 0, 0 is returned.\n"
1883 "\n"
b380b885
MD		 1884 "@lisp\n"
1885 "(logcount #b10101010)\n"
ca46fb90
RB		 1886 " @result{} 4\n"
1887 "(logcount 0)\n"
1888 " @result{} 0\n"
1889 "(logcount -2)\n"
1890 " @result{} 1\n"
1891 "@end lisp")
1892#define FUNC_NAME s_scm_logcount
1893{
1894 if (SCM_INUMP (n))
f872b822 1895 {
ca46fb90
RB		 1896 unsigned long int c = 0;
1897 long int nn = SCM_INUM (n);
1898 if (nn < 0)
1899 nn = -1 - nn;
1900 while (nn)
1901 {
1902 c += scm_logtab[15 & nn];
1903 nn >>= 4;
1904 }
1905 return SCM_MAKINUM (c);
f872b822 1906 }
ca46fb90 1907 else if (SCM_BIGP (n))
f872b822 1908 {
ca46fb90 1909 unsigned long count;
713a4259
KR		 1910 if (mpz_sgn (SCM_I_BIG_MPZ (n)) >= 0)
1911 count = mpz_popcount (SCM_I_BIG_MPZ (n));
ca46fb90 1912 else
713a4259
KR		 1913 count = mpz_hamdist (SCM_I_BIG_MPZ (n), z_negative_one);
1914 scm_remember_upto_here_1 (n);
ca46fb90 1915 return SCM_MAKINUM (count);
f872b822 1916 }
ca46fb90
RB		 1917 else
1918 SCM_WRONG_TYPE_ARG (SCM_ARG1, n);
0f2d19dd 1919}
ca46fb90 1920#undef FUNC_NAME
0f2d19dd
JB		 1921
1922
ca46fb90
RB		 1923static const char scm_ilentab[] = {
1924 0, 1, 2, 2, 3, 3, 3, 3, 4, 4, 4, 4, 4, 4, 4, 4
1925};
1926
0f2d19dd 1927
ca46fb90
RB		 1928SCM_DEFINE (scm_integer_length, "integer-length", 1, 0, 0,
1929 (SCM n),
1930 "Return the number of bits necessary to represent @var{n}.\n"
1931 "\n"
1932 "@lisp\n"
1933 "(integer-length #b10101010)\n"
1934 " @result{} 8\n"
1935 "(integer-length 0)\n"
1936 " @result{} 0\n"
1937 "(integer-length #b1111)\n"
1938 " @result{} 4\n"
1939 "@end lisp")
1940#define FUNC_NAME s_scm_integer_length
1941{
0aacf84e
MD		 1942 if (SCM_INUMP (n))
1943 {
1944 unsigned long int c = 0;
1945 unsigned int l = 4;
1946 long int nn = SCM_INUM (n);
1947 if (nn < 0)
1948 nn = -1 - nn;
1949 while (nn)
1950 {
1951 c += 4;
1952 l = scm_ilentab [15 & nn];
1953 nn >>= 4;
1954 }
1955 return SCM_MAKINUM (c - 4 + l);
1956 }
1957 else if (SCM_BIGP (n))
1958 {
1959 /* mpz_sizeinbase looks at the absolute value of negatives, whereas we
1960 want a ones-complement. If n is ...111100..00 then mpz_sizeinbase is
1961 1 too big, so check for that and adjust. */
1962 size_t size = mpz_sizeinbase (SCM_I_BIG_MPZ (n), 2);
1963 if (mpz_sgn (SCM_I_BIG_MPZ (n)) < 0
1964 && mpz_scan0 (SCM_I_BIG_MPZ (n), /* no 0 bits above the lowest 1 */
1965 mpz_scan1 (SCM_I_BIG_MPZ (n), 0)) == ULONG_MAX)
1966 size--;
1967 scm_remember_upto_here_1 (n);
1968 return SCM_MAKINUM (size);
1969 }
1970 else
ca46fb90 1971 SCM_WRONG_TYPE_ARG (SCM_ARG1, n);
ca46fb90
RB		 1972}
1973#undef FUNC_NAME
0f2d19dd
JB		 1974
1975/*** NUMBERS -> STRINGS ***/
0b799eea
MV		 1976#define SCM_MAX_DBL_PREC 60
1977#define SCM_MAX_DBL_RADIX 36
1978
1979/* this is an array starting with radix 2, and ending with radix SCM_MAX_DBL_RADIX */
1980static int scm_dblprec[SCM_MAX_DBL_RADIX - 1];
1981static double fx_per_radix[SCM_MAX_DBL_RADIX - 1][SCM_MAX_DBL_PREC];
1982
1983static
1984void init_dblprec(int *prec, int radix) {
1985 /* determine floating point precision by adding successively
1986 smaller increments to 1.0 until it is considered == 1.0 */
1987 double f = ((double)1.0)/radix;
1988 double fsum = 1.0 + f;
1989
1990 *prec = 0;
1991 while (fsum != 1.0)
1992 {
1993 if (++(*prec) > SCM_MAX_DBL_PREC)
1994 fsum = 1.0;
1995 else
1996 {
1997 f /= radix;
1998 fsum = f + 1.0;
1999 }
2000 }
2001 (*prec) -= 1;
2002}
2003
2004static
2005void init_fx_radix(double *fx_list, int radix)
2006{
2007 /* initialize a per-radix list of tolerances. When added
2008 to a number < 1.0, we can determine if we should raund
2009 up and quit converting a number to a string. */
2010 int i;
2011 fx_list[0] = 0.0;
2012 fx_list[1] = 0.5;
2013 for( i = 2 ; i < SCM_MAX_DBL_PREC; ++i )
2014 fx_list[i] = (fx_list[i-1] / radix);
2015}
2016
2017/* use this array as a way to generate a single digit */
2018static const char*number_chars="0123456789ABCDEFGHIJKLMNOPQRSTUVWXYZ";
0f2d19dd 2019
1be6b49c 2020static size_t
0b799eea 2021idbl2str (double f, char *a, int radix)
0f2d19dd 2022{
0b799eea
MV		 2023 int efmt, dpt, d, i, wp;
2024 double *fx;
2025#ifdef DBL_MIN_10_EXP
2026 double f_cpy;
2027 int exp_cpy;
2028#endif /* DBL_MIN_10_EXP */
2029 size_t ch = 0;
2030 int exp = 0;
2031
2032 if(radix < 2 ||
2033 radix > SCM_MAX_DBL_RADIX)
2034 {
2035 /* revert to existing behavior */
2036 radix = 10;
2037 }
2038
2039 wp = scm_dblprec[radix-2];
2040 fx = fx_per_radix[radix-2];
0f2d19dd 2041
f872b822 2042 if (f == 0.0)
abb7e44d
MV		 2043 {
2044#ifdef HAVE_COPYSIGN
2045 double sgn = copysign (1.0, f);
2046
2047 if (sgn < 0.0)
2048 a[ch++] = '-';
2049#endif
abb7e44d
MV		 2050 goto zero; /*{a[0]='0'; a[1]='.'; a[2]='0'; return 3;} */
2051 }
7351e207
MV		 2052
2053 if (xisinf (f))
2054 {
2055 if (f < 0)
2056 strcpy (a, "-inf.0");
2057 else
2058 strcpy (a, "+inf.0");
2059 return ch+6;
2060 }
2061 else if (xisnan (f))
2062 {
2063 strcpy (a, "+nan.0");
2064 return ch+6;
2065 }
2066
f872b822
MD		 2067 if (f < 0.0)
2068 {
2069 f = -f;
2070 a[ch++] = '-';
2071 }
7351e207 2072
f872b822
MD		 2073#ifdef DBL_MIN_10_EXP /* Prevent unnormalized values, as from
2074 make-uniform-vector, from causing infinite loops. */
0b799eea
MV		 2075 /* just do the checking...if it passes, we do the conversion for our
2076 radix again below */
2077 f_cpy = f;
2078 exp_cpy = exp;
2079
2080 while (f_cpy < 1.0)
f872b822 2081 {
0b799eea
MV		 2082 f_cpy *= 10.0;
2083 if (exp_cpy-- < DBL_MIN_10_EXP)
7351e207
MV		 2084 {
2085 a[ch++] = '#';
2086 a[ch++] = '.';
2087 a[ch++] = '#';
2088 return ch;
2089 }
f872b822 2090 }
0b799eea 2091 while (f_cpy > 10.0)
f872b822 2092 {
0b799eea
MV		 2093 f_cpy *= 0.10;
2094 if (exp_cpy++ > DBL_MAX_10_EXP)
7351e207
MV		 2095 {
2096 a[ch++] = '#';
2097 a[ch++] = '.';
2098 a[ch++] = '#';
2099 return ch;
2100 }
f872b822 2101 }
0b799eea
MV		 2102#endif
2103
f872b822
MD		 2104 while (f < 1.0)
2105 {
0b799eea 2106 f *= radix;
f872b822
MD		 2107 exp--;
2108 }
0b799eea 2109 while (f > radix)
f872b822 2110 {
0b799eea 2111 f /= radix;
f872b822
MD		 2112 exp++;
2113 }
0b799eea
MV		 2114
2115 if (f + fx[wp] >= radix)
f872b822
MD		 2116 {
2117 f = 1.0;
2118 exp++;
2119 }
0f2d19dd 2120 zero:
0b799eea
MV		 2121#ifdef ENGNOT
2122 /* adding 9999 makes this equivalent to abs(x) % 3 */
f872b822 2123 dpt = (exp + 9999) % 3;
0f2d19dd
JB		 2124 exp -= dpt++;
2125 efmt = 1;
f872b822
MD		 2126#else
2127 efmt = (exp < -3) || (exp > wp + 2);
0f2d19dd 2128 if (!efmt)
cda139a7
MD		 2129 {
2130 if (exp < 0)
2131 {
2132 a[ch++] = '0';
2133 a[ch++] = '.';
2134 dpt = exp;
f872b822
MD		 2135 while (++dpt)
2136 a[ch++] = '0';
cda139a7
MD		 2137 }
2138 else
f872b822 2139 dpt = exp + 1;
cda139a7 2140 }
0f2d19dd
JB		 2141 else
2142 dpt = 1;
f872b822
MD		 2143#endif
2144
2145 do
2146 {
2147 d = f;
2148 f -= d;
0b799eea 2149 a[ch++] = number_chars[d];
f872b822
MD		 2150 if (f < fx[wp])
2151 break;
2152 if (f + fx[wp] >= 1.0)
2153 {
0b799eea 2154 a[ch - 1] = number_chars[d+1];
f872b822
MD		 2155 break;
2156 }
0b799eea 2157 f *= radix;
f872b822
MD		 2158 if (!(--dpt))
2159 a[ch++] = '.';
0f2d19dd 2160 }
f872b822 2161 while (wp--);
0f2d19dd
JB		 2162
2163 if (dpt > 0)
cda139a7 2164 {
f872b822 2165#ifndef ENGNOT
cda139a7
MD		 2166 if ((dpt > 4) && (exp > 6))
2167 {
f872b822 2168 d = (a[0] == '-' ? 2 : 1);
cda139a7 2169 for (i = ch++; i > d; i--)
f872b822 2170 a[i] = a[i - 1];
cda139a7
MD		 2171 a[d] = '.';
2172 efmt = 1;
2173 }
2174 else
f872b822 2175#endif
cda139a7 2176 {
f872b822
MD		 2177 while (--dpt)
2178 a[ch++] = '0';
cda139a7
MD		 2179 a[ch++] = '.';
2180 }
2181 }
f872b822
MD		 2182 if (a[ch - 1] == '.')
2183 a[ch++] = '0'; /* trailing zero */
2184 if (efmt && exp)
2185 {
2186 a[ch++] = 'e';
2187 if (exp < 0)
2188 {
2189 exp = -exp;
2190 a[ch++] = '-';
2191 }
0b799eea
MV		 2192 for (i = radix; i <= exp; i *= radix);
2193 for (i /= radix; i; i /= radix)
f872b822 2194 {
0b799eea 2195 a[ch++] = number_chars[exp / i];
f872b822
MD		 2196 exp %= i;
2197 }
0f2d19dd 2198 }
0f2d19dd
JB		 2199 return ch;
2200}
2201
1be6b49c 2202static size_t
0b799eea 2203iflo2str (SCM flt, char *str, int radix)
0f2d19dd 2204{
1be6b49c 2205 size_t i;
3c9a524f 2206 if (SCM_REALP (flt))
0b799eea 2207 i = idbl2str (SCM_REAL_VALUE (flt), str, radix);
0f2d19dd 2208 else
f872b822 2209 {
0b799eea 2210 i = idbl2str (SCM_COMPLEX_REAL (flt), str, radix);
f3ae5d60
MD		 2211 if (SCM_COMPLEX_IMAG (flt) != 0.0)
2212 {
7351e207
MV		 2213 double imag = SCM_COMPLEX_IMAG (flt);
2214 /* Don't output a '+' for negative numbers or for Inf and
2215 NaN. They will provide their own sign. */
2216 if (0 <= imag && !xisinf (imag) && !xisnan (imag))
f3ae5d60 2217 str[i++] = '+';
0b799eea 2218 i += idbl2str (imag, &str[i], radix);
f3ae5d60
MD		 2219 str[i++] = 'i';
2220 }
f872b822 2221 }
0f2d19dd
JB		 2222 return i;
2223}
0f2d19dd 2224
5c11cc9d 2225/* convert a long to a string (unterminated). returns the number of
1bbd0b84
GB		 2226 characters in the result.
2227 rad is output base
2228 p is destination: worst case (base 2) is SCM_INTBUFLEN */
1be6b49c 2229size_t
1bbd0b84 2230scm_iint2str (long num, int rad, char *p)
0f2d19dd 2231{
1be6b49c
ML		 2232 size_t j = 1;
2233 size_t i;
5c11cc9d
GH		 2234 unsigned long n = (num < 0) ? -num : num;
2235
f872b822 2236 for (n /= rad; n > 0; n /= rad)
5c11cc9d
GH		 2237 j++;
2238
2239 i = j;
2240 if (num < 0)
f872b822 2241 {
f872b822 2242 *p++ = '-';
5c11cc9d
GH		 2243 j++;
2244 n = -num;
f872b822 2245 }
5c11cc9d
GH		 2246 else
2247 n = num;
f872b822
MD		 2248 while (i--)
2249 {
5c11cc9d
GH		 2250 int d = n % rad;
2251
f872b822
MD		 2252 n /= rad;
2253 p[i] = d + ((d < 10) ? '0' : 'a' - 10);
2254 }
0f2d19dd
JB		 2255 return j;
2256}
2257
a1ec6916 2258SCM_DEFINE (scm_number_to_string, "number->string", 1, 1, 0,
bb628794
DH		 2259 (SCM n, SCM radix),
2260 "Return a string holding the external representation of the\n"
942e5b91
MG		 2261 "number @var{n} in the given @var{radix}. If @var{n} is\n"
2262 "inexact, a radix of 10 will be used.")
1bbd0b84 2263#define FUNC_NAME s_scm_number_to_string
0f2d19dd 2264{
1bbd0b84 2265 int base;
98cb6e75 2266
0aacf84e 2267 if (SCM_UNBNDP (radix))
98cb6e75 2268 base = 10;
0aacf84e
MD		 2269 else
2270 {
2271 SCM_VALIDATE_INUM (2, radix);
2272 base = SCM_INUM (radix);
2273 /* FIXME: ask if range limit was OK, and if so, document */
2274 SCM_ASSERT_RANGE (2, radix, (base >= 2) && (base <= 36));
2275 }
98cb6e75 2276
0aacf84e
MD		 2277 if (SCM_INUMP (n))
2278 {
2279 char num_buf [SCM_INTBUFLEN];
2280 size_t length = scm_iint2str (SCM_INUM (n), base, num_buf);
2281 return scm_mem2string (num_buf, length);
2282 }
2283 else if (SCM_BIGP (n))
2284 {
2285 char *str = mpz_get_str (NULL, base, SCM_I_BIG_MPZ (n));
2286 scm_remember_upto_here_1 (n);
2287 return scm_take0str (str);
2288 }
f92e85f7
MV		 2289 else if (SCM_FRACTIONP (n))
2290 {
2291 scm_i_fraction_reduce (n);
2292 return scm_string_append (scm_list_3 (scm_number_to_string (SCM_FRACTION_NUMERATOR (n), radix),
2293 scm_mem2string ("/", 1),
2294 scm_number_to_string (SCM_FRACTION_DENOMINATOR (n), radix)));
2295 }
0aacf84e
MD		 2296 else if (SCM_INEXACTP (n))
2297 {
2298 char num_buf [FLOBUFLEN];
0b799eea 2299 return scm_mem2string (num_buf, iflo2str (n, num_buf, base));
0aacf84e
MD		 2300 }
2301 else
bb628794 2302 SCM_WRONG_TYPE_ARG (1, n);
0f2d19dd 2303}
1bbd0b84 2304#undef FUNC_NAME
0f2d19dd
JB		 2305
2306
ca46fb90
RB		 2307/* These print routines used to be stubbed here so that scm_repl.c
2308 wouldn't need SCM_BIGDIG conditionals (pre GMP) */
1cc91f1b 2309
0f2d19dd 2310int
e81d98ec 2311scm_print_real (SCM sexp, SCM port, scm_print_state *pstate SCM_UNUSED)
0f2d19dd 2312{
56e55ac7 2313 char num_buf[FLOBUFLEN];
0b799eea 2314 scm_lfwrite (num_buf, iflo2str (sexp, num_buf, 10), port);
0f2d19dd
JB		 2315 return !0;
2316}
2317
f3ae5d60 2318int
e81d98ec 2319scm_print_complex (SCM sexp, SCM port, scm_print_state *pstate SCM_UNUSED)
f92e85f7 2320
f3ae5d60 2321{
56e55ac7 2322 char num_buf[FLOBUFLEN];
0b799eea 2323 scm_lfwrite (num_buf, iflo2str (sexp, num_buf, 10), port);
f3ae5d60
MD		 2324 return !0;
2325}
1cc91f1b 2326
f92e85f7
MV		 2327int
2328scm_i_print_fraction (SCM sexp, SCM port, scm_print_state *pstate SCM_UNUSED)
2329{
2330 SCM str;
2331 scm_i_fraction_reduce (sexp);
2332 str = scm_number_to_string (sexp, SCM_UNDEFINED);
2333 scm_lfwrite (SCM_STRING_CHARS (str), SCM_STRING_LENGTH (str), port);
2334 scm_remember_upto_here_1 (str);
2335 return !0;
2336}
2337
0f2d19dd 2338int
e81d98ec 2339scm_bigprint (SCM exp, SCM port, scm_print_state *pstate SCM_UNUSED)
0f2d19dd 2340{
ca46fb90
RB		 2341 char *str = mpz_get_str (NULL, 10, SCM_I_BIG_MPZ (exp));
2342 scm_remember_upto_here_1 (exp);
2343 scm_lfwrite (str, (size_t) strlen (str), port);
2344 free (str);
0f2d19dd
JB		 2345 return !0;
2346}
2347/*** END nums->strs ***/
2348
3c9a524f 2349
0f2d19dd 2350/*** STRINGS -> NUMBERS ***/
2a8fecee 2351
3c9a524f
DH		 2352/* The following functions implement the conversion from strings to numbers.
2353 * The implementation somehow follows the grammar for numbers as it is given
2354 * in R5RS. Thus, the functions resemble syntactic units (<ureal R>,
2355 * <uinteger R>, ...) that are used to build up numbers in the grammar. Some
2356 * points should be noted about the implementation:
2357 * * Each function keeps a local index variable 'idx' that points at the
2358 * current position within the parsed string. The global index is only
2359 * updated if the function could parse the corresponding syntactic unit
2360 * successfully.
2361 * * Similarly, the functions keep track of indicators of inexactness ('#',
2362 * '.' or exponents) using local variables ('hash_seen', 'x'). Again, the
2363 * global exactness information is only updated after each part has been
2364 * successfully parsed.
2365 * * Sequences of digits are parsed into temporary variables holding fixnums.
2366 * Only if these fixnums would overflow, the result variables are updated
2367 * using the standard functions scm_add, scm_product, scm_divide etc. Then,
2368 * the temporary variables holding the fixnums are cleared, and the process
2369 * starts over again. If for example fixnums were able to store five decimal
2370 * digits, a number 1234567890 would be parsed in two parts 12345 and 67890,
2371 * and the result was computed as 12345 * 100000 + 67890. In other words,
2372 * only every five digits two bignum operations were performed.
2373 */
2374
2375enum t_exactness {NO_EXACTNESS, INEXACT, EXACT};
2376
2377/* R5RS, section 7.1.1, lexical structure of numbers: <uinteger R>. */
2378
2379/* In non ASCII-style encodings the following macro might not work. */
71df73ac
KR		 2380#define XDIGIT2UINT(d) \
2381 (isdigit ((int) (unsigned char) d) \
2382 ? (d) - '0' \
2383 : tolower ((int) (unsigned char) d) - 'a' + 10)
3c9a524f 2384
2a8fecee 2385static SCM
3c9a524f
DH		 2386mem2uinteger (const char* mem, size_t len, unsigned int *p_idx,
2387 unsigned int radix, enum t_exactness *p_exactness)
2a8fecee 2388{
3c9a524f
DH		 2389 unsigned int idx = *p_idx;
2390 unsigned int hash_seen = 0;
2391 scm_t_bits shift = 1;
2392 scm_t_bits add = 0;
2393 unsigned int digit_value;
2394 SCM result;
2395 char c;
2396
2397 if (idx == len)
2398 return SCM_BOOL_F;
2a8fecee 2399
3c9a524f 2400 c = mem[idx];
71df73ac 2401 if (!isxdigit ((int) (unsigned char) c))
3c9a524f
DH		 2402 return SCM_BOOL_F;
2403 digit_value = XDIGIT2UINT (c);
2404 if (digit_value >= radix)
2405 return SCM_BOOL_F;
2406
2407 idx++;
2408 result = SCM_MAKINUM (digit_value);
2409 while (idx != len)
f872b822 2410 {
3c9a524f 2411 char c = mem[idx];
71df73ac 2412 if (isxdigit ((int) (unsigned char) c))
f872b822 2413 {
3c9a524f 2414 if (hash_seen)
1fe5e088 2415 break;
3c9a524f
DH		 2416 digit_value = XDIGIT2UINT (c);
2417 if (digit_value >= radix)
1fe5e088 2418 break;
f872b822 2419 }
3c9a524f
DH		 2420 else if (c == '#')
2421 {
2422 hash_seen = 1;
2423 digit_value = 0;
2424 }
2425 else
2426 break;
2427
2428 idx++;
2429 if (SCM_MOST_POSITIVE_FIXNUM / radix < shift)
2430 {
2431 result = scm_product (result, SCM_MAKINUM (shift));
2432 if (add > 0)
2433 result = scm_sum (result, SCM_MAKINUM (add));
2434
2435 shift = radix;
2436 add = digit_value;
2437 }
2438 else
2439 {
2440 shift = shift * radix;
2441 add = add * radix + digit_value;
2442 }
2443 };
2444
2445 if (shift > 1)
2446 result = scm_product (result, SCM_MAKINUM (shift));
2447 if (add > 0)
2448 result = scm_sum (result, SCM_MAKINUM (add));
2449
2450 *p_idx = idx;
2451 if (hash_seen)
2452 *p_exactness = INEXACT;
2453
2454 return result;
2a8fecee
JB		 2455}
2456
2457
3c9a524f
DH		 2458/* R5RS, section 7.1.1, lexical structure of numbers: <decimal 10>. Only
2459 * covers the parts of the rules that start at a potential point. The value
2460 * of the digits up to the point have been parsed by the caller and are given
79d34f68
DH		 2461 * in variable result. The content of *p_exactness indicates, whether a hash
2462 * has already been seen in the digits before the point.
3c9a524f 2463 */
1cc91f1b 2464
3c9a524f
DH		 2465/* In non ASCII-style encodings the following macro might not work. */
2466#define DIGIT2UINT(d) ((d) - '0')
2467
2468static SCM
79d34f68 2469mem2decimal_from_point (SCM result, const char* mem, size_t len,
3c9a524f 2470 unsigned int *p_idx, enum t_exactness *p_exactness)
0f2d19dd 2471{
3c9a524f
DH		 2472 unsigned int idx = *p_idx;
2473 enum t_exactness x = *p_exactness;
3c9a524f
DH		 2474
2475 if (idx == len)
79d34f68 2476 return result;
3c9a524f
DH		 2477
2478 if (mem[idx] == '.')
2479 {
2480 scm_t_bits shift = 1;
2481 scm_t_bits add = 0;
2482 unsigned int digit_value;
79d34f68 2483 SCM big_shift = SCM_MAKINUM (1);
3c9a524f
DH		 2484
2485 idx++;
2486 while (idx != len)
2487 {
2488 char c = mem[idx];
71df73ac 2489 if (isdigit ((int) (unsigned char) c))
3c9a524f
DH		 2490 {
2491 if (x == INEXACT)
2492 return SCM_BOOL_F;
2493 else
2494 digit_value = DIGIT2UINT (c);
2495 }
2496 else if (c == '#')
2497 {
2498 x = INEXACT;
2499 digit_value = 0;
2500 }
2501 else
2502 break;
2503
2504 idx++;
2505 if (SCM_MOST_POSITIVE_FIXNUM / 10 < shift)
2506 {
2507 big_shift = scm_product (big_shift, SCM_MAKINUM (shift));
79d34f68 2508 result = scm_product (result, SCM_MAKINUM (shift));
3c9a524f 2509 if (add > 0)
79d34f68 2510 result = scm_sum (result, SCM_MAKINUM (add));
3c9a524f
DH		 2511
2512 shift = 10;
2513 add = digit_value;
2514 }
2515 else
2516 {
2517 shift = shift * 10;
2518 add = add * 10 + digit_value;
2519 }
2520 };
2521
2522 if (add > 0)
2523 {
2524 big_shift = scm_product (big_shift, SCM_MAKINUM (shift));
79d34f68
DH		 2525 result = scm_product (result, SCM_MAKINUM (shift));
2526 result = scm_sum (result, SCM_MAKINUM (add));
3c9a524f
DH		 2527 }
2528
d8592269 2529 result = scm_divide (result, big_shift);
79d34f68 2530
3c9a524f
DH		 2531 /* We've seen a decimal point, thus the value is implicitly inexact. */
2532 x = INEXACT;
f872b822 2533 }
3c9a524f 2534
3c9a524f 2535 if (idx != len)
f872b822 2536 {
3c9a524f
DH		 2537 int sign = 1;
2538 unsigned int start;
2539 char c;
2540 int exponent;
2541 SCM e;
2542
2543 /* R5RS, section 7.1.1, lexical structure of numbers: <suffix> */
2544
2545 switch (mem[idx])
f872b822 2546 {
3c9a524f
DH		 2547 case 'd': case 'D':
2548 case 'e': case 'E':
2549 case 'f': case 'F':
2550 case 'l': case 'L':
2551 case 's': case 'S':
2552 idx++;
2553 start = idx;
2554 c = mem[idx];
2555 if (c == '-')
2556 {
2557 idx++;
2558 sign = -1;
2559 c = mem[idx];
2560 }
2561 else if (c == '+')
2562 {
2563 idx++;
2564 sign = 1;
2565 c = mem[idx];
2566 }
2567 else
2568 sign = 1;
2569
71df73ac 2570 if (!isdigit ((int) (unsigned char) c))
3c9a524f
DH		 2571 return SCM_BOOL_F;
2572
2573 idx++;
2574 exponent = DIGIT2UINT (c);
2575 while (idx != len)
f872b822 2576 {
3c9a524f 2577 char c = mem[idx];
71df73ac 2578 if (isdigit ((int) (unsigned char) c))
3c9a524f
DH		 2579 {
2580 idx++;
2581 if (exponent <= SCM_MAXEXP)
2582 exponent = exponent * 10 + DIGIT2UINT (c);
2583 }
2584 else
2585 break;
f872b822 2586 }
3c9a524f
DH		 2587
2588 if (exponent > SCM_MAXEXP)
f872b822 2589 {
3c9a524f
DH		 2590 size_t exp_len = idx - start;
2591 SCM exp_string = scm_mem2string (&mem[start], exp_len);
2592 SCM exp_num = scm_string_to_number (exp_string, SCM_UNDEFINED);
2593 scm_out_of_range ("string->number", exp_num);
f872b822 2594 }
3c9a524f
DH		 2595
2596 e = scm_integer_expt (SCM_MAKINUM (10), SCM_MAKINUM (exponent));
2597 if (sign == 1)
2598 result = scm_product (result, e);
2599 else
f92e85f7 2600 result = scm_divide2real (result, e);
3c9a524f
DH		 2601
2602 /* We've seen an exponent, thus the value is implicitly inexact. */
2603 x = INEXACT;
2604
f872b822 2605 break;
3c9a524f 2606
f872b822 2607 default:
3c9a524f 2608 break;
f872b822 2609 }
0f2d19dd 2610 }
3c9a524f
DH		 2611
2612 *p_idx = idx;
2613 if (x == INEXACT)
2614 *p_exactness = x;
2615
2616 return result;
0f2d19dd 2617}
0f2d19dd 2618
3c9a524f
DH		 2619
2620/* R5RS, section 7.1.1, lexical structure of numbers: <ureal R> */
2621
2622static SCM
2623mem2ureal (const char* mem, size_t len, unsigned int *p_idx,
2624 unsigned int radix, enum t_exactness *p_exactness)
0f2d19dd 2625{
3c9a524f 2626 unsigned int idx = *p_idx;
164d2481 2627 SCM result;
3c9a524f
DH		 2628
2629 if (idx == len)
2630 return SCM_BOOL_F;
2631
7351e207
MV		 2632 if (idx+5 <= len && !strncmp (mem+idx, "inf.0", 5))
2633 {
2634 *p_idx = idx+5;
2635 return scm_inf ();
2636 }
2637
2638 if (idx+4 < len && !strncmp (mem+idx, "nan.", 4))
2639 {
2640 enum t_exactness x = EXACT;
2641
d8592269
MV		 2642 /* Cobble up the fractional part. We might want to set the
2643 NaN's mantissa from it. */
7351e207
MV		 2644 idx += 4;
2645 mem2uinteger (mem, len, &idx, 10, &x);
2646 *p_idx = idx;
2647 return scm_nan ();
2648 }
2649
3c9a524f
DH		 2650 if (mem[idx] == '.')
2651 {
2652 if (radix != 10)
2653 return SCM_BOOL_F;
2654 else if (idx + 1 == len)
2655 return SCM_BOOL_F;
71df73ac 2656 else if (!isdigit ((int) (unsigned char) mem[idx + 1]))
3c9a524f
DH		 2657 return SCM_BOOL_F;
2658 else
164d2481
MV		 2659 result = mem2decimal_from_point (SCM_MAKINUM (0), mem, len,
2660 p_idx, p_exactness);
f872b822 2661 }
3c9a524f
DH		 2662 else
2663 {
2664 enum t_exactness x = EXACT;
2665 SCM uinteger;
3c9a524f
DH		 2666
2667 uinteger = mem2uinteger (mem, len, &idx, radix, &x);
2668 if (SCM_FALSEP (uinteger))
2669 return SCM_BOOL_F;
2670
2671 if (idx == len)
2672 result = uinteger;
2673 else if (mem[idx] == '/')
f872b822 2674 {
3c9a524f
DH		 2675 SCM divisor;
2676
2677 idx++;
2678
2679 divisor = mem2uinteger (mem, len, &idx, radix, &x);
2680 if (SCM_FALSEP (divisor))
2681 return SCM_BOOL_F;
2682
f92e85f7
MV		 2683 /* both are int/big here, I assume */
2684 result = scm_make_ratio (uinteger, divisor);
f872b822 2685 }
3c9a524f
DH		 2686 else if (radix == 10)
2687 {
2688 result = mem2decimal_from_point (uinteger, mem, len, &idx, &x);
2689 if (SCM_FALSEP (result))
2690 return SCM_BOOL_F;
2691 }
2692 else
2693 result = uinteger;
2694
2695 *p_idx = idx;
2696 if (x == INEXACT)
2697 *p_exactness = x;
f872b822 2698 }
164d2481
MV		 2699
2700 /* When returning an inexact zero, make sure it is represented as a
2701 floating point value so that we can change its sign.
2702 */
2703 if (SCM_EQ_P (result, SCM_MAKINUM(0)) && *p_exactness == INEXACT)
2704 result = scm_make_real (0.0);
2705
2706 return result;
3c9a524f 2707}
0f2d19dd 2708
0f2d19dd 2709
3c9a524f 2710/* R5RS, section 7.1.1, lexical structure of numbers: <complex R> */
0f2d19dd 2711
3c9a524f
DH		 2712static SCM
2713mem2complex (const char* mem, size_t len, unsigned int idx,
2714 unsigned int radix, enum t_exactness *p_exactness)
2715{
2716 char c;
2717 int sign = 0;
2718 SCM ureal;
2719
2720 if (idx == len)
2721 return SCM_BOOL_F;
2722
2723 c = mem[idx];
2724 if (c == '+')
2725 {
2726 idx++;
2727 sign = 1;
2728 }
2729 else if (c == '-')
2730 {
2731 idx++;
2732 sign = -1;
0f2d19dd 2733 }
0f2d19dd 2734
3c9a524f
DH		 2735 if (idx == len)
2736 return SCM_BOOL_F;
2737
2738 ureal = mem2ureal (mem, len, &idx, radix, p_exactness);
2739 if (SCM_FALSEP (ureal))
f872b822 2740 {
3c9a524f
DH		 2741 /* input must be either +i or -i */
2742
2743 if (sign == 0)
2744 return SCM_BOOL_F;
2745
2746 if (mem[idx] == 'i' || mem[idx] == 'I')
f872b822 2747 {
3c9a524f
DH		 2748 idx++;
2749 if (idx != len)
2750 return SCM_BOOL_F;
2751
2752 return scm_make_rectangular (SCM_MAKINUM (0), SCM_MAKINUM (sign));
f872b822 2753 }
3c9a524f
DH		 2754 else
2755 return SCM_BOOL_F;
0f2d19dd 2756 }
3c9a524f
DH		 2757 else
2758 {
fc194577 2759 if (sign == -1 && SCM_FALSEP (scm_nan_p (ureal)))
3c9a524f 2760 ureal = scm_difference (ureal, SCM_UNDEFINED);
f872b822 2761
3c9a524f
DH		 2762 if (idx == len)
2763 return ureal;
2764
2765 c = mem[idx];
2766 switch (c)
f872b822 2767 {
3c9a524f
DH		 2768 case 'i': case 'I':
2769 /* either +<ureal>i or -<ureal>i */
2770
2771 idx++;
2772 if (sign == 0)
2773 return SCM_BOOL_F;
2774 if (idx != len)
2775 return SCM_BOOL_F;
2776 return scm_make_rectangular (SCM_MAKINUM (0), ureal);
2777
2778 case '@':
2779 /* polar input: <real>@<real>. */
2780
2781 idx++;
2782 if (idx == len)
2783 return SCM_BOOL_F;
2784 else
f872b822 2785 {
3c9a524f
DH		 2786 int sign;
2787 SCM angle;
2788 SCM result;
2789
2790 c = mem[idx];
2791 if (c == '+')
2792 {
2793 idx++;
2794 sign = 1;
2795 }
2796 else if (c == '-')
2797 {
2798 idx++;
2799 sign = -1;
2800 }
2801 else
2802 sign = 1;
2803
2804 angle = mem2ureal (mem, len, &idx, radix, p_exactness);
2805 if (SCM_FALSEP (angle))
2806 return SCM_BOOL_F;
2807 if (idx != len)
2808 return SCM_BOOL_F;
2809
fc194577 2810 if (sign == -1 && SCM_FALSEP (scm_nan_p (ureal)))
3c9a524f
DH		 2811 angle = scm_difference (angle, SCM_UNDEFINED);
2812
2813 result = scm_make_polar (ureal, angle);
2814 return result;
f872b822 2815 }
3c9a524f
DH		 2816 case '+':
2817 case '-':
2818 /* expecting input matching <real>[+-]<ureal>?i */
0f2d19dd 2819
3c9a524f
DH		 2820 idx++;
2821 if (idx == len)
2822 return SCM_BOOL_F;
2823 else
2824 {
2825 int sign = (c == '+') ? 1 : -1;
2826 SCM imag = mem2ureal (mem, len, &idx, radix, p_exactness);
0f2d19dd 2827
3c9a524f
DH		 2828 if (SCM_FALSEP (imag))
2829 imag = SCM_MAKINUM (sign);
fc194577 2830 else if (sign == -1 && SCM_FALSEP (scm_nan_p (ureal)))
1fe5e088 2831 imag = scm_difference (imag, SCM_UNDEFINED);
0f2d19dd 2832
3c9a524f
DH		 2833 if (idx == len)
2834 return SCM_BOOL_F;
2835 if (mem[idx] != 'i' && mem[idx] != 'I')
2836 return SCM_BOOL_F;
0f2d19dd 2837
3c9a524f
DH		 2838 idx++;
2839 if (idx != len)
2840 return SCM_BOOL_F;
0f2d19dd 2841
1fe5e088 2842 return scm_make_rectangular (ureal, imag);
3c9a524f
DH		 2843 }
2844 default:
2845 return SCM_BOOL_F;
2846 }
2847 }
0f2d19dd 2848}
0f2d19dd
JB		 2849
2850
3c9a524f
DH		 2851/* R5RS, section 7.1.1, lexical structure of numbers: <number> */
2852
2853enum t_radix {NO_RADIX=0, DUAL=2, OCT=8, DEC=10, HEX=16};
1cc91f1b 2854
0f2d19dd 2855SCM
3c9a524f 2856scm_i_mem2number (const char* mem, size_t len, unsigned int default_radix)
0f2d19dd 2857{
3c9a524f
DH		 2858 unsigned int idx = 0;
2859 unsigned int radix = NO_RADIX;
2860 enum t_exactness forced_x = NO_EXACTNESS;
2861 enum t_exactness implicit_x = EXACT;
2862 SCM result;
2863
2864 /* R5RS, section 7.1.1, lexical structure of numbers: <prefix R> */
2865 while (idx + 2 < len && mem[idx] == '#')
2866 {
2867 switch (mem[idx + 1])
2868 {
2869 case 'b': case 'B':
2870 if (radix != NO_RADIX)
2871 return SCM_BOOL_F;
2872 radix = DUAL;
2873 break;
2874 case 'd': case 'D':
2875 if (radix != NO_RADIX)
2876 return SCM_BOOL_F;
2877 radix = DEC;
2878 break;
2879 case 'i': case 'I':
2880 if (forced_x != NO_EXACTNESS)
2881 return SCM_BOOL_F;
2882 forced_x = INEXACT;
2883 break;
2884 case 'e': case 'E':
2885 if (forced_x != NO_EXACTNESS)
2886 return SCM_BOOL_F;
2887 forced_x = EXACT;
2888 break;
2889 case 'o': case 'O':
2890 if (radix != NO_RADIX)
2891 return SCM_BOOL_F;
2892 radix = OCT;
2893 break;
2894 case 'x': case 'X':
2895 if (radix != NO_RADIX)
2896 return SCM_BOOL_F;
2897 radix = HEX;
2898 break;
2899 default:
f872b822 2900 return SCM_BOOL_F;
3c9a524f
DH		 2901 }
2902 idx += 2;
2903 }
2904
2905 /* R5RS, section 7.1.1, lexical structure of numbers: <complex R> */
2906 if (radix == NO_RADIX)
2907 result = mem2complex (mem, len, idx, default_radix, &implicit_x);
2908 else
2909 result = mem2complex (mem, len, idx, (unsigned int) radix, &implicit_x);
2910
2911 if (SCM_FALSEP (result))
2912 return SCM_BOOL_F;
f872b822 2913
3c9a524f 2914 switch (forced_x)
f872b822 2915 {
3c9a524f
DH		 2916 case EXACT:
2917 if (SCM_INEXACTP (result))
3c9a524f
DH		 2918 return scm_inexact_to_exact (result);
2919 else
2920 return result;
2921 case INEXACT:
2922 if (SCM_INEXACTP (result))
2923 return result;
2924 else
2925 return scm_exact_to_inexact (result);
2926 case NO_EXACTNESS:
2927 default:
2928 if (implicit_x == INEXACT)
2929 {
2930 if (SCM_INEXACTP (result))
2931 return result;
2932 else
2933 return scm_exact_to_inexact (result);
2934 }
2935 else
2936 return result;
f872b822 2937 }
0f2d19dd
JB		 2938}
2939
2940
a1ec6916 2941SCM_DEFINE (scm_string_to_number, "string->number", 1, 1, 0,
bb628794 2942 (SCM string, SCM radix),
1e6808ea 2943 "Return a number of the maximally precise representation\n"
942e5b91 2944 "expressed by the given @var{string}. @var{radix} must be an\n"
5352393c
MG		 2945 "exact integer, either 2, 8, 10, or 16. If supplied, @var{radix}\n"
2946 "is a default radix that may be overridden by an explicit radix\n"
2947 "prefix in @var{string} (e.g. \"#o177\"). If @var{radix} is not\n"
2948 "supplied, then the default radix is 10. If string is not a\n"
2949 "syntactically valid notation for a number, then\n"
2950 "@code{string->number} returns @code{#f}.")
1bbd0b84 2951#define FUNC_NAME s_scm_string_to_number
0f2d19dd
JB		 2952{
2953 SCM answer;
1bbd0b84 2954 int base;
a6d9e5ab 2955 SCM_VALIDATE_STRING (1, string);
34d19ef6 2956 SCM_VALIDATE_INUM_MIN_DEF_COPY (2, radix,2,10, base);
3c9a524f 2957 answer = scm_i_mem2number (SCM_STRING_CHARS (string),
d8592269
MV		 2958 SCM_STRING_LENGTH (string),
2959 base);
bb628794 2960 return scm_return_first (answer, string);
0f2d19dd 2961}
1bbd0b84 2962#undef FUNC_NAME
3c9a524f
DH		 2963
2964
0f2d19dd
JB		 2965/*** END strs->nums ***/
2966
5986c47d 2967
0f2d19dd 2968SCM
f3ae5d60 2969scm_make_real (double x)
0f2d19dd 2970{
3553e1d1
GH		 2971 SCM z = scm_double_cell (scm_tc16_real, 0, 0, 0);
2972
3a9809df 2973 SCM_REAL_VALUE (z) = x;
0f2d19dd
JB		 2974 return z;
2975}
0f2d19dd 2976
5986c47d 2977
f3ae5d60
MD		 2978SCM
2979scm_make_complex (double x, double y)
2980{
0aacf84e 2981 if (y == 0.0)
3a9809df 2982 return scm_make_real (x);
0aacf84e
MD		 2983 else
2984 {
2985 SCM z;
29c4382a 2986 SCM_NEWSMOB (z, scm_tc16_complex, scm_gc_malloc (sizeof (scm_t_complex),
0aacf84e
MD		 2987 "complex"));
2988 SCM_COMPLEX_REAL (z) = x;
2989 SCM_COMPLEX_IMAG (z) = y;
2990 return z;
2991 }
f3ae5d60 2992}
1cc91f1b 2993
5986c47d 2994
0f2d19dd 2995SCM
1bbd0b84 2996scm_bigequal (SCM x, SCM y)
0f2d19dd 2997{
47ae1f0e 2998 int result = mpz_cmp (SCM_I_BIG_MPZ (x), SCM_I_BIG_MPZ (y));
ca46fb90
RB		 2999 scm_remember_upto_here_2 (x, y);
3000 return SCM_BOOL (0 == result);
0f2d19dd
JB		 3001}
3002
0f2d19dd 3003SCM
f3ae5d60 3004scm_real_equalp (SCM x, SCM y)
0f2d19dd 3005{
f3ae5d60 3006 return SCM_BOOL (SCM_REAL_VALUE (x) == SCM_REAL_VALUE (y));
0f2d19dd
JB		 3007}
3008
f3ae5d60
MD		 3009SCM
3010scm_complex_equalp (SCM x, SCM y)
3011{
3012 return SCM_BOOL (SCM_COMPLEX_REAL (x) == SCM_COMPLEX_REAL (y)
3013 && SCM_COMPLEX_IMAG (x) == SCM_COMPLEX_IMAG (y));
3014}
0f2d19dd 3015
f92e85f7
MV		 3016SCM
3017scm_i_fraction_equalp (SCM x, SCM y)
3018{
3019 scm_i_fraction_reduce (x);
3020 scm_i_fraction_reduce (y);
02164269
MV		 3021 if (SCM_FALSEP (scm_equal_p (SCM_FRACTION_NUMERATOR (x),
3022 SCM_FRACTION_NUMERATOR (y)))
3023 || SCM_FALSEP (scm_equal_p (SCM_FRACTION_DENOMINATOR (x),
3024 SCM_FRACTION_DENOMINATOR (y))))
3025 return SCM_BOOL_F;
3026 else
3027 return SCM_BOOL_T;
f92e85f7 3028}
0f2d19dd
JB		 3029
3030
1bbd0b84 3031SCM_REGISTER_PROC (s_number_p, "number?", 1, 0, 0, scm_number_p);
942e5b91
MG		 3032/* "Return @code{#t} if @var{x} is a number, @code{#f}\n"
3033 * "else. Note that the sets of complex, real, rational and\n"
3034 * "integer values form subsets of the set of numbers, i. e. the\n"
3035 * "predicate will be fulfilled for any number."
3036 */
a1ec6916 3037SCM_DEFINE (scm_number_p, "complex?", 1, 0, 0,
1bbd0b84 3038 (SCM x),
942e5b91 3039 "Return @code{#t} if @var{x} is a complex number, @code{#f}\n"
bb2c02f2 3040 "otherwise. Note that the sets of real, rational and integer\n"
942e5b91
MG		 3041 "values form subsets of the set of complex numbers, i. e. the\n"
3042 "predicate will also be fulfilled if @var{x} is a real,\n"
3043 "rational or integer number.")
1bbd0b84 3044#define FUNC_NAME s_scm_number_p
0f2d19dd 3045{
bb628794 3046 return SCM_BOOL (SCM_NUMBERP (x));
0f2d19dd 3047}
1bbd0b84 3048#undef FUNC_NAME
0f2d19dd
JB		 3049
3050
f92e85f7
MV		 3051SCM_DEFINE (scm_real_p, "real?", 1, 0, 0,
3052 (SCM x),
3053 "Return @code{#t} if @var{x} is a real number, @code{#f}\n"
3054 "otherwise. Note that the set of integer values forms a subset of\n"
3055 "the set of real numbers, i. e. the predicate will also be\n"
3056 "fulfilled if @var{x} is an integer number.")
3057#define FUNC_NAME s_scm_real_p
3058{
3059 /* we can't represent irrational numbers. */
3060 return scm_rational_p (x);
3061}
3062#undef FUNC_NAME
3063
3064SCM_DEFINE (scm_rational_p, "rational?", 1, 0, 0,
1bbd0b84 3065 (SCM x),
942e5b91 3066 "Return @code{#t} if @var{x} is a rational number, @code{#f}\n"
bb2c02f2 3067 "otherwise. Note that the set of integer values forms a subset of\n"
942e5b91 3068 "the set of rational numbers, i. e. the predicate will also be\n"
f92e85f7
MV		 3069 "fulfilled if @var{x} is an integer number.")
3070#define FUNC_NAME s_scm_rational_p
0f2d19dd 3071{
0aacf84e 3072 if (SCM_INUMP (x))
0f2d19dd 3073 return SCM_BOOL_T;
0aacf84e 3074 else if (SCM_IMP (x))
0f2d19dd 3075 return SCM_BOOL_F;
0aacf84e 3076 else if (SCM_BIGP (x))
0f2d19dd 3077 return SCM_BOOL_T;
f92e85f7
MV		 3078 else if (SCM_FRACTIONP (x))
3079 return SCM_BOOL_T;
3080 else if (SCM_REALP (x))
3081 /* due to their limited precision, all floating point numbers are
3082 rational as well. */
3083 return SCM_BOOL_T;
0aacf84e 3084 else
bb628794 3085 return SCM_BOOL_F;
0f2d19dd 3086}
1bbd0b84 3087#undef FUNC_NAME
0f2d19dd
JB		 3088
3089
a1ec6916 3090SCM_DEFINE (scm_integer_p, "integer?", 1, 0, 0,
1bbd0b84 3091 (SCM x),
942e5b91
MG		 3092 "Return @code{#t} if @var{x} is an integer number, @code{#f}\n"
3093 "else.")
1bbd0b84 3094#define FUNC_NAME s_scm_integer_p
0f2d19dd
JB		 3095{
3096 double r;
f872b822
MD		 3097 if (SCM_INUMP (x))
3098 return SCM_BOOL_T;
3099 if (SCM_IMP (x))
3100 return SCM_BOOL_F;
f872b822
MD		 3101 if (SCM_BIGP (x))
3102 return SCM_BOOL_T;
3c9a524f 3103 if (!SCM_INEXACTP (x))
f872b822 3104 return SCM_BOOL_F;
3c9a524f 3105 if (SCM_COMPLEXP (x))
f872b822 3106 return SCM_BOOL_F;
5986c47d 3107 r = SCM_REAL_VALUE (x);
f872b822
MD		 3108 if (r == floor (r))
3109 return SCM_BOOL_T;
0f2d19dd
JB		 3110 return SCM_BOOL_F;
3111}
1bbd0b84 3112#undef FUNC_NAME
0f2d19dd
JB		 3113
3114
a1ec6916 3115SCM_DEFINE (scm_inexact_p, "inexact?", 1, 0, 0,
1bbd0b84 3116 (SCM x),
942e5b91
MG		 3117 "Return @code{#t} if @var{x} is an inexact number, @code{#f}\n"
3118 "else.")
1bbd0b84 3119#define FUNC_NAME s_scm_inexact_p
0f2d19dd 3120{
eb927cb9
MV		 3121 if (SCM_INEXACTP (x))
3122 return SCM_BOOL_T;
3123 if (SCM_NUMBERP (x))
3124 return SCM_BOOL_F;
3125 SCM_WRONG_TYPE_ARG (1, x);
0f2d19dd 3126}
1bbd0b84 3127#undef FUNC_NAME
0f2d19dd
JB		 3128
3129
152f82bf 3130SCM_GPROC1 (s_eq_p, "=", scm_tc7_rpsubr, scm_num_eq_p, g_eq_p);
942e5b91 3131/* "Return @code{#t} if all parameters are numerically equal." */
0f2d19dd 3132SCM
6e8d25a6 3133scm_num_eq_p (SCM x, SCM y)
0f2d19dd 3134{
d8b95e27 3135 again:
0aacf84e
MD		 3136 if (SCM_INUMP (x))
3137 {
3138 long xx = SCM_INUM (x);
3139 if (SCM_INUMP (y))
3140 {
3141 long yy = SCM_INUM (y);
3142 return SCM_BOOL (xx == yy);
3143 }
3144 else if (SCM_BIGP (y))
3145 return SCM_BOOL_F;
3146 else if (SCM_REALP (y))
3147 return SCM_BOOL ((double) xx == SCM_REAL_VALUE (y));
3148 else if (SCM_COMPLEXP (y))
3149 return SCM_BOOL (((double) xx == SCM_COMPLEX_REAL (y))
3150 && (0.0 == SCM_COMPLEX_IMAG (y)));
f92e85f7
MV		 3151 else if (SCM_FRACTIONP (y))
3152 return SCM_BOOL_F;
0aacf84e
MD		 3153 else
3154 SCM_WTA_DISPATCH_2 (g_eq_p, x, y, SCM_ARGn, s_eq_p);
f872b822 3155 }
0aacf84e
MD		 3156 else if (SCM_BIGP (x))
3157 {
3158 if (SCM_INUMP (y))
3159 return SCM_BOOL_F;
3160 else if (SCM_BIGP (y))
3161 {
3162 int cmp = mpz_cmp (SCM_I_BIG_MPZ (x), SCM_I_BIG_MPZ (y));
3163 scm_remember_upto_here_2 (x, y);
3164 return SCM_BOOL (0 == cmp);
3165 }
3166 else if (SCM_REALP (y))
3167 {
3168 int cmp;
3169 if (xisnan (SCM_REAL_VALUE (y)))
3170 return SCM_BOOL_F;
3171 cmp = xmpz_cmp_d (SCM_I_BIG_MPZ (x), SCM_REAL_VALUE (y));
3172 scm_remember_upto_here_1 (x);
3173 return SCM_BOOL (0 == cmp);
3174 }
3175 else if (SCM_COMPLEXP (y))
3176 {
3177 int cmp;
3178 if (0.0 != SCM_COMPLEX_IMAG (y))
3179 return SCM_BOOL_F;
3180 if (xisnan (SCM_COMPLEX_REAL (y)))
3181 return SCM_BOOL_F;
3182 cmp = xmpz_cmp_d (SCM_I_BIG_MPZ (x), SCM_COMPLEX_REAL (y));
3183 scm_remember_upto_here_1 (x);
3184 return SCM_BOOL (0 == cmp);
3185 }
f92e85f7
MV		 3186 else if (SCM_FRACTIONP (y))
3187 return SCM_BOOL_F;
0aacf84e
MD		 3188 else
3189 SCM_WTA_DISPATCH_2 (g_eq_p, x, y, SCM_ARGn, s_eq_p);
f4c627b3 3190 }
0aacf84e
MD		 3191 else if (SCM_REALP (x))
3192 {
3193 if (SCM_INUMP (y))
3194 return SCM_BOOL (SCM_REAL_VALUE (x) == (double) SCM_INUM (y));
3195 else if (SCM_BIGP (y))
3196 {
3197 int cmp;
3198 if (xisnan (SCM_REAL_VALUE (x)))
3199 return SCM_BOOL_F;
3200 cmp = xmpz_cmp_d (SCM_I_BIG_MPZ (y), SCM_REAL_VALUE (x));
3201 scm_remember_upto_here_1 (y);
3202 return SCM_BOOL (0 == cmp);
3203 }
3204 else if (SCM_REALP (y))
3205 return SCM_BOOL (SCM_REAL_VALUE (x) == SCM_REAL_VALUE (y));
3206 else if (SCM_COMPLEXP (y))
3207 return SCM_BOOL ((SCM_REAL_VALUE (x) == SCM_COMPLEX_REAL (y))
3208 && (0.0 == SCM_COMPLEX_IMAG (y)));
f92e85f7 3209 else if (SCM_FRACTIONP (y))
d8b95e27
KR		 3210 {
3211 double xx = SCM_REAL_VALUE (x);
3212 if (xisnan (xx))
3213 return SCM_BOOL_F;
3214 if (xisinf (xx))
3215 return SCM_BOOL (xx < 0.0);
3216 x = scm_inexact_to_exact (x); /* with x as frac or int */
3217 goto again;
3218 }
0aacf84e
MD		 3219 else
3220 SCM_WTA_DISPATCH_2 (g_eq_p, x, y, SCM_ARGn, s_eq_p);
f872b822 3221 }
0aacf84e
MD		 3222 else if (SCM_COMPLEXP (x))
3223 {
3224 if (SCM_INUMP (y))
3225 return SCM_BOOL ((SCM_COMPLEX_REAL (x) == (double) SCM_INUM (y))
3226 && (SCM_COMPLEX_IMAG (x) == 0.0));
3227 else if (SCM_BIGP (y))
3228 {
3229 int cmp;
3230 if (0.0 != SCM_COMPLEX_IMAG (x))
3231 return SCM_BOOL_F;
3232 if (xisnan (SCM_COMPLEX_REAL (x)))
3233 return SCM_BOOL_F;
3234 cmp = xmpz_cmp_d (SCM_I_BIG_MPZ (y), SCM_COMPLEX_REAL (x));
3235 scm_remember_upto_here_1 (y);
3236 return SCM_BOOL (0 == cmp);
3237 }
3238 else if (SCM_REALP (y))
3239 return SCM_BOOL ((SCM_COMPLEX_REAL (x) == SCM_REAL_VALUE (y))
3240 && (SCM_COMPLEX_IMAG (x) == 0.0));
3241 else if (SCM_COMPLEXP (y))
3242 return SCM_BOOL ((SCM_COMPLEX_REAL (x) == SCM_COMPLEX_REAL (y))
3243 && (SCM_COMPLEX_IMAG (x) == SCM_COMPLEX_IMAG (y)));
f92e85f7 3244 else if (SCM_FRACTIONP (y))
d8b95e27
KR		 3245 {
3246 double xx;
3247 if (SCM_COMPLEX_IMAG (x) != 0.0)
3248 return SCM_BOOL_F;
3249 xx = SCM_COMPLEX_REAL (x);
3250 if (xisnan (xx))
3251 return SCM_BOOL_F;
3252 if (xisinf (xx))
3253 return SCM_BOOL (xx < 0.0);
3254 x = scm_inexact_to_exact (x); /* with x as frac or int */
3255 goto again;
3256 }
f92e85f7
MV		 3257 else
3258 SCM_WTA_DISPATCH_2 (g_eq_p, x, y, SCM_ARGn, s_eq_p);
3259 }
3260 else if (SCM_FRACTIONP (x))
3261 {
3262 if (SCM_INUMP (y))
3263 return SCM_BOOL_F;
3264 else if (SCM_BIGP (y))
3265 return SCM_BOOL_F;
3266 else if (SCM_REALP (y))
d8b95e27
KR		 3267 {
3268 double yy = SCM_REAL_VALUE (y);
3269 if (xisnan (yy))
3270 return SCM_BOOL_F;
3271 if (xisinf (yy))
3272 return SCM_BOOL (0.0 < yy);
3273 y = scm_inexact_to_exact (y); /* with y as frac or int */
3274 goto again;
3275 }
f92e85f7 3276 else if (SCM_COMPLEXP (y))
d8b95e27
KR		 3277 {
3278 double yy;
3279 if (SCM_COMPLEX_IMAG (y) != 0.0)
3280 return SCM_BOOL_F;
3281 yy = SCM_COMPLEX_REAL (y);
3282 if (xisnan (yy))
3283 return SCM_BOOL_F;
3284 if (xisinf (yy))
3285 return SCM_BOOL (0.0 < yy);
3286 y = scm_inexact_to_exact (y); /* with y as frac or int */
3287 goto again;
3288 }
f92e85f7
MV		 3289 else if (SCM_FRACTIONP (y))
3290 return scm_i_fraction_equalp (x, y);
0aacf84e
MD		 3291 else
3292 SCM_WTA_DISPATCH_2 (g_eq_p, x, y, SCM_ARGn, s_eq_p);
f4c627b3 3293 }
0aacf84e 3294 else
f4c627b3 3295 SCM_WTA_DISPATCH_2 (g_eq_p, x, y, SCM_ARG1, s_eq_p);
0f2d19dd
JB		 3296}
3297
3298
a5f0b599
KR		 3299/* OPTIMIZE-ME: For int/frac and frac/frac compares, the multiplications
3300 done are good for inums, but for bignums an answer can almost always be
3301 had by just examining a few high bits of the operands, as done by GMP in
3302 mpq_cmp. flonum/frac compares likewise, but with the slight complication
3303 of the float exponent to take into account. */
3304
152f82bf 3305SCM_GPROC1 (s_less_p, "<", scm_tc7_rpsubr, scm_less_p, g_less_p);
942e5b91
MG		 3306/* "Return @code{#t} if the list of parameters is monotonically\n"
3307 * "increasing."
3308 */
0f2d19dd 3309SCM
6e8d25a6 3310scm_less_p (SCM x, SCM y)
0f2d19dd 3311{
a5f0b599 3312 again:
0aacf84e
MD		 3313 if (SCM_INUMP (x))
3314 {
3315 long xx = SCM_INUM (x);
3316 if (SCM_INUMP (y))
3317 {
3318 long yy = SCM_INUM (y);
3319 return SCM_BOOL (xx < yy);
3320 }
3321 else if (SCM_BIGP (y))
3322 {
3323 int sgn = mpz_sgn (SCM_I_BIG_MPZ (y));
3324 scm_remember_upto_here_1 (y);
3325 return SCM_BOOL (sgn > 0);
3326 }
3327 else if (SCM_REALP (y))
3328 return SCM_BOOL ((double) xx < SCM_REAL_VALUE (y));
f92e85f7 3329 else if (SCM_FRACTIONP (y))
a5f0b599
KR		 3330 {
3331 /* "x < a/b" becomes "x*b < a" */
3332 int_frac:
3333 x = scm_product (x, SCM_FRACTION_DENOMINATOR (y));
3334 y = SCM_FRACTION_NUMERATOR (y);
3335 goto again;
3336 }
0aacf84e
MD		 3337 else
3338 SCM_WTA_DISPATCH_2 (g_less_p, x, y, SCM_ARGn, s_less_p);
f872b822 3339 }
0aacf84e
MD		 3340 else if (SCM_BIGP (x))
3341 {
3342 if (SCM_INUMP (y))
3343 {
3344 int sgn = mpz_sgn (SCM_I_BIG_MPZ (x));
3345 scm_remember_upto_here_1 (x);
3346 return SCM_BOOL (sgn < 0);
3347 }
3348 else if (SCM_BIGP (y))
3349 {
3350 int cmp = mpz_cmp (SCM_I_BIG_MPZ (x), SCM_I_BIG_MPZ (y));
3351 scm_remember_upto_here_2 (x, y);
3352 return SCM_BOOL (cmp < 0);
3353 }
3354 else if (SCM_REALP (y))
3355 {
3356 int cmp;
3357 if (xisnan (SCM_REAL_VALUE (y)))
3358 return SCM_BOOL_F;
3359 cmp = xmpz_cmp_d (SCM_I_BIG_MPZ (x), SCM_REAL_VALUE (y));
3360 scm_remember_upto_here_1 (x);
3361 return SCM_BOOL (cmp < 0);
3362 }
f92e85f7 3363 else if (SCM_FRACTIONP (y))
a5f0b599 3364 goto int_frac;
0aacf84e
MD		 3365 else
3366 SCM_WTA_DISPATCH_2 (g_less_p, x, y, SCM_ARGn, s_less_p);
f4c627b3 3367 }
0aacf84e
MD		 3368 else if (SCM_REALP (x))
3369 {
3370 if (SCM_INUMP (y))
3371 return SCM_BOOL (SCM_REAL_VALUE (x) < (double) SCM_INUM (y));
3372 else if (SCM_BIGP (y))
3373 {
3374 int cmp;
3375 if (xisnan (SCM_REAL_VALUE (x)))
3376 return SCM_BOOL_F;
3377 cmp = xmpz_cmp_d (SCM_I_BIG_MPZ (y), SCM_REAL_VALUE (x));
3378 scm_remember_upto_here_1 (y);
3379 return SCM_BOOL (cmp > 0);
3380 }
3381 else if (SCM_REALP (y))
3382 return SCM_BOOL (SCM_REAL_VALUE (x) < SCM_REAL_VALUE (y));
f92e85f7 3383 else if (SCM_FRACTIONP (y))
a5f0b599
KR		 3384 {
3385 double xx = SCM_REAL_VALUE (x);
3386 if (xisnan (xx))
3387 return SCM_BOOL_F;
3388 if (xisinf (xx))
3389 return SCM_BOOL (xx < 0.0);
3390 x = scm_inexact_to_exact (x); /* with x as frac or int */
3391 goto again;
3392 }
f92e85f7
MV		 3393 else
3394 SCM_WTA_DISPATCH_2 (g_less_p, x, y, SCM_ARGn, s_less_p);
3395 }
3396 else if (SCM_FRACTIONP (x))
3397 {
a5f0b599
KR		 3398 if (SCM_INUMP (y) || SCM_BIGP (y))
3399 {
3400 /* "a/b < y" becomes "a < y*b" */
3401 y = scm_product (y, SCM_FRACTION_DENOMINATOR (x));
3402 x = SCM_FRACTION_NUMERATOR (x);
3403 goto again;
3404 }
f92e85f7 3405 else if (SCM_REALP (y))
a5f0b599
KR		 3406 {
3407 double yy = SCM_REAL_VALUE (y);
3408 if (xisnan (yy))
3409 return SCM_BOOL_F;
3410 if (xisinf (yy))
3411 return SCM_BOOL (0.0 < yy);
3412 y = scm_inexact_to_exact (y); /* with y as frac or int */
3413 goto again;
3414 }
f92e85f7 3415 else if (SCM_FRACTIONP (y))
a5f0b599
KR		 3416 {
3417 /* "a/b < c/d" becomes "a*d < c*b" */
3418 SCM new_x = scm_product (SCM_FRACTION_NUMERATOR (x),
3419 SCM_FRACTION_DENOMINATOR (y));
3420 SCM new_y = scm_product (SCM_FRACTION_NUMERATOR (y),
3421 SCM_FRACTION_DENOMINATOR (x));
3422 x = new_x;
3423 y = new_y;
3424 goto again;
3425 }
0aacf84e
MD		 3426 else
3427 SCM_WTA_DISPATCH_2 (g_less_p, x, y, SCM_ARGn, s_less_p);
f872b822 3428 }
0aacf84e 3429 else
f4c627b3 3430 SCM_WTA_DISPATCH_2 (g_less_p, x, y, SCM_ARG1, s_less_p);
0f2d19dd
JB		 3431}
3432
3433
c76b1eaf 3434SCM_GPROC1 (s_scm_gr_p, ">", scm_tc7_rpsubr, scm_gr_p, g_gr_p);
942e5b91
MG		 3435/* "Return @code{#t} if the list of parameters is monotonically\n"
3436 * "decreasing."
c76b1eaf 3437 */
1bbd0b84 3438#define FUNC_NAME s_scm_gr_p
c76b1eaf
MD		 3439SCM
3440scm_gr_p (SCM x, SCM y)
0f2d19dd 3441{
c76b1eaf
MD		 3442 if (!SCM_NUMBERP (x))
3443 SCM_WTA_DISPATCH_2 (g_gr_p, x, y, SCM_ARG1, FUNC_NAME);
3444 else if (!SCM_NUMBERP (y))
3445 SCM_WTA_DISPATCH_2 (g_gr_p, x, y, SCM_ARG2, FUNC_NAME);
3446 else
3447 return scm_less_p (y, x);
0f2d19dd 3448}
1bbd0b84 3449#undef FUNC_NAME
0f2d19dd
JB		 3450
3451
c76b1eaf 3452SCM_GPROC1 (s_scm_leq_p, "<=", scm_tc7_rpsubr, scm_leq_p, g_leq_p);
942e5b91 3453/* "Return @code{#t} if the list of parameters is monotonically\n"
c76b1eaf
MD		 3454 * "non-decreasing."
3455 */
1bbd0b84 3456#define FUNC_NAME s_scm_leq_p
c76b1eaf
MD		 3457SCM
3458scm_leq_p (SCM x, SCM y)
0f2d19dd 3459{
c76b1eaf
MD		 3460 if (!SCM_NUMBERP (x))
3461 SCM_WTA_DISPATCH_2 (g_leq_p, x, y, SCM_ARG1, FUNC_NAME);
3462 else if (!SCM_NUMBERP (y))
3463 SCM_WTA_DISPATCH_2 (g_leq_p, x, y, SCM_ARG2, FUNC_NAME);
fc194577
MV		 3464 else if (SCM_NFALSEP (scm_nan_p (x)) || SCM_NFALSEP (scm_nan_p (y)))
3465 return SCM_BOOL_F;
c76b1eaf
MD		 3466 else
3467 return SCM_BOOL_NOT (scm_less_p (y, x));
0f2d19dd 3468}
1bbd0b84 3469#undef FUNC_NAME
0f2d19dd
JB		 3470
3471
c76b1eaf 3472SCM_GPROC1 (s_scm_geq_p, ">=", scm_tc7_rpsubr, scm_geq_p, g_geq_p);
942e5b91 3473/* "Return @code{#t} if the list of parameters is monotonically\n"
c76b1eaf
MD		 3474 * "non-increasing."
3475 */
1bbd0b84 3476#define FUNC_NAME s_scm_geq_p
c76b1eaf
MD		 3477SCM
3478scm_geq_p (SCM x, SCM y)
0f2d19dd 3479{
c76b1eaf
MD		 3480 if (!SCM_NUMBERP (x))
3481 SCM_WTA_DISPATCH_2 (g_geq_p, x, y, SCM_ARG1, FUNC_NAME);
3482 else if (!SCM_NUMBERP (y))
3483 SCM_WTA_DISPATCH_2 (g_geq_p, x, y, SCM_ARG2, FUNC_NAME);
fc194577
MV		 3484 else if (SCM_NFALSEP (scm_nan_p (x)) || SCM_NFALSEP (scm_nan_p (y)))
3485 return SCM_BOOL_F;
c76b1eaf 3486 else
fc194577 3487 return SCM_BOOL_NOT (scm_less_p (x, y));
0f2d19dd 3488}
1bbd0b84 3489#undef FUNC_NAME
0f2d19dd
JB		 3490
3491
152f82bf 3492SCM_GPROC (s_zero_p, "zero?", 1, 0, 0, scm_zero_p, g_zero_p);
942e5b91
MG		 3493/* "Return @code{#t} if @var{z} is an exact or inexact number equal to\n"
3494 * "zero."
3495 */
0f2d19dd 3496SCM
6e8d25a6 3497scm_zero_p (SCM z)
0f2d19dd 3498{
0aacf84e 3499 if (SCM_INUMP (z))
c2ff8ab0 3500 return SCM_BOOL (SCM_EQ_P (z, SCM_INUM0));
0aacf84e 3501 else if (SCM_BIGP (z))
c2ff8ab0 3502 return SCM_BOOL_F;
0aacf84e 3503 else if (SCM_REALP (z))
c2ff8ab0 3504 return SCM_BOOL (SCM_REAL_VALUE (z) == 0.0);
0aacf84e 3505 else if (SCM_COMPLEXP (z))
c2ff8ab0
DH		 3506 return SCM_BOOL (SCM_COMPLEX_REAL (z) == 0.0
3507 && SCM_COMPLEX_IMAG (z) == 0.0);
f92e85f7
MV		 3508 else if (SCM_FRACTIONP (z))
3509 return SCM_BOOL_F;
0aacf84e 3510 else
c2ff8ab0 3511 SCM_WTA_DISPATCH_1 (g_zero_p, z, SCM_ARG1, s_zero_p);
0f2d19dd
JB		 3512}
3513
3514
152f82bf 3515SCM_GPROC (s_positive_p, "positive?", 1, 0, 0, scm_positive_p, g_positive_p);
942e5b91
MG		 3516/* "Return @code{#t} if @var{x} is an exact or inexact number greater than\n"
3517 * "zero."
3518 */
0f2d19dd 3519SCM
6e8d25a6 3520scm_positive_p (SCM x)
0f2d19dd 3521{
0aacf84e 3522 if (SCM_INUMP (x))
c2ff8ab0 3523 return SCM_BOOL (SCM_INUM (x) > 0);
0aacf84e
MD		 3524 else if (SCM_BIGP (x))
3525 {
3526 int sgn = mpz_sgn (SCM_I_BIG_MPZ (x));
3527 scm_remember_upto_here_1 (x);
3528 return SCM_BOOL (sgn > 0);
3529 }
3530 else if (SCM_REALP (x))
c2ff8ab0 3531 return SCM_BOOL(SCM_REAL_VALUE (x) > 0.0);
f92e85f7
MV		 3532 else if (SCM_FRACTIONP (x))
3533 return scm_positive_p (SCM_FRACTION_NUMERATOR (x));
0aacf84e 3534 else
c2ff8ab0 3535 SCM_WTA_DISPATCH_1 (g_positive_p, x, SCM_ARG1, s_positive_p);
0f2d19dd
JB		 3536}
3537
3538
152f82bf 3539SCM_GPROC (s_negative_p, "negative?", 1, 0, 0, scm_negative_p, g_negative_p);
942e5b91
MG		 3540/* "Return @code{#t} if @var{x} is an exact or inexact number less than\n"
3541 * "zero."
3542 */
0f2d19dd 3543SCM
6e8d25a6 3544scm_negative_p (SCM x)
0f2d19dd 3545{
0aacf84e 3546 if (SCM_INUMP (x))
c2ff8ab0 3547 return SCM_BOOL (SCM_INUM (x) < 0);
0aacf84e
MD		 3548 else if (SCM_BIGP (x))
3549 {
3550 int sgn = mpz_sgn (SCM_I_BIG_MPZ (x));
3551 scm_remember_upto_here_1 (x);
3552 return SCM_BOOL (sgn < 0);
3553 }
3554 else if (SCM_REALP (x))
c2ff8ab0 3555 return SCM_BOOL(SCM_REAL_VALUE (x) < 0.0);
f92e85f7
MV		 3556 else if (SCM_FRACTIONP (x))
3557 return scm_negative_p (SCM_FRACTION_NUMERATOR (x));
0aacf84e 3558 else
c2ff8ab0 3559 SCM_WTA_DISPATCH_1 (g_negative_p, x, SCM_ARG1, s_negative_p);
0f2d19dd
JB		 3560}
3561
3562
2a06f791
KR		 3563/* scm_min and scm_max return an inexact when either argument is inexact, as
3564 required by r5rs. On that basis, for exact/inexact combinations the
3565 exact is converted to inexact to compare and possibly return. This is
3566 unlike scm_less_p above which takes some trouble to preserve all bits in
3567 its test, such trouble is not required for min and max. */
3568
9de33deb 3569SCM_GPROC1 (s_max, "max", scm_tc7_asubr, scm_max, g_max);
942e5b91
MG		 3570/* "Return the maximum of all parameter values."
3571 */
0f2d19dd 3572SCM
6e8d25a6 3573scm_max (SCM x, SCM y)
0f2d19dd 3574{
0aacf84e
MD		 3575 if (SCM_UNBNDP (y))
3576 {
3577 if (SCM_UNBNDP (x))
3578 SCM_WTA_DISPATCH_0 (g_max, s_max);
dab4e67a 3579 else if (SCM_INUMP(x) || SCM_BIGP(x) || SCM_REALP(x) || SCM_FRACTIONP(x))
0aacf84e
MD		 3580 return x;
3581 else
3582 SCM_WTA_DISPATCH_1 (g_max, x, SCM_ARG1, s_max);
f872b822 3583 }
f4c627b3 3584
0aacf84e
MD		 3585 if (SCM_INUMP (x))
3586 {
3587 long xx = SCM_INUM (x);
3588 if (SCM_INUMP (y))
3589 {
3590 long yy = SCM_INUM (y);
3591 return (xx < yy) ? y : x;
3592 }
3593 else if (SCM_BIGP (y))
3594 {
3595 int sgn = mpz_sgn (SCM_I_BIG_MPZ (y));
3596 scm_remember_upto_here_1 (y);
3597 return (sgn < 0) ? x : y;
3598 }
3599 else if (SCM_REALP (y))
3600 {
3601 double z = xx;
3602 /* if y==NaN then ">" is false and we return NaN */
3603 return (z > SCM_REAL_VALUE (y)) ? scm_make_real (z) : y;
3604 }
f92e85f7
MV		 3605 else if (SCM_FRACTIONP (y))
3606 {
e4bc5d6c
KR		 3607 use_less:
3608 return (SCM_FALSEP (scm_less_p (x, y)) ? x : y);
f92e85f7 3609 }
0aacf84e
MD		 3610 else
3611 SCM_WTA_DISPATCH_2 (g_max, x, y, SCM_ARGn, s_max);
f872b822 3612 }
0aacf84e
MD		 3613 else if (SCM_BIGP (x))
3614 {
3615 if (SCM_INUMP (y))
3616 {
3617 int sgn = mpz_sgn (SCM_I_BIG_MPZ (x));
3618 scm_remember_upto_here_1 (x);
3619 return (sgn < 0) ? y : x;
3620 }
3621 else if (SCM_BIGP (y))
3622 {
3623 int cmp = mpz_cmp (SCM_I_BIG_MPZ (x), SCM_I_BIG_MPZ (y));
3624 scm_remember_upto_here_2 (x, y);
3625 return (cmp > 0) ? x : y;
3626 }
3627 else if (SCM_REALP (y))
3628 {
2a06f791
KR		 3629 /* if y==NaN then xx>yy is false, so we return the NaN y */
3630 double xx, yy;
3631 big_real:
3632 xx = scm_i_big2dbl (x);
3633 yy = SCM_REAL_VALUE (y);
3634 return (xx > yy ? scm_make_real (xx) : y);
0aacf84e 3635 }
f92e85f7
MV		 3636 else if (SCM_FRACTIONP (y))
3637 {
e4bc5d6c 3638 goto use_less;
f92e85f7 3639 }
0aacf84e
MD		 3640 else
3641 SCM_WTA_DISPATCH_2 (g_max, x, y, SCM_ARGn, s_max);
f4c627b3 3642 }
0aacf84e
MD		 3643 else if (SCM_REALP (x))
3644 {
3645 if (SCM_INUMP (y))
3646 {
3647 double z = SCM_INUM (y);
3648 /* if x==NaN then "<" is false and we return NaN */
3649 return (SCM_REAL_VALUE (x) < z) ? scm_make_real (z) : x;
3650 }
3651 else if (SCM_BIGP (y))
3652 {
b6f8f763 3653 SCM_SWAP (x, y);
2a06f791 3654 goto big_real;
0aacf84e
MD		 3655 }
3656 else if (SCM_REALP (y))
3657 {
3658 /* if x==NaN then our explicit check means we return NaN
3659 if y==NaN then ">" is false and we return NaN
3660 calling isnan is unavoidable, since it's the only way to know
3661 which of x or y causes any compares to be false */
3662 double xx = SCM_REAL_VALUE (x);
3663 return (xisnan (xx) || xx > SCM_REAL_VALUE (y)) ? x : y;
3664 }
f92e85f7
MV		 3665 else if (SCM_FRACTIONP (y))
3666 {
3667 double yy = scm_i_fraction2double (y);
3668 double xx = SCM_REAL_VALUE (x);
3669 return (xx < yy) ? scm_make_real (yy) : x;
3670 }
3671 else
3672 SCM_WTA_DISPATCH_2 (g_max, x, y, SCM_ARGn, s_max);
3673 }
3674 else if (SCM_FRACTIONP (x))
3675 {
3676 if (SCM_INUMP (y))
3677 {
e4bc5d6c 3678 goto use_less;
f92e85f7
MV		 3679 }
3680 else if (SCM_BIGP (y))
3681 {
e4bc5d6c 3682 goto use_less;
f92e85f7
MV		 3683 }
3684 else if (SCM_REALP (y))
3685 {
3686 double xx = scm_i_fraction2double (x);
3687 return (xx < SCM_REAL_VALUE (y)) ? y : scm_make_real (xx);
3688 }
3689 else if (SCM_FRACTIONP (y))
3690 {
e4bc5d6c 3691 goto use_less;
f92e85f7 3692 }
0aacf84e
MD		 3693 else
3694 SCM_WTA_DISPATCH_2 (g_max, x, y, SCM_ARGn, s_max);
f872b822 3695 }
0aacf84e 3696 else
f4c627b3 3697 SCM_WTA_DISPATCH_2 (g_max, x, y, SCM_ARG1, s_max);
0f2d19dd
JB		 3698}
3699
3700
9de33deb 3701SCM_GPROC1 (s_min, "min", scm_tc7_asubr, scm_min, g_min);
942e5b91
MG		 3702/* "Return the minium of all parameter values."
3703 */
0f2d19dd 3704SCM
6e8d25a6 3705scm_min (SCM x, SCM y)
0f2d19dd 3706{
0aacf84e
MD		 3707 if (SCM_UNBNDP (y))
3708 {
3709 if (SCM_UNBNDP (x))
3710 SCM_WTA_DISPATCH_0 (g_min, s_min);
dab4e67a 3711 else if (SCM_INUMP(x) || SCM_BIGP(x) || SCM_REALP(x) || SCM_FRACTIONP(x))
0aacf84e
MD		 3712 return x;
3713 else
3714 SCM_WTA_DISPATCH_1 (g_min, x, SCM_ARG1, s_min);
f872b822 3715 }
f4c627b3 3716
0aacf84e
MD		 3717 if (SCM_INUMP (x))
3718 {
3719 long xx = SCM_INUM (x);
3720 if (SCM_INUMP (y))
3721 {
3722 long yy = SCM_INUM (y);
3723 return (xx < yy) ? x : y;
3724 }
3725 else if (SCM_BIGP (y))
3726 {
3727 int sgn = mpz_sgn (SCM_I_BIG_MPZ (y));
3728 scm_remember_upto_here_1 (y);
3729 return (sgn < 0) ? y : x;
3730 }
3731 else if (SCM_REALP (y))
3732 {
3733 double z = xx;
3734 /* if y==NaN then "<" is false and we return NaN */
3735 return (z < SCM_REAL_VALUE (y)) ? scm_make_real (z) : y;
3736 }
f92e85f7
MV		 3737 else if (SCM_FRACTIONP (y))
3738 {
e4bc5d6c
KR		 3739 use_less:
3740 return (SCM_FALSEP (scm_less_p (x, y)) ? y : x);
f92e85f7 3741 }
0aacf84e
MD		 3742 else
3743 SCM_WTA_DISPATCH_2 (g_min, x, y, SCM_ARGn, s_min);
f872b822 3744 }
0aacf84e
MD		 3745 else if (SCM_BIGP (x))
3746 {
3747 if (SCM_INUMP (y))
3748 {
3749 int sgn = mpz_sgn (SCM_I_BIG_MPZ (x));
3750 scm_remember_upto_here_1 (x);
3751 return (sgn < 0) ? x : y;
3752 }
3753 else if (SCM_BIGP (y))
3754 {
3755 int cmp = mpz_cmp (SCM_I_BIG_MPZ (x), SCM_I_BIG_MPZ (y));
3756 scm_remember_upto_here_2 (x, y);
3757 return (cmp > 0) ? y : x;
3758 }
3759 else if (SCM_REALP (y))
3760 {
2a06f791
KR		 3761 /* if y==NaN then xx<yy is false, so we return the NaN y */
3762 double xx, yy;
3763 big_real:
3764 xx = scm_i_big2dbl (x);
3765 yy = SCM_REAL_VALUE (y);
3766 return (xx < yy ? scm_make_real (xx) : y);
0aacf84e 3767 }
f92e85f7
MV		 3768 else if (SCM_FRACTIONP (y))
3769 {
e4bc5d6c 3770 goto use_less;
f92e85f7 3771 }
0aacf84e
MD		 3772 else
3773 SCM_WTA_DISPATCH_2 (g_min, x, y, SCM_ARGn, s_min);
f4c627b3 3774 }
0aacf84e
MD		 3775 else if (SCM_REALP (x))
3776 {
3777 if (SCM_INUMP (y))
3778 {
3779 double z = SCM_INUM (y);
3780 /* if x==NaN then "<" is false and we return NaN */
3781 return (z < SCM_REAL_VALUE (x)) ? scm_make_real (z) : x;
3782 }
3783 else if (SCM_BIGP (y))
3784 {
b6f8f763 3785 SCM_SWAP (x, y);
2a06f791 3786 goto big_real;
0aacf84e
MD		 3787 }
3788 else if (SCM_REALP (y))
3789 {
3790 /* if x==NaN then our explicit check means we return NaN
3791 if y==NaN then "<" is false and we return NaN
3792 calling isnan is unavoidable, since it's the only way to know
3793 which of x or y causes any compares to be false */
3794 double xx = SCM_REAL_VALUE (x);
3795 return (xisnan (xx) || xx < SCM_REAL_VALUE (y)) ? x : y;
3796 }
f92e85f7
MV		 3797 else if (SCM_FRACTIONP (y))
3798 {
3799 double yy = scm_i_fraction2double (y);
3800 double xx = SCM_REAL_VALUE (x);
3801 return (yy < xx) ? scm_make_real (yy) : x;
3802 }
0aacf84e
MD		 3803 else
3804 SCM_WTA_DISPATCH_2 (g_min, x, y, SCM_ARGn, s_min);
f872b822 3805 }
f92e85f7
MV		 3806 else if (SCM_FRACTIONP (x))
3807 {
3808 if (SCM_INUMP (y))
3809 {
e4bc5d6c 3810 goto use_less;
f92e85f7
MV		 3811 }
3812 else if (SCM_BIGP (y))
3813 {
e4bc5d6c 3814 goto use_less;
f92e85f7
MV		 3815 }
3816 else if (SCM_REALP (y))
3817 {
3818 double xx = scm_i_fraction2double (x);
3819 return (SCM_REAL_VALUE (y) < xx) ? y : scm_make_real (xx);
3820 }
3821 else if (SCM_FRACTIONP (y))
3822 {
e4bc5d6c 3823 goto use_less;
f92e85f7
MV		 3824 }
3825 else
3826 SCM_WTA_DISPATCH_2 (g_max, x, y, SCM_ARGn, s_max);
3827 }
0aacf84e 3828 else
f4c627b3 3829 SCM_WTA_DISPATCH_2 (g_min, x, y, SCM_ARG1, s_min);
0f2d19dd
JB		 3830}
3831
3832
9de33deb 3833SCM_GPROC1 (s_sum, "+", scm_tc7_asubr, scm_sum, g_sum);
942e5b91
MG		 3834/* "Return the sum of all parameter values. Return 0 if called without\n"
3835 * "any parameters."
3836 */
0f2d19dd 3837SCM
6e8d25a6 3838scm_sum (SCM x, SCM y)
0f2d19dd 3839{
ca46fb90
RB		 3840 if (SCM_UNBNDP (y))
3841 {
3842 if (SCM_NUMBERP (x)) return x;
3843 if (SCM_UNBNDP (x)) return SCM_INUM0;
98cb6e75 3844 SCM_WTA_DISPATCH_1 (g_sum, x, SCM_ARG1, s_sum);
f872b822 3845 }
c209c88e 3846
ca46fb90
RB		 3847 if (SCM_INUMP (x))
3848 {
3849 if (SCM_INUMP (y))
3850 {
3851 long xx = SCM_INUM (x);
3852 long yy = SCM_INUM (y);
3853 long int z = xx + yy;
3854 return SCM_FIXABLE (z) ? SCM_MAKINUM (z) : scm_i_long2big (z);
3855 }
3856 else if (SCM_BIGP (y))
3857 {
3858 SCM_SWAP (x, y);
3859 goto add_big_inum;
3860 }
3861 else if (SCM_REALP (y))
3862 {
3863 long int xx = SCM_INUM (x);
3864 return scm_make_real (xx + SCM_REAL_VALUE (y));
3865 }
3866 else if (SCM_COMPLEXP (y))
3867 {
3868 long int xx = SCM_INUM (x);
3869 return scm_make_complex (xx + SCM_COMPLEX_REAL (y),
3870 SCM_COMPLEX_IMAG (y));
3871 }
f92e85f7
MV		 3872 else if (SCM_FRACTIONP (y))
3873 return scm_make_ratio (scm_sum (SCM_FRACTION_NUMERATOR (y),
3874 scm_product (x, SCM_FRACTION_DENOMINATOR (y))),
3875 SCM_FRACTION_DENOMINATOR (y));
ca46fb90
RB		 3876 else
3877 SCM_WTA_DISPATCH_2 (g_sum, x, y, SCM_ARGn, s_sum);
0aacf84e
MD		 3878 } else if (SCM_BIGP (x))
3879 {
3880 if (SCM_INUMP (y))
3881 {
3882 long int inum;
3883 int bigsgn;
3884 add_big_inum:
3885 inum = SCM_INUM (y);
3886 if (inum == 0)
3887 return x;
3888 bigsgn = mpz_sgn (SCM_I_BIG_MPZ (x));
3889 if (inum < 0)
3890 {
3891 SCM result = scm_i_mkbig ();
3892 mpz_sub_ui (SCM_I_BIG_MPZ (result), SCM_I_BIG_MPZ (x), - inum);
3893 scm_remember_upto_here_1 (x);
3894 /* we know the result will have to be a bignum */
3895 if (bigsgn == -1)
3896 return result;
3897 return scm_i_normbig (result);
3898 }
3899 else
3900 {
3901 SCM result = scm_i_mkbig ();
3902 mpz_add_ui (SCM_I_BIG_MPZ (result), SCM_I_BIG_MPZ (x), inum);
3903 scm_remember_upto_here_1 (x);
3904 /* we know the result will have to be a bignum */
3905 if (bigsgn == 1)
3906 return result;
3907 return scm_i_normbig (result);
3908 }
3909 }
3910 else if (SCM_BIGP (y))
3911 {
3912 SCM result = scm_i_mkbig ();
3913 int sgn_x = mpz_sgn (SCM_I_BIG_MPZ (x));
3914 int sgn_y = mpz_sgn (SCM_I_BIG_MPZ (y));
3915 mpz_add (SCM_I_BIG_MPZ (result),
3916 SCM_I_BIG_MPZ (x),
3917 SCM_I_BIG_MPZ (y));
3918 scm_remember_upto_here_2 (x, y);
3919 /* we know the result will have to be a bignum */
3920 if (sgn_x == sgn_y)
3921 return result;
3922 return scm_i_normbig (result);
3923 }
3924 else if (SCM_REALP (y))
3925 {
3926 double result = mpz_get_d (SCM_I_BIG_MPZ (x)) + SCM_REAL_VALUE (y);
3927 scm_remember_upto_here_1 (x);
3928 return scm_make_real (result);
3929 }
3930 else if (SCM_COMPLEXP (y))
3931 {
3932 double real_part = (mpz_get_d (SCM_I_BIG_MPZ (x))
3933 + SCM_COMPLEX_REAL (y));
3934 scm_remember_upto_here_1 (x);
3935 return scm_make_complex (real_part, SCM_COMPLEX_IMAG (y));
3936 }
f92e85f7
MV		 3937 else if (SCM_FRACTIONP (y))
3938 return scm_make_ratio (scm_sum (SCM_FRACTION_NUMERATOR (y),
3939 scm_product (x, SCM_FRACTION_DENOMINATOR (y))),
3940 SCM_FRACTION_DENOMINATOR (y));
0aacf84e
MD		 3941 else
3942 SCM_WTA_DISPATCH_2 (g_sum, x, y, SCM_ARGn, s_sum);
0f2d19dd 3943 }
0aacf84e
MD		 3944 else if (SCM_REALP (x))
3945 {
3946 if (SCM_INUMP (y))
3947 return scm_make_real (SCM_REAL_VALUE (x) + SCM_INUM (y));
3948 else if (SCM_BIGP (y))
3949 {
3950 double result = mpz_get_d (SCM_I_BIG_MPZ (y)) + SCM_REAL_VALUE (x);
3951 scm_remember_upto_here_1 (y);
3952 return scm_make_real (result);
3953 }
3954 else if (SCM_REALP (y))
3955 return scm_make_real (SCM_REAL_VALUE (x) + SCM_REAL_VALUE (y));
3956 else if (SCM_COMPLEXP (y))
3957 return scm_make_complex (SCM_REAL_VALUE (x) + SCM_COMPLEX_REAL (y),
3958 SCM_COMPLEX_IMAG (y));
f92e85f7
MV		 3959 else if (SCM_FRACTIONP (y))
3960 return scm_make_real (SCM_REAL_VALUE (x) + scm_i_fraction2double (y));
0aacf84e
MD		 3961 else
3962 SCM_WTA_DISPATCH_2 (g_sum, x, y, SCM_ARGn, s_sum);
f872b822 3963 }
0aacf84e
MD		 3964 else if (SCM_COMPLEXP (x))
3965 {
3966 if (SCM_INUMP (y))
3967 return scm_make_complex (SCM_COMPLEX_REAL (x) + SCM_INUM (y),
3968 SCM_COMPLEX_IMAG (x));
3969 else if (SCM_BIGP (y))
3970 {
3971 double real_part = (mpz_get_d (SCM_I_BIG_MPZ (y))
3972 + SCM_COMPLEX_REAL (x));
3973 scm_remember_upto_here_1 (y);
3974 return scm_make_complex (real_part, SCM_COMPLEX_IMAG (x));
3975 }
3976 else if (SCM_REALP (y))
3977 return scm_make_complex (SCM_COMPLEX_REAL (x) + SCM_REAL_VALUE (y),
3978 SCM_COMPLEX_IMAG (x));
3979 else if (SCM_COMPLEXP (y))
3980 return scm_make_complex (SCM_COMPLEX_REAL (x) + SCM_COMPLEX_REAL (y),
3981 SCM_COMPLEX_IMAG (x) + SCM_COMPLEX_IMAG (y));
f92e85f7
MV		 3982 else if (SCM_FRACTIONP (y))
3983 return scm_make_complex (SCM_COMPLEX_REAL (x) + scm_i_fraction2double (y),
3984 SCM_COMPLEX_IMAG (x));
3985 else
3986 SCM_WTA_DISPATCH_2 (g_sum, x, y, SCM_ARGn, s_sum);
3987 }
3988 else if (SCM_FRACTIONP (x))
3989 {
3990 if (SCM_INUMP (y))
3991 return scm_make_ratio (scm_sum (SCM_FRACTION_NUMERATOR (x),
3992 scm_product (y, SCM_FRACTION_DENOMINATOR (x))),
3993 SCM_FRACTION_DENOMINATOR (x));
3994 else if (SCM_BIGP (y))
3995 return scm_make_ratio (scm_sum (SCM_FRACTION_NUMERATOR (x),
3996 scm_product (y, SCM_FRACTION_DENOMINATOR (x))),
3997 SCM_FRACTION_DENOMINATOR (x));
3998 else if (SCM_REALP (y))
3999 return scm_make_real (SCM_REAL_VALUE (y) + scm_i_fraction2double (x));
4000 else if (SCM_COMPLEXP (y))
4001 return scm_make_complex (SCM_COMPLEX_REAL (y) + scm_i_fraction2double (x),
4002 SCM_COMPLEX_IMAG (y));
4003 else if (SCM_FRACTIONP (y))
4004 /* a/b + c/d = (ad + bc) / bd */
4005 return scm_make_ratio (scm_sum (scm_product (SCM_FRACTION_NUMERATOR (x), SCM_FRACTION_DENOMINATOR (y)),
4006 scm_product (SCM_FRACTION_NUMERATOR (y), SCM_FRACTION_DENOMINATOR (x))),
4007 scm_product (SCM_FRACTION_DENOMINATOR (x), SCM_FRACTION_DENOMINATOR (y)));
0aacf84e
MD		 4008 else
4009 SCM_WTA_DISPATCH_2 (g_sum, x, y, SCM_ARGn, s_sum);
98cb6e75 4010 }
0aacf84e 4011 else
98cb6e75 4012 SCM_WTA_DISPATCH_2 (g_sum, x, y, SCM_ARG1, s_sum);
0f2d19dd
JB		 4013}
4014
4015
9de33deb 4016SCM_GPROC1 (s_difference, "-", scm_tc7_asubr, scm_difference, g_difference);
609c3d30
MG		 4017/* If called with one argument @var{z1}, -@var{z1} returned. Otherwise
4018 * the sum of all but the first argument are subtracted from the first
4019 * argument. */
c05e97b7 4020#define FUNC_NAME s_difference
0f2d19dd 4021SCM
6e8d25a6 4022scm_difference (SCM x, SCM y)
0f2d19dd 4023{
ca46fb90
RB		 4024 if (SCM_UNBNDP (y))
4025 {
4026 if (SCM_UNBNDP (x))
4027 SCM_WTA_DISPATCH_0 (g_difference, s_difference);
4028 else
4029 if (SCM_INUMP (x))
4030 {
4031 long xx = -SCM_INUM (x);
4032 if (SCM_FIXABLE (xx))
4033 return SCM_MAKINUM (xx);
4034 else
4035 return scm_i_long2big (xx);
4036 }
4037 else if (SCM_BIGP (x))
4038 /* FIXME: do we really need to normalize here? */
4039 return scm_i_normbig (scm_i_clonebig (x, 0));
4040 else if (SCM_REALP (x))
4041 return scm_make_real (-SCM_REAL_VALUE (x));
4042 else if (SCM_COMPLEXP (x))
4043 return scm_make_complex (-SCM_COMPLEX_REAL (x),
4044 -SCM_COMPLEX_IMAG (x));
f92e85f7
MV		 4045 else if (SCM_FRACTIONP (x))
4046 return scm_make_ratio (scm_difference (SCM_FRACTION_NUMERATOR (x), SCM_UNDEFINED),
4047 SCM_FRACTION_DENOMINATOR (x));
ca46fb90
RB		 4048 else
4049 SCM_WTA_DISPATCH_1 (g_difference, x, SCM_ARG1, s_difference);
f872b822 4050 }
ca46fb90 4051
0aacf84e
MD		 4052 if (SCM_INUMP (x))
4053 {
4054 if (SCM_INUMP (y))
4055 {
4056 long int xx = SCM_INUM (x);
4057 long int yy = SCM_INUM (y);
4058 long int z = xx - yy;
4059 if (SCM_FIXABLE (z))
4060 return SCM_MAKINUM (z);
4061 else
4062 return scm_i_long2big (z);
4063 }
4064 else if (SCM_BIGP (y))
4065 {
4066 /* inum-x - big-y */
4067 long xx = SCM_INUM (x);
ca46fb90 4068
0aacf84e
MD		 4069 if (xx == 0)
4070 return scm_i_clonebig (y, 0);
4071 else
4072 {
4073 int sgn_y = mpz_sgn (SCM_I_BIG_MPZ (y));
4074 SCM result = scm_i_mkbig ();
ca46fb90 4075
0aacf84e
MD		 4076 if (xx >= 0)
4077 mpz_ui_sub (SCM_I_BIG_MPZ (result), xx, SCM_I_BIG_MPZ (y));
4078 else
4079 {
4080 /* x - y == -(y + -x) */
4081 mpz_add_ui (SCM_I_BIG_MPZ (result), SCM_I_BIG_MPZ (y), -xx);
4082 mpz_neg (SCM_I_BIG_MPZ (result), SCM_I_BIG_MPZ (result));
4083 }
4084 scm_remember_upto_here_1 (y);
ca46fb90 4085
0aacf84e
MD		 4086 if ((xx < 0 && (sgn_y > 0)) || ((xx > 0) && sgn_y < 0))
4087 /* we know the result will have to be a bignum */
4088 return result;
4089 else
4090 return scm_i_normbig (result);
4091 }
4092 }
4093 else if (SCM_REALP (y))
4094 {
4095 long int xx = SCM_INUM (x);
4096 return scm_make_real (xx - SCM_REAL_VALUE (y));
4097 }
4098 else if (SCM_COMPLEXP (y))
4099 {
4100 long int xx = SCM_INUM (x);
4101 return scm_make_complex (xx - SCM_COMPLEX_REAL (y),
4102 - SCM_COMPLEX_IMAG (y));
4103 }
f92e85f7
MV		 4104 else if (SCM_FRACTIONP (y))
4105 /* a - b/c = (ac - b) / c */
4106 return scm_make_ratio (scm_difference (scm_product (x, SCM_FRACTION_DENOMINATOR (y)),
4107 SCM_FRACTION_NUMERATOR (y)),
4108 SCM_FRACTION_DENOMINATOR (y));
0aacf84e
MD		 4109 else
4110 SCM_WTA_DISPATCH_2 (g_difference, x, y, SCM_ARGn, s_difference);
f872b822 4111 }
0aacf84e
MD		 4112 else if (SCM_BIGP (x))
4113 {
4114 if (SCM_INUMP (y))
4115 {
4116 /* big-x - inum-y */
4117 long yy = SCM_INUM (y);
4118 int sgn_x = mpz_sgn (SCM_I_BIG_MPZ (x));
ca46fb90 4119
0aacf84e
MD		 4120 scm_remember_upto_here_1 (x);
4121 if (sgn_x == 0)
4122 return SCM_FIXABLE (-yy) ? SCM_MAKINUM (-yy) : scm_long2num (-yy);
4123 else
4124 {
4125 SCM result = scm_i_mkbig ();
ca46fb90 4126
708f22c6
KR		 4127 if (yy >= 0)
4128 mpz_sub_ui (SCM_I_BIG_MPZ (result), SCM_I_BIG_MPZ (x), yy);
4129 else
4130 mpz_add_ui (SCM_I_BIG_MPZ (result), SCM_I_BIG_MPZ (x), -yy);
0aacf84e 4131 scm_remember_upto_here_1 (x);
ca46fb90 4132
0aacf84e
MD		 4133 if ((sgn_x < 0 && (yy > 0)) || ((sgn_x > 0) && yy < 0))
4134 /* we know the result will have to be a bignum */
4135 return result;
4136 else
4137 return scm_i_normbig (result);
4138 }
4139 }
4140 else if (SCM_BIGP (y))
4141 {
4142 int sgn_x = mpz_sgn (SCM_I_BIG_MPZ (x));
4143 int sgn_y = mpz_sgn (SCM_I_BIG_MPZ (y));
4144 SCM result = scm_i_mkbig ();
4145 mpz_sub (SCM_I_BIG_MPZ (result),
4146 SCM_I_BIG_MPZ (x),
4147 SCM_I_BIG_MPZ (y));
4148 scm_remember_upto_here_2 (x, y);
4149 /* we know the result will have to be a bignum */
4150 if ((sgn_x == 1) && (sgn_y == -1))
4151 return result;
4152 if ((sgn_x == -1) && (sgn_y == 1))
4153 return result;
4154 return scm_i_normbig (result);
4155 }
4156 else if (SCM_REALP (y))
4157 {
4158 double result = mpz_get_d (SCM_I_BIG_MPZ (x)) - SCM_REAL_VALUE (y);
4159 scm_remember_upto_here_1 (x);
4160 return scm_make_real (result);
4161 }
4162 else if (SCM_COMPLEXP (y))
4163 {
4164 double real_part = (mpz_get_d (SCM_I_BIG_MPZ (x))
4165 - SCM_COMPLEX_REAL (y));
4166 scm_remember_upto_here_1 (x);
4167 return scm_make_complex (real_part, - SCM_COMPLEX_IMAG (y));
4168 }
f92e85f7
MV		 4169 else if (SCM_FRACTIONP (y))
4170 return scm_make_ratio (scm_difference (scm_product (x, SCM_FRACTION_DENOMINATOR (y)),
4171 SCM_FRACTION_NUMERATOR (y)),
4172 SCM_FRACTION_DENOMINATOR (y));
0aacf84e 4173 else SCM_WTA_DISPATCH_2 (g_difference, x, y, SCM_ARGn, s_difference);
ca46fb90 4174 }
0aacf84e
MD		 4175 else if (SCM_REALP (x))
4176 {
4177 if (SCM_INUMP (y))
4178 return scm_make_real (SCM_REAL_VALUE (x) - SCM_INUM (y));
4179 else if (SCM_BIGP (y))
4180 {
4181 double result = SCM_REAL_VALUE (x) - mpz_get_d (SCM_I_BIG_MPZ (y));
4182 scm_remember_upto_here_1 (x);
4183 return scm_make_real (result);
4184 }
4185 else if (SCM_REALP (y))
4186 return scm_make_real (SCM_REAL_VALUE (x) - SCM_REAL_VALUE (y));
4187 else if (SCM_COMPLEXP (y))
4188 return scm_make_complex (SCM_REAL_VALUE (x) - SCM_COMPLEX_REAL (y),
4189 -SCM_COMPLEX_IMAG (y));
f92e85f7
MV		 4190 else if (SCM_FRACTIONP (y))
4191 return scm_make_real (SCM_REAL_VALUE (x) - scm_i_fraction2double (y));
0aacf84e
MD		 4192 else
4193 SCM_WTA_DISPATCH_2 (g_difference, x, y, SCM_ARGn, s_difference);
98cb6e75 4194 }
0aacf84e
MD		 4195 else if (SCM_COMPLEXP (x))
4196 {
4197 if (SCM_INUMP (y))
4198 return scm_make_complex (SCM_COMPLEX_REAL (x) - SCM_INUM (y),
4199 SCM_COMPLEX_IMAG (x));
4200 else if (SCM_BIGP (y))
4201 {
4202 double real_part = (SCM_COMPLEX_REAL (x)
4203 - mpz_get_d (SCM_I_BIG_MPZ (y)));
4204 scm_remember_upto_here_1 (x);
4205 return scm_make_complex (real_part, SCM_COMPLEX_IMAG (y));
4206 }
4207 else if (SCM_REALP (y))
4208 return scm_make_complex (SCM_COMPLEX_REAL (x) - SCM_REAL_VALUE (y),
4209 SCM_COMPLEX_IMAG (x));
4210 else if (SCM_COMPLEXP (y))
4211 return scm_make_complex (SCM_COMPLEX_REAL (x) - SCM_COMPLEX_REAL (y),
4212 SCM_COMPLEX_IMAG (x) - SCM_COMPLEX_IMAG (y));
f92e85f7
MV		 4213 else if (SCM_FRACTIONP (y))
4214 return scm_make_complex (SCM_COMPLEX_REAL (x) - scm_i_fraction2double (y),
4215 SCM_COMPLEX_IMAG (x));
4216 else
4217 SCM_WTA_DISPATCH_2 (g_difference, x, y, SCM_ARGn, s_difference);
4218 }
4219 else if (SCM_FRACTIONP (x))
4220 {
4221 if (SCM_INUMP (y))
4222 /* a/b - c = (a - cb) / b */
4223 return scm_make_ratio (scm_difference (SCM_FRACTION_NUMERATOR (x),
4224 scm_product(y, SCM_FRACTION_DENOMINATOR (x))),
4225 SCM_FRACTION_DENOMINATOR (x));
4226 else if (SCM_BIGP (y))
4227 return scm_make_ratio (scm_difference (SCM_FRACTION_NUMERATOR (x),
4228 scm_product(y, SCM_FRACTION_DENOMINATOR (x))),
4229 SCM_FRACTION_DENOMINATOR (x));
4230 else if (SCM_REALP (y))
4231 return scm_make_real (scm_i_fraction2double (x) - SCM_REAL_VALUE (y));
4232 else if (SCM_COMPLEXP (y))
4233 return scm_make_complex (scm_i_fraction2double (x) - SCM_COMPLEX_REAL (y),
4234 -SCM_COMPLEX_IMAG (y));
4235 else if (SCM_FRACTIONP (y))
4236 /* a/b - c/d = (ad - bc) / bd */
4237 return scm_make_ratio (scm_difference (scm_product (SCM_FRACTION_NUMERATOR (x), SCM_FRACTION_DENOMINATOR (y)),
4238 scm_product (SCM_FRACTION_NUMERATOR (y), SCM_FRACTION_DENOMINATOR (x))),
4239 scm_product (SCM_FRACTION_DENOMINATOR (x), SCM_FRACTION_DENOMINATOR (y)));
0aacf84e
MD		 4240 else
4241 SCM_WTA_DISPATCH_2 (g_difference, x, y, SCM_ARGn, s_difference);
98cb6e75 4242 }
0aacf84e 4243 else
98cb6e75 4244 SCM_WTA_DISPATCH_2 (g_difference, x, y, SCM_ARG1, s_difference);
0f2d19dd 4245}
c05e97b7 4246#undef FUNC_NAME
0f2d19dd 4247
ca46fb90 4248
9de33deb 4249SCM_GPROC1 (s_product, "*", scm_tc7_asubr, scm_product, g_product);
942e5b91
MG		 4250/* "Return the product of all arguments. If called without arguments,\n"
4251 * "1 is returned."
4252 */
0f2d19dd 4253SCM
6e8d25a6 4254scm_product (SCM x, SCM y)
0f2d19dd 4255{
0aacf84e
MD		 4256 if (SCM_UNBNDP (y))
4257 {
4258 if (SCM_UNBNDP (x))
4259 return SCM_MAKINUM (1L);
4260 else if (SCM_NUMBERP (x))
4261 return x;
4262 else
4263 SCM_WTA_DISPATCH_1 (g_product, x, SCM_ARG1, s_product);
f872b822 4264 }
ca46fb90 4265
0aacf84e
MD		 4266 if (SCM_INUMP (x))
4267 {
4268 long xx;
f4c627b3 4269
0aacf84e
MD		 4270 intbig:
4271 xx = SCM_INUM (x);
f4c627b3 4272
0aacf84e
MD		 4273 switch (xx)
4274 {
ca46fb90
RB		 4275 case 0: return x; break;
4276 case 1: return y; break;
0aacf84e 4277 }
f4c627b3 4278
0aacf84e
MD		 4279 if (SCM_INUMP (y))
4280 {
4281 long yy = SCM_INUM (y);
4282 long kk = xx * yy;
4283 SCM k = SCM_MAKINUM (kk);
4284 if ((kk == SCM_INUM (k)) && (kk / xx == yy))
4285 return k;
4286 else
4287 {
4288 SCM result = scm_i_long2big (xx);
4289 mpz_mul_si (SCM_I_BIG_MPZ (result), SCM_I_BIG_MPZ (result), yy);
4290 return scm_i_normbig (result);
4291 }
4292 }
4293 else if (SCM_BIGP (y))
4294 {
4295 SCM result = scm_i_mkbig ();
4296 mpz_mul_si (SCM_I_BIG_MPZ (result), SCM_I_BIG_MPZ (y), xx);
4297 scm_remember_upto_here_1 (y);
4298 return result;
4299 }
4300 else if (SCM_REALP (y))
4301 return scm_make_real (xx * SCM_REAL_VALUE (y));
4302 else if (SCM_COMPLEXP (y))
4303 return scm_make_complex (xx * SCM_COMPLEX_REAL (y),
4304 xx * SCM_COMPLEX_IMAG (y));
f92e85f7
MV		 4305 else if (SCM_FRACTIONP (y))
4306 return scm_make_ratio (scm_product (x, SCM_FRACTION_NUMERATOR (y)),
4307 SCM_FRACTION_DENOMINATOR (y));
0aacf84e
MD		 4308 else
4309 SCM_WTA_DISPATCH_2 (g_product, x, y, SCM_ARGn, s_product);
f4c627b3 4310 }
0aacf84e
MD		 4311 else if (SCM_BIGP (x))
4312 {
4313 if (SCM_INUMP (y))
4314 {
4315 SCM_SWAP (x, y);
4316 goto intbig;
4317 }
4318 else if (SCM_BIGP (y))
4319 {
4320 SCM result = scm_i_mkbig ();
4321 mpz_mul (SCM_I_BIG_MPZ (result),
4322 SCM_I_BIG_MPZ (x),
4323 SCM_I_BIG_MPZ (y));
4324 scm_remember_upto_here_2 (x, y);
4325 return result;
4326 }
4327 else if (SCM_REALP (y))
4328 {
4329 double result = mpz_get_d (SCM_I_BIG_MPZ (x)) * SCM_REAL_VALUE (y);
4330 scm_remember_upto_here_1 (x);
4331 return scm_make_real (result);
4332 }
4333 else if (SCM_COMPLEXP (y))
4334 {
4335 double z = mpz_get_d (SCM_I_BIG_MPZ (x));
4336 scm_remember_upto_here_1 (x);
4337 return scm_make_complex (z * SCM_COMPLEX_REAL (y),
4338 z * SCM_COMPLEX_IMAG (y));
4339 }
f92e85f7
MV		 4340 else if (SCM_FRACTIONP (y))
4341 return scm_make_ratio (scm_product (x, SCM_FRACTION_NUMERATOR (y)),
4342 SCM_FRACTION_DENOMINATOR (y));
0aacf84e
MD		 4343 else
4344 SCM_WTA_DISPATCH_2 (g_product, x, y, SCM_ARGn, s_product);
f4c627b3 4345 }
0aacf84e
MD		 4346 else if (SCM_REALP (x))
4347 {
4348 if (SCM_INUMP (y))
4349 return scm_make_real (SCM_INUM (y) * SCM_REAL_VALUE (x));
4350 else if (SCM_BIGP (y))
4351 {
4352 double result = mpz_get_d (SCM_I_BIG_MPZ (y)) * SCM_REAL_VALUE (x);
4353 scm_remember_upto_here_1 (y);
4354 return scm_make_real (result);
4355 }
4356 else if (SCM_REALP (y))
4357 return scm_make_real (SCM_REAL_VALUE (x) * SCM_REAL_VALUE (y));
4358 else if (SCM_COMPLEXP (y))
4359 return scm_make_complex (SCM_REAL_VALUE (x) * SCM_COMPLEX_REAL (y),
4360 SCM_REAL_VALUE (x) * SCM_COMPLEX_IMAG (y));
f92e85f7
MV		 4361 else if (SCM_FRACTIONP (y))
4362 return scm_make_real (SCM_REAL_VALUE (x) * scm_i_fraction2double (y));
0aacf84e
MD		 4363 else
4364 SCM_WTA_DISPATCH_2 (g_product, x, y, SCM_ARGn, s_product);
f4c627b3 4365 }
0aacf84e
MD		 4366 else if (SCM_COMPLEXP (x))
4367 {
4368 if (SCM_INUMP (y))
4369 return scm_make_complex (SCM_INUM (y) * SCM_COMPLEX_REAL (x),
4370 SCM_INUM (y) * SCM_COMPLEX_IMAG (x));
4371 else if (SCM_BIGP (y))
4372 {
4373 double z = mpz_get_d (SCM_I_BIG_MPZ (y));
4374 scm_remember_upto_here_1 (y);
76506335
KR		 4375 return scm_make_complex (z * SCM_COMPLEX_REAL (x),
4376 z * SCM_COMPLEX_IMAG (x));
0aacf84e
MD		 4377 }
4378 else if (SCM_REALP (y))
4379 return scm_make_complex (SCM_REAL_VALUE (y) * SCM_COMPLEX_REAL (x),
4380 SCM_REAL_VALUE (y) * SCM_COMPLEX_IMAG (x));
4381 else if (SCM_COMPLEXP (y))
4382 {
4383 return scm_make_complex (SCM_COMPLEX_REAL (x) * SCM_COMPLEX_REAL (y)
4384 - SCM_COMPLEX_IMAG (x) * SCM_COMPLEX_IMAG (y),
4385 SCM_COMPLEX_REAL (x) * SCM_COMPLEX_IMAG (y)
4386 + SCM_COMPLEX_IMAG (x) * SCM_COMPLEX_REAL (y));
4387 }
f92e85f7
MV		 4388 else if (SCM_FRACTIONP (y))
4389 {
4390 double yy = scm_i_fraction2double (y);
4391 return scm_make_complex (yy * SCM_COMPLEX_REAL (x),
4392 yy * SCM_COMPLEX_IMAG (x));
4393 }
4394 else
4395 SCM_WTA_DISPATCH_2 (g_product, x, y, SCM_ARGn, s_product);
4396 }
4397 else if (SCM_FRACTIONP (x))
4398 {
4399 if (SCM_INUMP (y))
4400 return scm_make_ratio (scm_product (y, SCM_FRACTION_NUMERATOR (x)),
4401 SCM_FRACTION_DENOMINATOR (x));
4402 else if (SCM_BIGP (y))
4403 return scm_make_ratio (scm_product (y, SCM_FRACTION_NUMERATOR (x)),
4404 SCM_FRACTION_DENOMINATOR (x));
4405 else if (SCM_REALP (y))
4406 return scm_make_real (scm_i_fraction2double (x) * SCM_REAL_VALUE (y));
4407 else if (SCM_COMPLEXP (y))
4408 {
4409 double xx = scm_i_fraction2double (x);
4410 return scm_make_complex (xx * SCM_COMPLEX_REAL (y),
4411 xx * SCM_COMPLEX_IMAG (y));
4412 }
4413 else if (SCM_FRACTIONP (y))
4414 /* a/b * c/d = ac / bd */
4415 return scm_make_ratio (scm_product (SCM_FRACTION_NUMERATOR (x),
4416 SCM_FRACTION_NUMERATOR (y)),
4417 scm_product (SCM_FRACTION_DENOMINATOR (x),
4418 SCM_FRACTION_DENOMINATOR (y)));
0aacf84e
MD		 4419 else
4420 SCM_WTA_DISPATCH_2 (g_product, x, y, SCM_ARGn, s_product);
f4c627b3 4421 }
0aacf84e 4422 else
f4c627b3 4423 SCM_WTA_DISPATCH_2 (g_product, x, y, SCM_ARG1, s_product);
0f2d19dd
JB		 4424}
4425
0f2d19dd 4426double
6e8d25a6 4427scm_num2dbl (SCM a, const char *why)
f4c627b3 4428#define FUNC_NAME why
0f2d19dd 4429{
0aacf84e 4430 if (SCM_INUMP (a))
0f2d19dd 4431 return (double) SCM_INUM (a);
0aacf84e
MD		 4432 else if (SCM_BIGP (a))
4433 {
4434 double result = mpz_get_d (SCM_I_BIG_MPZ (a));
4435 scm_remember_upto_here_1 (a);
4436 return result;
4437 }
4438 else if (SCM_REALP (a))
f4c627b3 4439 return (SCM_REAL_VALUE (a));
f92e85f7
MV		 4440 else if (SCM_FRACTIONP (a))
4441 return scm_i_fraction2double (a);
0aacf84e 4442 else
f4c627b3 4443 SCM_WRONG_TYPE_ARG (SCM_ARGn, a);
0f2d19dd 4444}
f4c627b3 4445#undef FUNC_NAME
0f2d19dd 4446
7351e207
MV		 4447#if ((defined (HAVE_ISINF) && defined (HAVE_ISNAN)) \
4448 || (defined (HAVE_FINITE) && defined (HAVE_ISNAN)))
4449#define ALLOW_DIVIDE_BY_ZERO
4450/* #define ALLOW_DIVIDE_BY_EXACT_ZERO */
4451#endif
0f2d19dd 4452
ba74ef4e
MV		 4453/* The code below for complex division is adapted from the GNU
4454 libstdc++, which adapted it from f2c's libF77, and is subject to
4455 this copyright: */
4456
4457/****************************************************************
4458Copyright 1990, 1991, 1992, 1993 by AT&T Bell Laboratories and Bellcore.
4459
4460Permission to use, copy, modify, and distribute this software
4461and its documentation for any purpose and without fee is hereby
4462granted, provided that the above copyright notice appear in all
4463copies and that both that the copyright notice and this
4464permission notice and warranty disclaimer appear in supporting
4465documentation, and that the names of AT&T Bell Laboratories or
4466Bellcore or any of their entities not be used in advertising or
4467publicity pertaining to distribution of the software without
4468specific, written prior permission.
4469
4470AT&T and Bellcore disclaim all warranties with regard to this
4471software, including all implied warranties of merchantability
4472and fitness. In no event shall AT&T or Bellcore be liable for
4473any special, indirect or consequential damages or any damages
4474whatsoever resulting from loss of use, data or profits, whether
4475in an action of contract, negligence or other tortious action,
4476arising out of or in connection with the use or performance of
4477this software.
4478****************************************************************/
4479
9de33deb 4480SCM_GPROC1 (s_divide, "/", scm_tc7_asubr, scm_divide, g_divide);
609c3d30
MG		 4481/* Divide the first argument by the product of the remaining
4482 arguments. If called with one argument @var{z1}, 1/@var{z1} is
4483 returned. */
c05e97b7 4484#define FUNC_NAME s_divide
f92e85f7
MV		 4485static SCM
4486scm_i_divide (SCM x, SCM y, int inexact)
0f2d19dd 4487{
f8de44c1
DH		 4488 double a;
4489
0aacf84e
MD		 4490 if (SCM_UNBNDP (y))
4491 {
4492 if (SCM_UNBNDP (x))
4493 SCM_WTA_DISPATCH_0 (g_divide, s_divide);
4494 else if (SCM_INUMP (x))
4495 {
4496 long xx = SCM_INUM (x);
4497 if (xx == 1 || xx == -1)
4498 return x;
7351e207 4499#ifndef ALLOW_DIVIDE_BY_EXACT_ZERO
0aacf84e
MD		 4500 else if (xx == 0)
4501 scm_num_overflow (s_divide);
7351e207 4502#endif
0aacf84e 4503 else
f92e85f7
MV		 4504 {
4505 if (inexact)
4506 return scm_make_real (1.0 / (double) xx);
4507 else return scm_make_ratio (SCM_MAKINUM(1), x);
4508 }
0aacf84e
MD		 4509 }
4510 else if (SCM_BIGP (x))
f92e85f7
MV		 4511 {
4512 if (inexact)
4513 return scm_make_real (1.0 / scm_i_big2dbl (x));
4514 else return scm_make_ratio (SCM_MAKINUM(1), x);
4515 }
0aacf84e
MD		 4516 else if (SCM_REALP (x))
4517 {
4518 double xx = SCM_REAL_VALUE (x);
7351e207 4519#ifndef ALLOW_DIVIDE_BY_ZERO
0aacf84e
MD		 4520 if (xx == 0.0)
4521 scm_num_overflow (s_divide);
4522 else
7351e207 4523#endif
0aacf84e
MD		 4524 return scm_make_real (1.0 / xx);
4525 }
4526 else if (SCM_COMPLEXP (x))
4527 {
4528 double r = SCM_COMPLEX_REAL (x);
4529 double i = SCM_COMPLEX_IMAG (x);
4530 if (r <= i)
4531 {
4532 double t = r / i;
4533 double d = i * (1.0 + t * t);
4534 return scm_make_complex (t / d, -1.0 / d);
4535 }
4536 else
4537 {
4538 double t = i / r;
4539 double d = r * (1.0 + t * t);
4540 return scm_make_complex (1.0 / d, -t / d);
4541 }
4542 }
f92e85f7
MV		 4543 else if (SCM_FRACTIONP (x))
4544 return scm_make_ratio (SCM_FRACTION_DENOMINATOR (x),
4545 SCM_FRACTION_NUMERATOR (x));
0aacf84e
MD		 4546 else
4547 SCM_WTA_DISPATCH_1 (g_divide, x, SCM_ARG1, s_divide);
f8de44c1 4548 }
f8de44c1 4549
0aacf84e
MD		 4550 if (SCM_INUMP (x))
4551 {
4552 long xx = SCM_INUM (x);
4553 if (SCM_INUMP (y))
4554 {
4555 long yy = SCM_INUM (y);
4556 if (yy == 0)
4557 {
7351e207 4558#ifndef ALLOW_DIVIDE_BY_EXACT_ZERO
0aacf84e 4559 scm_num_overflow (s_divide);
7351e207 4560#else
0aacf84e 4561 return scm_make_real ((double) xx / (double) yy);
7351e207 4562#endif
0aacf84e
MD		 4563 }
4564 else if (xx % yy != 0)
f92e85f7
MV		 4565 {
4566 if (inexact)
4567 return scm_make_real ((double) xx / (double) yy);
4568 else return scm_make_ratio (x, y);
4569 }
0aacf84e
MD		 4570 else
4571 {
4572 long z = xx / yy;
4573 if (SCM_FIXABLE (z))
4574 return SCM_MAKINUM (z);
4575 else
4576 return scm_i_long2big (z);
4577 }
f872b822 4578 }
0aacf84e 4579 else if (SCM_BIGP (y))
f92e85f7
MV		 4580 {
4581 if (inexact)
4582 return scm_make_real ((double) xx / scm_i_big2dbl (y));
4583 else return scm_make_ratio (x, y);
4584 }
0aacf84e
MD		 4585 else if (SCM_REALP (y))
4586 {
4587 double yy = SCM_REAL_VALUE (y);
7351e207 4588#ifndef ALLOW_DIVIDE_BY_ZERO
0aacf84e
MD		 4589 if (yy == 0.0)
4590 scm_num_overflow (s_divide);
4591 else
7351e207 4592#endif
0aacf84e 4593 return scm_make_real ((double) xx / yy);
ba74ef4e 4594 }
0aacf84e
MD		 4595 else if (SCM_COMPLEXP (y))
4596 {
4597 a = xx;
4598 complex_div: /* y _must_ be a complex number */
4599 {
4600 double r = SCM_COMPLEX_REAL (y);
4601 double i = SCM_COMPLEX_IMAG (y);
4602 if (r <= i)
4603 {
4604 double t = r / i;
4605 double d = i * (1.0 + t * t);
4606 return scm_make_complex ((a * t) / d, -a / d);
4607 }
4608 else
4609 {
4610 double t = i / r;
4611 double d = r * (1.0 + t * t);
4612 return scm_make_complex (a / d, -(a * t) / d);
4613 }
4614 }
4615 }
f92e85f7
MV		 4616 else if (SCM_FRACTIONP (y))
4617 /* a / b/c = ac / b */
4618 return scm_make_ratio (scm_product (x, SCM_FRACTION_DENOMINATOR (y)),
4619 SCM_FRACTION_NUMERATOR (y));
0aacf84e
MD		 4620 else
4621 SCM_WTA_DISPATCH_2 (g_divide, x, y, SCM_ARGn, s_divide);
f8de44c1 4622 }
0aacf84e
MD		 4623 else if (SCM_BIGP (x))
4624 {
4625 if (SCM_INUMP (y))
4626 {
4627 long int yy = SCM_INUM (y);
4628 if (yy == 0)
4629 {
7351e207 4630#ifndef ALLOW_DIVIDE_BY_EXACT_ZERO
0aacf84e 4631 scm_num_overflow (s_divide);
7351e207 4632#else
0aacf84e
MD		 4633 int sgn = mpz_sgn (SCM_I_BIG_MPZ (x));
4634 scm_remember_upto_here_1 (x);
4635 return (sgn == 0) ? scm_nan () : scm_inf ();
7351e207 4636#endif
0aacf84e
MD		 4637 }
4638 else if (yy == 1)
4639 return x;
4640 else
4641 {
4642 /* FIXME: HMM, what are the relative performance issues here?
4643 We need to test. Is it faster on average to test
4644 divisible_p, then perform whichever operation, or is it
4645 faster to perform the integer div opportunistically and
4646 switch to real if there's a remainder? For now we take the
4647 middle ground: test, then if divisible, use the faster div
4648 func. */
4649
4650 long abs_yy = yy < 0 ? -yy : yy;
4651 int divisible_p = mpz_divisible_ui_p (SCM_I_BIG_MPZ (x), abs_yy);
4652
4653 if (divisible_p)
4654 {
4655 SCM result = scm_i_mkbig ();
4656 mpz_divexact_ui (SCM_I_BIG_MPZ (result), SCM_I_BIG_MPZ (x), abs_yy);
4657 scm_remember_upto_here_1 (x);
4658 if (yy < 0)
4659 mpz_neg (SCM_I_BIG_MPZ (result), SCM_I_BIG_MPZ (result));
4660 return scm_i_normbig (result);
4661 }
4662 else
f92e85f7
MV		 4663 {
4664 if (inexact)
4665 return scm_make_real (scm_i_big2dbl (x) / (double) yy);
4666 else return scm_make_ratio (x, y);
4667 }
0aacf84e
MD		 4668 }
4669 }
4670 else if (SCM_BIGP (y))
4671 {
4672 int y_is_zero = (mpz_sgn (SCM_I_BIG_MPZ (y)) == 0);
4673 if (y_is_zero)
4674 {
ca46fb90 4675#ifndef ALLOW_DIVIDE_BY_EXACT_ZERO
0aacf84e 4676 scm_num_overflow (s_divide);
f872b822 4677#else
0aacf84e
MD		 4678 int sgn = mpz_sgn (SCM_I_BIG_MPZ (x));
4679 scm_remember_upto_here_1 (x);
4680 return (sgn == 0) ? scm_nan () : scm_inf ();
f872b822 4681#endif
0aacf84e
MD		 4682 }
4683 else
4684 {
4685 /* big_x / big_y */
4686 int divisible_p = mpz_divisible_p (SCM_I_BIG_MPZ (x),
4687 SCM_I_BIG_MPZ (y));
4688 if (divisible_p)
4689 {
4690 SCM result = scm_i_mkbig ();
4691 mpz_divexact (SCM_I_BIG_MPZ (result),
4692 SCM_I_BIG_MPZ (x),
4693 SCM_I_BIG_MPZ (y));
4694 scm_remember_upto_here_2 (x, y);
4695 return scm_i_normbig (result);
4696 }
4697 else
4698 {
f92e85f7
MV		 4699 if (inexact)
4700 {
4701 double dbx = mpz_get_d (SCM_I_BIG_MPZ (x));
4702 double dby = mpz_get_d (SCM_I_BIG_MPZ (y));
4703 scm_remember_upto_here_2 (x, y);
4704 return scm_make_real (dbx / dby);
4705 }
4706 else return scm_make_ratio (x, y);
0aacf84e
MD		 4707 }
4708 }
4709 }
4710 else if (SCM_REALP (y))
4711 {
4712 double yy = SCM_REAL_VALUE (y);
7351e207 4713#ifndef ALLOW_DIVIDE_BY_ZERO
0aacf84e
MD		 4714 if (yy == 0.0)
4715 scm_num_overflow (s_divide);
4716 else
7351e207 4717#endif
0aacf84e
MD		 4718 return scm_make_real (scm_i_big2dbl (x) / yy);
4719 }
4720 else if (SCM_COMPLEXP (y))
4721 {
4722 a = scm_i_big2dbl (x);
4723 goto complex_div;
4724 }
f92e85f7
MV		 4725 else if (SCM_FRACTIONP (y))
4726 return scm_make_ratio (scm_product (x, SCM_FRACTION_DENOMINATOR (y)),
4727 SCM_FRACTION_NUMERATOR (y));
0aacf84e
MD		 4728 else
4729 SCM_WTA_DISPATCH_2 (g_divide, x, y, SCM_ARGn, s_divide);
f872b822 4730 }
0aacf84e
MD		 4731 else if (SCM_REALP (x))
4732 {
4733 double rx = SCM_REAL_VALUE (x);
4734 if (SCM_INUMP (y))
4735 {
4736 long int yy = SCM_INUM (y);
7351e207 4737#ifndef ALLOW_DIVIDE_BY_EXACT_ZERO
0aacf84e
MD		 4738 if (yy == 0)
4739 scm_num_overflow (s_divide);
4740 else
7351e207 4741#endif
0aacf84e
MD		 4742 return scm_make_real (rx / (double) yy);
4743 }
4744 else if (SCM_BIGP (y))
4745 {
4746 double dby = mpz_get_d (SCM_I_BIG_MPZ (y));
4747 scm_remember_upto_here_1 (y);
4748 return scm_make_real (rx / dby);
4749 }
4750 else if (SCM_REALP (y))
4751 {
4752 double yy = SCM_REAL_VALUE (y);
7351e207 4753#ifndef ALLOW_DIVIDE_BY_ZERO
0aacf84e
MD		 4754 if (yy == 0.0)
4755 scm_num_overflow (s_divide);
4756 else
7351e207 4757#endif
0aacf84e
MD		 4758 return scm_make_real (rx / yy);
4759 }
4760 else if (SCM_COMPLEXP (y))
4761 {
4762 a = rx;
4763 goto complex_div;
4764 }
f92e85f7
MV		 4765 else if (SCM_FRACTIONP (y))
4766 return scm_make_real (rx / scm_i_fraction2double (y));
0aacf84e
MD		 4767 else
4768 SCM_WTA_DISPATCH_2 (g_divide, x, y, SCM_ARGn, s_divide);
f872b822 4769 }
0aacf84e
MD		 4770 else if (SCM_COMPLEXP (x))
4771 {
4772 double rx = SCM_COMPLEX_REAL (x);
4773 double ix = SCM_COMPLEX_IMAG (x);
4774 if (SCM_INUMP (y))
4775 {
4776 long int yy = SCM_INUM (y);
7351e207 4777#ifndef ALLOW_DIVIDE_BY_EXACT_ZERO
0aacf84e
MD		 4778 if (yy == 0)
4779 scm_num_overflow (s_divide);
4780 else
7351e207 4781#endif
0aacf84e
MD		 4782 {
4783 double d = yy;
4784 return scm_make_complex (rx / d, ix / d);
4785 }
4786 }
4787 else if (SCM_BIGP (y))
4788 {
4789 double dby = mpz_get_d (SCM_I_BIG_MPZ (y));
4790 scm_remember_upto_here_1 (y);
4791 return scm_make_complex (rx / dby, ix / dby);
4792 }
4793 else if (SCM_REALP (y))
4794 {
4795 double yy = SCM_REAL_VALUE (y);
7351e207 4796#ifndef ALLOW_DIVIDE_BY_ZERO
0aacf84e
MD		 4797 if (yy == 0.0)
4798 scm_num_overflow (s_divide);
4799 else
7351e207 4800#endif
0aacf84e
MD		 4801 return scm_make_complex (rx / yy, ix / yy);
4802 }
4803 else if (SCM_COMPLEXP (y))
4804 {
4805 double ry = SCM_COMPLEX_REAL (y);
4806 double iy = SCM_COMPLEX_IMAG (y);
4807 if (ry <= iy)
4808 {
4809 double t = ry / iy;
4810 double d = iy * (1.0 + t * t);
4811 return scm_make_complex ((rx * t + ix) / d, (ix * t - rx) / d);
4812 }
4813 else
4814 {
4815 double t = iy / ry;
4816 double d = ry * (1.0 + t * t);
4817 return scm_make_complex ((rx + ix * t) / d, (ix - rx * t) / d);
4818 }
4819 }
f92e85f7
MV		 4820 else if (SCM_FRACTIONP (y))
4821 {
4822 double yy = scm_i_fraction2double (y);
4823 return scm_make_complex (rx / yy, ix / yy);
4824 }
0aacf84e
MD		 4825 else
4826 SCM_WTA_DISPATCH_2 (g_divide, x, y, SCM_ARGn, s_divide);
f8de44c1 4827 }
f92e85f7
MV		 4828 else if (SCM_FRACTIONP (x))
4829 {
4830 if (SCM_INUMP (y))
4831 {
4832 long int yy = SCM_INUM (y);
4833#ifndef ALLOW_DIVIDE_BY_EXACT_ZERO
4834 if (yy == 0)
4835 scm_num_overflow (s_divide);
4836 else
4837#endif
4838 return scm_make_ratio (SCM_FRACTION_NUMERATOR (x),
4839 scm_product (SCM_FRACTION_DENOMINATOR (x), y));
4840 }
4841 else if (SCM_BIGP (y))
4842 {
4843 return scm_make_ratio (SCM_FRACTION_NUMERATOR (x),
4844 scm_product (SCM_FRACTION_DENOMINATOR (x), y));
4845 }
4846 else if (SCM_REALP (y))
4847 {
4848 double yy = SCM_REAL_VALUE (y);
4849#ifndef ALLOW_DIVIDE_BY_ZERO
4850 if (yy == 0.0)
4851 scm_num_overflow (s_divide);
4852 else
4853#endif
4854 return scm_make_real (scm_i_fraction2double (x) / yy);
4855 }
4856 else if (SCM_COMPLEXP (y))
4857 {
4858 a = scm_i_fraction2double (x);
4859 goto complex_div;
4860 }
4861 else if (SCM_FRACTIONP (y))
4862 return scm_make_ratio (scm_product (SCM_FRACTION_NUMERATOR (x), SCM_FRACTION_DENOMINATOR (y)),
4863 scm_product (SCM_FRACTION_NUMERATOR (y), SCM_FRACTION_DENOMINATOR (x)));
4864 else
4865 SCM_WTA_DISPATCH_2 (g_divide, x, y, SCM_ARGn, s_divide);
4866 }
0aacf84e 4867 else
f8de44c1 4868 SCM_WTA_DISPATCH_2 (g_divide, x, y, SCM_ARG1, s_divide);
0f2d19dd 4869}
f92e85f7
MV		 4870
4871SCM
4872scm_divide (SCM x, SCM y)
4873{
4874 return scm_i_divide (x, y, 0);
4875}
4876
4877static SCM scm_divide2real (SCM x, SCM y)
4878{
4879 return scm_i_divide (x, y, 1);
4880}
c05e97b7 4881#undef FUNC_NAME
0f2d19dd 4882
fa605590 4883
0f2d19dd 4884double
6e8d25a6 4885scm_asinh (double x)
0f2d19dd 4886{
fa605590
KR		 4887#if HAVE_ASINH
4888 return asinh (x);
4889#else
4890#define asinh scm_asinh
f872b822 4891 return log (x + sqrt (x * x + 1));
fa605590 4892#endif
0f2d19dd 4893}
fa605590
KR		 4894SCM_GPROC1 (s_asinh, "$asinh", scm_tc7_dsubr, (SCM (*)()) asinh, g_asinh);
4895/* "Return the inverse hyperbolic sine of @var{x}."
4896 */
0f2d19dd
JB		 4897
4898
0f2d19dd 4899double
6e8d25a6 4900scm_acosh (double x)
0f2d19dd 4901{
fa605590
KR		 4902#if HAVE_ACOSH
4903 return acosh (x);
4904#else
4905#define acosh scm_acosh
f872b822 4906 return log (x + sqrt (x * x - 1));
fa605590 4907#endif
0f2d19dd 4908}
fa605590
KR		 4909SCM_GPROC1 (s_acosh, "$acosh", scm_tc7_dsubr, (SCM (*)()) acosh, g_acosh);
4910/* "Return the inverse hyperbolic cosine of @var{x}."
4911 */
0f2d19dd
JB		 4912
4913
0f2d19dd 4914double
6e8d25a6 4915scm_atanh (double x)
0f2d19dd 4916{
fa605590
KR		 4917#if HAVE_ATANH
4918 return atanh (x);
4919#else
4920#define atanh scm_atanh
f872b822 4921 return 0.5 * log ((1 + x) / (1 - x));
fa605590 4922#endif
0f2d19dd 4923}
fa605590
KR		 4924SCM_GPROC1 (s_atanh, "$atanh", scm_tc7_dsubr, (SCM (*)()) atanh, g_atanh);
4925/* "Return the inverse hyperbolic tangent of @var{x}."
4926 */
0f2d19dd
JB		 4927
4928
f92e85f7
MV		 4929/* XXX - eventually, we should remove this definition of scm_round and
4930 rename scm_round_number to scm_round. Likewise for scm_truncate
4931 and scm_truncate_number.
4932 */
4933
0f2d19dd 4934double
6e8d25a6 4935scm_truncate (double x)
0f2d19dd 4936{
fa605590
KR		 4937#if HAVE_TRUNC
4938 return trunc (x);
4939#else
4940#define trunc scm_truncate
f872b822
MD		 4941 if (x < 0.0)
4942 return -floor (-x);
4943 return floor (x);
fa605590 4944#endif
0f2d19dd 4945}
0f2d19dd 4946
6187f48b
KR		 4947/* scm_round is done using floor(x+0.5) to round to nearest and with
4948 half-way case (ie. when x is an integer plus 0.5) going upwards. Then
4949 half-way cases are identified and adjusted down if the round-upwards
4950 didn't give the desired even integer.
4951
4952 "plus_half == result" identifies a half-way case. If plus_half, which is
4953 x + 0.5, is an integer then x must be an integer plus 0.5.
4954
4955 An odd "result" value is identified with result/2 != floor(result/2).
4956 This is done with plus_half, since that value is ready for use sooner in
4957 a pipelined cpu, and we're already requiring plus_half == result.
4958
4959 Note however that we need to be careful when x is big and already an
4960 integer. In that case "x+0.5" may round to an adjacent integer, causing
4961 us to return such a value, incorrectly. For instance if the hardware is
4962 in the usual default nearest-even rounding, then for x = 0x1FFFFFFFFFFFFF
4963 (ie. 53 one bits) we will have x+0.5 = 0x20000000000000 and that value
4964 returned. Or if the hardware is in round-upwards mode, then other bigger
4965 values like say x == 2^128 will see x+0.5 rounding up to the next higher
4966 representable value, 2^128+2^76 (or whatever), again incorrect.
4967
4968 These bad roundings of x+0.5 are avoided by testing at the start whether
4969 x is already an integer. If it is then clearly that's the desired result
4970 already. And if it's not then the exponent must be small enough to allow
4971 an 0.5 to be represented, and hence added without a bad rounding. */
4972
0f2d19dd 4973double
6e8d25a6 4974scm_round (double x)
0f2d19dd 4975{
6187f48b
KR		 4976 double plus_half, result;
4977
4978 if (x == floor (x))
4979 return x;
4980
4981 plus_half = x + 0.5;
4982 result = floor (plus_half);
0f2d19dd 4983 /* Adjust so that the scm_round is towards even. */
0aacf84e
MD		 4984 return ((plus_half == result && plus_half / 2 != floor (plus_half / 2))
4985 ? result - 1
4986 : result);
0f2d19dd
JB		 4987}
4988
f92e85f7
MV		 4989SCM_DEFINE (scm_truncate_number, "truncate", 1, 0, 0,
4990 (SCM x),
4991 "Round the number @var{x} towards zero.")
4992#define FUNC_NAME s_scm_truncate_number
4993{
4994 if (SCM_FALSEP (scm_negative_p (x)))
4995 return scm_floor (x);
4996 else
4997 return scm_ceiling (x);
4998}
4999#undef FUNC_NAME
5000
5001static SCM exactly_one_half;
5002
5003SCM_DEFINE (scm_round_number, "round", 1, 0, 0,
5004 (SCM x),
5005 "Round the number @var{x} towards the nearest integer. "
5006 "When it is exactly halfway between two integers, "
5007 "round towards the even one.")
5008#define FUNC_NAME s_scm_round_number
5009{
5010 SCM plus_half = scm_sum (x, exactly_one_half);
5011 SCM result = scm_floor (plus_half);
5012 /* Adjust so that the scm_round is towards even. */
5013 if (!SCM_FALSEP (scm_num_eq_p (plus_half, result))
5014 && !SCM_FALSEP (scm_odd_p (result)))
5015 return scm_difference (result, SCM_MAKINUM (1));
5016 else
5017 return result;
5018}
5019#undef FUNC_NAME
5020
5021SCM_PRIMITIVE_GENERIC (scm_floor, "floor", 1, 0, 0,
5022 (SCM x),
5023 "Round the number @var{x} towards minus infinity.")
5024#define FUNC_NAME s_scm_floor
5025{
5026 if (SCM_INUMP (x) || SCM_BIGP (x))
5027 return x;
5028 else if (SCM_REALP (x))
5029 return scm_make_real (floor (SCM_REAL_VALUE (x)));
5030 else if (SCM_FRACTIONP (x))
5031 {
5032 SCM q = scm_quotient (SCM_FRACTION_NUMERATOR (x),
5033 SCM_FRACTION_DENOMINATOR (x));
5034 if (SCM_FALSEP (scm_negative_p (x)))
5035 {
5036 /* For positive x, rounding towards zero is correct. */
5037 return q;
5038 }
5039 else
5040 {
5041 /* For negative x, we need to return q-1 unless x is an
5042 integer. But fractions are never integer, per our
5043 assumptions. */
5044 return scm_difference (q, SCM_MAKINUM (1));
5045 }
5046 }
5047 else
5048 SCM_WTA_DISPATCH_1 (g_scm_floor, x, 1, s_scm_floor);
5049}
5050#undef FUNC_NAME
5051
5052SCM_PRIMITIVE_GENERIC (scm_ceiling, "ceiling", 1, 0, 0,
5053 (SCM x),
5054 "Round the number @var{x} towards infinity.")
5055#define FUNC_NAME s_scm_ceiling
5056{
5057 if (SCM_INUMP (x) || SCM_BIGP (x))
5058 return x;
5059 else if (SCM_REALP (x))
5060 return scm_make_real (ceil (SCM_REAL_VALUE (x)));
5061 else if (SCM_FRACTIONP (x))
5062 {
5063 SCM q = scm_quotient (SCM_FRACTION_NUMERATOR (x),
5064 SCM_FRACTION_DENOMINATOR (x));
5065 if (SCM_FALSEP (scm_positive_p (x)))
5066 {
5067 /* For negative x, rounding towards zero is correct. */
5068 return q;
5069 }
5070 else
5071 {
5072 /* For positive x, we need to return q+1 unless x is an
5073 integer. But fractions are never integer, per our
5074 assumptions. */
5075 return scm_sum (q, SCM_MAKINUM (1));
5076 }
5077 }
5078 else
5079 SCM_WTA_DISPATCH_1 (g_scm_ceiling, x, 1, s_scm_ceiling);
5080}
5081#undef FUNC_NAME
0f2d19dd 5082
14b18ed6 5083SCM_GPROC1 (s_i_sqrt, "$sqrt", scm_tc7_dsubr, (SCM (*)()) sqrt, g_i_sqrt);
942e5b91
MG		 5084/* "Return the square root of the real number @var{x}."
5085 */
14b18ed6 5086SCM_GPROC1 (s_i_abs, "$abs", scm_tc7_dsubr, (SCM (*)()) fabs, g_i_abs);
942e5b91
MG		 5087/* "Return the absolute value of the real number @var{x}."
5088 */
14b18ed6 5089SCM_GPROC1 (s_i_exp, "$exp", scm_tc7_dsubr, (SCM (*)()) exp, g_i_exp);
942e5b91
MG		 5090/* "Return the @var{x}th power of e."
5091 */
14b18ed6 5092SCM_GPROC1 (s_i_log, "$log", scm_tc7_dsubr, (SCM (*)()) log, g_i_log);
b3fcac34 5093/* "Return the natural logarithm of the real number @var{x}."
942e5b91 5094 */
14b18ed6 5095SCM_GPROC1 (s_i_sin, "$sin", scm_tc7_dsubr, (SCM (*)()) sin, g_i_sin);
942e5b91
MG		 5096/* "Return the sine of the real number @var{x}."
5097 */
14b18ed6 5098SCM_GPROC1 (s_i_cos, "$cos", scm_tc7_dsubr, (SCM (*)()) cos, g_i_cos);
942e5b91
MG		 5099/* "Return the cosine of the real number @var{x}."
5100 */
14b18ed6 5101SCM_GPROC1 (s_i_tan, "$tan", scm_tc7_dsubr, (SCM (*)()) tan, g_i_tan);
942e5b91
MG		 5102/* "Return the tangent of the real number @var{x}."
5103 */
14b18ed6 5104SCM_GPROC1 (s_i_asin, "$asin", scm_tc7_dsubr, (SCM (*)()) asin, g_i_asin);
942e5b91
MG		 5105/* "Return the arc sine of the real number @var{x}."
5106 */
14b18ed6 5107SCM_GPROC1 (s_i_acos, "$acos", scm_tc7_dsubr, (SCM (*)()) acos, g_i_acos);
942e5b91
MG		 5108/* "Return the arc cosine of the real number @var{x}."
5109 */
14b18ed6 5110SCM_GPROC1 (s_i_atan, "$atan", scm_tc7_dsubr, (SCM (*)()) atan, g_i_atan);
942e5b91
MG		 5111/* "Return the arc tangent of the real number @var{x}."
5112 */
14b18ed6 5113SCM_GPROC1 (s_i_sinh, "$sinh", scm_tc7_dsubr, (SCM (*)()) sinh, g_i_sinh);
942e5b91
MG		 5114/* "Return the hyperbolic sine of the real number @var{x}."
5115 */
14b18ed6 5116SCM_GPROC1 (s_i_cosh, "$cosh", scm_tc7_dsubr, (SCM (*)()) cosh, g_i_cosh);
942e5b91
MG		 5117/* "Return the hyperbolic cosine of the real number @var{x}."
5118 */
14b18ed6 5119SCM_GPROC1 (s_i_tanh, "$tanh", scm_tc7_dsubr, (SCM (*)()) tanh, g_i_tanh);
942e5b91
MG		 5120/* "Return the hyperbolic tangent of the real number @var{x}."
5121 */
f872b822
MD		 5122
5123struct dpair
5124{
5125 double x, y;
5126};
5127
27c37006
NJ		 5128static void scm_two_doubles (SCM x,
5129 SCM y,
3eeba8d4
JB		 5130 const char *sstring,
5131 struct dpair * xy);
f872b822
MD		 5132
5133static void
27c37006
NJ		 5134scm_two_doubles (SCM x, SCM y, const char *sstring, struct dpair *xy)
5135{
0aacf84e 5136 if (SCM_INUMP (x))
27c37006 5137 xy->x = SCM_INUM (x);
0aacf84e 5138 else if (SCM_BIGP (x))
1be6b49c 5139 xy->x = scm_i_big2dbl (x);
0aacf84e 5140 else if (SCM_REALP (x))
27c37006 5141 xy->x = SCM_REAL_VALUE (x);
f92e85f7
MV		 5142 else if (SCM_FRACTIONP (x))
5143 xy->x = scm_i_fraction2double (x);
0aacf84e 5144 else
27c37006 5145 scm_wrong_type_arg (sstring, SCM_ARG1, x);
98cb6e75 5146
0aacf84e 5147 if (SCM_INUMP (y))
27c37006 5148 xy->y = SCM_INUM (y);
0aacf84e 5149 else if (SCM_BIGP (y))
1be6b49c 5150 xy->y = scm_i_big2dbl (y);
0aacf84e 5151 else if (SCM_REALP (y))
27c37006 5152 xy->y = SCM_REAL_VALUE (y);
f92e85f7
MV		 5153 else if (SCM_FRACTIONP (y))
5154 xy->y = scm_i_fraction2double (y);
0aacf84e 5155 else
27c37006 5156 scm_wrong_type_arg (sstring, SCM_ARG2, y);
0f2d19dd
JB		 5157}
5158
5159
a1ec6916 5160SCM_DEFINE (scm_sys_expt, "$expt", 2, 0, 0,
27c37006
NJ		 5161 (SCM x, SCM y),
5162 "Return @var{x} raised to the power of @var{y}. This\n"
0137a31b 5163 "procedure does not accept complex arguments.")
1bbd0b84 5164#define FUNC_NAME s_scm_sys_expt
0f2d19dd
JB		 5165{
5166 struct dpair xy;
27c37006 5167 scm_two_doubles (x, y, FUNC_NAME, &xy);
f8de44c1 5168 return scm_make_real (pow (xy.x, xy.y));
0f2d19dd 5169}
1bbd0b84 5170#undef FUNC_NAME
0f2d19dd
JB		 5171
5172
a1ec6916 5173SCM_DEFINE (scm_sys_atan2, "$atan2", 2, 0, 0,
27c37006
NJ		 5174 (SCM x, SCM y),
5175 "Return the arc tangent of the two arguments @var{x} and\n"
5176 "@var{y}. This is similar to calculating the arc tangent of\n"
5177 "@var{x} / @var{y}, except that the signs of both arguments\n"
0137a31b
MG		 5178 "are used to determine the quadrant of the result. This\n"
5179 "procedure does not accept complex arguments.")
1bbd0b84 5180#define FUNC_NAME s_scm_sys_atan2
0f2d19dd
JB		 5181{
5182 struct dpair xy;
27c37006 5183 scm_two_doubles (x, y, FUNC_NAME, &xy);
f8de44c1 5184 return scm_make_real (atan2 (xy.x, xy.y));
0f2d19dd 5185}
1bbd0b84 5186#undef FUNC_NAME
0f2d19dd
JB		 5187
5188
a1ec6916 5189SCM_DEFINE (scm_make_rectangular, "make-rectangular", 2, 0, 0,
bb628794 5190 (SCM real, SCM imaginary),
942e5b91
MG		 5191 "Return a complex number constructed of the given @var{real} and\n"
5192 "@var{imaginary} parts.")
1bbd0b84 5193#define FUNC_NAME s_scm_make_rectangular
0f2d19dd
JB		 5194{
5195 struct dpair xy;
bb628794 5196 scm_two_doubles (real, imaginary, FUNC_NAME, &xy);
f8de44c1 5197 return scm_make_complex (xy.x, xy.y);
0f2d19dd 5198}
1bbd0b84 5199#undef FUNC_NAME
0f2d19dd
JB		 5200
5201
5202
a1ec6916 5203SCM_DEFINE (scm_make_polar, "make-polar", 2, 0, 0,
27c37006 5204 (SCM x, SCM y),
942e5b91 5205 "Return the complex number @var{x} * e^(i * @var{y}).")
1bbd0b84 5206#define FUNC_NAME s_scm_make_polar
0f2d19dd
JB		 5207{
5208 struct dpair xy;
6efadd7c 5209 double s, c;
27c37006 5210 scm_two_doubles (x, y, FUNC_NAME, &xy);
6efadd7c
KR		 5211#if HAVE_SINCOS
5212 sincos (xy.y, &s, &c);
5213#else
5214 s = sin (xy.y);
5215 c = cos (xy.y);
5216#endif
5217 return scm_make_complex (xy.x * c, xy.x * s);
0f2d19dd 5218}
1bbd0b84 5219#undef FUNC_NAME
0f2d19dd
JB		 5220
5221
152f82bf 5222SCM_GPROC (s_real_part, "real-part", 1, 0, 0, scm_real_part, g_real_part);
942e5b91
MG		 5223/* "Return the real part of the number @var{z}."
5224 */
0f2d19dd 5225SCM
6e8d25a6 5226scm_real_part (SCM z)
0f2d19dd 5227{
0aacf84e 5228 if (SCM_INUMP (z))
c2ff8ab0 5229 return z;
0aacf84e 5230 else if (SCM_BIGP (z))
c2ff8ab0 5231 return z;
0aacf84e 5232 else if (SCM_REALP (z))
c2ff8ab0 5233 return z;
0aacf84e 5234 else if (SCM_COMPLEXP (z))
c2ff8ab0 5235 return scm_make_real (SCM_COMPLEX_REAL (z));
f92e85f7 5236 else if (SCM_FRACTIONP (z))
2fa2d879 5237 return z;
0aacf84e 5238 else
c2ff8ab0 5239 SCM_WTA_DISPATCH_1 (g_real_part, z, SCM_ARG1, s_real_part);
0f2d19dd
JB		 5240}
5241
5242
152f82bf 5243SCM_GPROC (s_imag_part, "imag-part", 1, 0, 0, scm_imag_part, g_imag_part);
942e5b91
MG		 5244/* "Return the imaginary part of the number @var{z}."
5245 */
0f2d19dd 5246SCM
6e8d25a6 5247scm_imag_part (SCM z)
0f2d19dd 5248{
0aacf84e 5249 if (SCM_INUMP (z))
f872b822 5250 return SCM_INUM0;
0aacf84e 5251 else if (SCM_BIGP (z))
f872b822 5252 return SCM_INUM0;
0aacf84e 5253 else if (SCM_REALP (z))
c2ff8ab0 5254 return scm_flo0;
0aacf84e 5255 else if (SCM_COMPLEXP (z))
c2ff8ab0 5256 return scm_make_real (SCM_COMPLEX_IMAG (z));
f92e85f7
MV		 5257 else if (SCM_FRACTIONP (z))
5258 return SCM_INUM0;
0aacf84e 5259 else
c2ff8ab0 5260 SCM_WTA_DISPATCH_1 (g_imag_part, z, SCM_ARG1, s_imag_part);
0f2d19dd
JB		 5261}
5262
f92e85f7
MV		 5263SCM_GPROC (s_numerator, "numerator", 1, 0, 0, scm_numerator, g_numerator);
5264/* "Return the numerator of the number @var{z}."
5265 */
5266SCM
5267scm_numerator (SCM z)
5268{
5269 if (SCM_INUMP (z))
5270 return z;
5271 else if (SCM_BIGP (z))
5272 return z;
5273 else if (SCM_FRACTIONP (z))
5274 {
5275 scm_i_fraction_reduce (z);
5276 return SCM_FRACTION_NUMERATOR (z);
5277 }
5278 else if (SCM_REALP (z))
5279 return scm_exact_to_inexact (scm_numerator (scm_inexact_to_exact (z)));
5280 else
5281 SCM_WTA_DISPATCH_1 (g_numerator, z, SCM_ARG1, s_numerator);
5282}
5283
5284
5285SCM_GPROC (s_denominator, "denominator", 1, 0, 0, scm_denominator, g_denominator);
5286/* "Return the denominator of the number @var{z}."
5287 */
5288SCM
5289scm_denominator (SCM z)
5290{
5291 if (SCM_INUMP (z))
5292 return SCM_MAKINUM (1);
5293 else if (SCM_BIGP (z))
5294 return SCM_MAKINUM (1);
5295 else if (SCM_FRACTIONP (z))
5296 {
5297 scm_i_fraction_reduce (z);
5298 return SCM_FRACTION_DENOMINATOR (z);
5299 }
5300 else if (SCM_REALP (z))
5301 return scm_exact_to_inexact (scm_denominator (scm_inexact_to_exact (z)));
5302 else
5303 SCM_WTA_DISPATCH_1 (g_denominator, z, SCM_ARG1, s_denominator);
5304}
0f2d19dd 5305
9de33deb 5306SCM_GPROC (s_magnitude, "magnitude", 1, 0, 0, scm_magnitude, g_magnitude);
942e5b91
MG		 5307/* "Return the magnitude of the number @var{z}. This is the same as\n"
5308 * "@code{abs} for real arguments, but also allows complex numbers."
5309 */
0f2d19dd 5310SCM
6e8d25a6 5311scm_magnitude (SCM z)
0f2d19dd 5312{
0aacf84e
MD		 5313 if (SCM_INUMP (z))
5314 {
5315 long int zz = SCM_INUM (z);
5316 if (zz >= 0)
5317 return z;
5318 else if (SCM_POSFIXABLE (-zz))
5319 return SCM_MAKINUM (-zz);
5320 else
5321 return scm_i_long2big (-zz);
5986c47d 5322 }
0aacf84e
MD		 5323 else if (SCM_BIGP (z))
5324 {
5325 int sgn = mpz_sgn (SCM_I_BIG_MPZ (z));
5326 scm_remember_upto_here_1 (z);
5327 if (sgn < 0)
5328 return scm_i_clonebig (z, 0);
5329 else
5330 return z;
5986c47d 5331 }
0aacf84e 5332 else if (SCM_REALP (z))
c2ff8ab0 5333 return scm_make_real (fabs (SCM_REAL_VALUE (z)));
0aacf84e 5334 else if (SCM_COMPLEXP (z))
6efadd7c 5335 return scm_make_real (hypot (SCM_COMPLEX_REAL (z), SCM_COMPLEX_IMAG (z)));
f92e85f7
MV		 5336 else if (SCM_FRACTIONP (z))
5337 {
5338 if (SCM_FALSEP (scm_negative_p (SCM_FRACTION_NUMERATOR (z))))
5339 return z;
5340 return scm_make_ratio (scm_difference (SCM_FRACTION_NUMERATOR (z), SCM_UNDEFINED),
5341 SCM_FRACTION_DENOMINATOR (z));
5342 }
0aacf84e 5343 else
c2ff8ab0 5344 SCM_WTA_DISPATCH_1 (g_magnitude, z, SCM_ARG1, s_magnitude);
0f2d19dd
JB		 5345}
5346
5347
9de33deb 5348SCM_GPROC (s_angle, "angle", 1, 0, 0, scm_angle, g_angle);
942e5b91
MG		 5349/* "Return the angle of the complex number @var{z}."
5350 */
0f2d19dd 5351SCM
6e8d25a6 5352scm_angle (SCM z)
0f2d19dd 5353{
c8ae173e
KR		 5354 /* atan(0,-1) is pi and it'd be possible to have that as a constant like
5355 scm_flo0 to save allocating a new flonum with scm_make_real each time.
5356 But if atan2 follows the floating point rounding mode, then the value
5357 is not a constant. Maybe it'd be close enough though. */
0aacf84e
MD		 5358 if (SCM_INUMP (z))
5359 {
5360 if (SCM_INUM (z) >= 0)
c8ae173e 5361 return scm_flo0;
0aacf84e
MD		 5362 else
5363 return scm_make_real (atan2 (0.0, -1.0));
f872b822 5364 }
0aacf84e
MD		 5365 else if (SCM_BIGP (z))
5366 {
5367 int sgn = mpz_sgn (SCM_I_BIG_MPZ (z));
5368 scm_remember_upto_here_1 (z);
5369 if (sgn < 0)
5370 return scm_make_real (atan2 (0.0, -1.0));
5371 else
c8ae173e 5372 return scm_flo0;
0f2d19dd 5373 }
0aacf84e 5374 else if (SCM_REALP (z))
c8ae173e
KR		 5375 {
5376 if (SCM_REAL_VALUE (z) >= 0)
5377 return scm_flo0;
5378 else
5379 return scm_make_real (atan2 (0.0, -1.0));
5380 }
0aacf84e 5381 else if (SCM_COMPLEXP (z))
f4c627b3 5382 return scm_make_real (atan2 (SCM_COMPLEX_IMAG (z), SCM_COMPLEX_REAL (z)));
f92e85f7
MV		 5383 else if (SCM_FRACTIONP (z))
5384 {
5385 if (SCM_FALSEP (scm_negative_p (SCM_FRACTION_NUMERATOR (z))))
5386 return scm_flo0;
5387 else return scm_make_real (atan2 (0.0, -1.0));
5388 }
0aacf84e 5389 else
f4c627b3 5390 SCM_WTA_DISPATCH_1 (g_angle, z, SCM_ARG1, s_angle);
0f2d19dd
JB		 5391}
5392
5393
3c9a524f
DH		 5394SCM_GPROC (s_exact_to_inexact, "exact->inexact", 1, 0, 0, scm_exact_to_inexact, g_exact_to_inexact);
5395/* Convert the number @var{x} to its inexact representation.\n"
5396 */
5397SCM
5398scm_exact_to_inexact (SCM z)
5399{
5400 if (SCM_INUMP (z))
5401 return scm_make_real ((double) SCM_INUM (z));
5402 else if (SCM_BIGP (z))
5403 return scm_make_real (scm_i_big2dbl (z));
f92e85f7
MV		 5404 else if (SCM_FRACTIONP (z))
5405 return scm_make_real (scm_i_fraction2double (z));
3c9a524f
DH		 5406 else if (SCM_INEXACTP (z))
5407 return z;
5408 else
5409 SCM_WTA_DISPATCH_1 (g_exact_to_inexact, z, 1, s_exact_to_inexact);
5410}
5411
5412
a1ec6916 5413SCM_DEFINE (scm_inexact_to_exact, "inexact->exact", 1, 0, 0,
1bbd0b84 5414 (SCM z),
1e6808ea 5415 "Return an exact number that is numerically closest to @var{z}.")
1bbd0b84 5416#define FUNC_NAME s_scm_inexact_to_exact
0f2d19dd 5417{
0aacf84e 5418 if (SCM_INUMP (z))
f872b822 5419 return z;
0aacf84e 5420 else if (SCM_BIGP (z))
f872b822 5421 return z;
0aacf84e
MD		 5422 else if (SCM_REALP (z))
5423 {
f92e85f7
MV		 5424 if (xisinf (SCM_REAL_VALUE (z)) || xisnan (SCM_REAL_VALUE (z)))
5425 SCM_OUT_OF_RANGE (1, z);
2be24db4 5426 else
f92e85f7
MV		 5427 {
5428 mpq_t frac;
5429 SCM q;
5430
5431 mpq_init (frac);
5432 mpq_set_d (frac, SCM_REAL_VALUE (z));
5433 q = scm_make_ratio (scm_i_mpz2num (mpq_numref (frac)),
5434 scm_i_mpz2num (mpq_denref (frac)));
5435
5436 /* When scm_make_ratio throws, we leak the memory allocated
5437 for frac...
5438 */
5439 mpq_clear (frac);
5440 return q;
5441 }
c2ff8ab0 5442 }
f92e85f7
MV		 5443 else if (SCM_FRACTIONP (z))
5444 return z;
0aacf84e 5445 else
c2ff8ab0 5446 SCM_WRONG_TYPE_ARG (1, z);
0f2d19dd 5447}
1bbd0b84 5448#undef FUNC_NAME
0f2d19dd 5449
f92e85f7
MV		 5450SCM_DEFINE (scm_rationalize, "rationalize", 2, 0, 0,
5451 (SCM x, SCM err),
5452 "Return an exact number that is within @var{err} of @var{x}.")
5453#define FUNC_NAME s_scm_rationalize
5454{
5455 if (SCM_INUMP (x))
5456 return x;
5457 else if (SCM_BIGP (x))
5458 return x;
5459 else if ((SCM_REALP (x)) || SCM_FRACTIONP (x))
5460 {
5461 /* Use continued fractions to find closest ratio. All
5462 arithmetic is done with exact numbers.
5463 */
5464
5465 SCM ex = scm_inexact_to_exact (x);
5466 SCM int_part = scm_floor (ex);
5467 SCM tt = SCM_MAKINUM (1);
5468 SCM a1 = SCM_MAKINUM (0), a2 = SCM_MAKINUM (1), a = SCM_MAKINUM (0);
5469 SCM b1 = SCM_MAKINUM (1), b2 = SCM_MAKINUM (0), b = SCM_MAKINUM (0);
5470 SCM rx;
5471 int i = 0;
5472
5473 if (!SCM_FALSEP (scm_num_eq_p (ex, int_part)))
5474 return ex;
5475
5476 ex = scm_difference (ex, int_part); /* x = x-int_part */
5477 rx = scm_divide (ex, SCM_UNDEFINED); /* rx = 1/x */
5478
5479 /* We stop after a million iterations just to be absolutely sure
5480 that we don't go into an infinite loop. The process normally
5481 converges after less than a dozen iterations.
5482 */
5483
5484 err = scm_abs (err);
5485 while (++i < 1000000)
5486 {
5487 a = scm_sum (scm_product (a1, tt), a2); /* a = a1*tt + a2 */
5488 b = scm_sum (scm_product (b1, tt), b2); /* b = b1*tt + b2 */
5489 if (SCM_FALSEP (scm_zero_p (b)) && /* b != 0 */
5490 SCM_FALSEP
5491 (scm_gr_p (scm_abs (scm_difference (ex, scm_divide (a, b))),
5492 err))) /* abs(x-a/b) <= err */
02164269
MV		 5493 {
5494 SCM res = scm_sum (int_part, scm_divide (a, b));
5495 if (SCM_FALSEP (scm_exact_p (x))
5496 || SCM_FALSEP (scm_exact_p (err)))
5497 return scm_exact_to_inexact (res);
5498 else
5499 return res;
5500 }
f92e85f7
MV		 5501 rx = scm_divide (scm_difference (rx, tt), /* rx = 1/(rx - tt) */
5502 SCM_UNDEFINED);
5503 tt = scm_floor (rx); /* tt = floor (rx) */
5504 a2 = a1;
5505 b2 = b1;
5506 a1 = a;
5507 b1 = b;
5508 }
5509 scm_num_overflow (s_scm_rationalize);
5510 }
5511 else
5512 SCM_WRONG_TYPE_ARG (1, x);
5513}
5514#undef FUNC_NAME
5515
87617347 5516/* if you need to change this, change test-num2integral.c as well */
ee33d62a 5517#if SCM_SIZEOF_LONG_LONG != 0
1be6b49c
ML		 5518# ifndef LLONG_MAX
5519# define ULLONG_MAX ((unsigned long long) (-1))
5520# define LLONG_MAX ((long long) (ULLONG_MAX >> 1))
5521# define LLONG_MIN (~LLONG_MAX)
5522# endif
f872b822 5523#endif
0f2d19dd 5524
3d2e8ceb
MV		 5525/* Parameters for creating integer conversion routines.
5526
5527 Define the following preprocessor macros before including
5528 "libguile/num2integral.i.c":
5529
5530 NUM2INTEGRAL - the name of the function for converting from a
ca46fb90
RB		 5531 Scheme object to the integral type. This function will be
5532 defined when including "num2integral.i.c".
3d2e8ceb
MV		 5533
5534 INTEGRAL2NUM - the name of the function for converting from the
ca46fb90 5535 integral type to a Scheme object. This function will be defined.
3d2e8ceb
MV		 5536
5537 INTEGRAL2BIG - the name of an internal function that createas a
ca46fb90
RB		 5538 bignum from the integral type. This function will be defined.
5539 The name should start with "scm_i_".
5540
5541 ITYPE - the name of the integral type.
5542
9dd023e1
MV		 5543 UNSIGNED - Define this to 1 when ITYPE is an unsigned type. Define
5544 it to 0 otherwise.
ca46fb90
RB		 5545
5546 UNSIGNED_ITYPE - the name of the the unsigned variant of the
5547 integral type. If you don't define this, it defaults to
5548 "unsigned ITYPE" for signed types and simply "ITYPE" for unsigned
5549 ones.
5550
5551 SIZEOF_ITYPE - an expression giving the size of the integral type
5552 in bytes. This expression must be computable by the
5553 preprocessor. (SIZEOF_FOO values are calculated by configure.in
5554 for common types).
5555
3d2e8ceb
MV		 5556*/
5557
1be6b49c
ML		 5558#define NUM2INTEGRAL scm_num2short
5559#define INTEGRAL2NUM scm_short2num
5560#define INTEGRAL2BIG scm_i_short2big
ca46fb90 5561#define UNSIGNED 0
1be6b49c 5562#define ITYPE short
3d2e8ceb 5563#define SIZEOF_ITYPE SIZEOF_SHORT
1be6b49c
ML		 5564#include "libguile/num2integral.i.c"
5565
5566#define NUM2INTEGRAL scm_num2ushort
5567#define INTEGRAL2NUM scm_ushort2num
5568#define INTEGRAL2BIG scm_i_ushort2big
ca46fb90 5569#define UNSIGNED 1
1be6b49c 5570#define ITYPE unsigned short
ca46fb90 5571#define SIZEOF_ITYPE SIZEOF_UNSIGNED_SHORT
1be6b49c
ML		 5572#include "libguile/num2integral.i.c"
5573
5574#define NUM2INTEGRAL scm_num2int
5575#define INTEGRAL2NUM scm_int2num
5576#define INTEGRAL2BIG scm_i_int2big
ca46fb90 5577#define UNSIGNED 0
1be6b49c 5578#define ITYPE int
3d2e8ceb 5579#define SIZEOF_ITYPE SIZEOF_INT
1be6b49c
ML		 5580#include "libguile/num2integral.i.c"
5581
5582#define NUM2INTEGRAL scm_num2uint
5583#define INTEGRAL2NUM scm_uint2num
5584#define INTEGRAL2BIG scm_i_uint2big
ca46fb90 5585#define UNSIGNED 1
1be6b49c 5586#define ITYPE unsigned int
ca46fb90 5587#define SIZEOF_ITYPE SIZEOF_UNSIGNED_INT
1be6b49c
ML		 5588#include "libguile/num2integral.i.c"
5589
5590#define NUM2INTEGRAL scm_num2long
5591#define INTEGRAL2NUM scm_long2num
5592#define INTEGRAL2BIG scm_i_long2big
ca46fb90 5593#define UNSIGNED 0
1be6b49c 5594#define ITYPE long
3d2e8ceb 5595#define SIZEOF_ITYPE SIZEOF_LONG
1be6b49c
ML		 5596#include "libguile/num2integral.i.c"
5597
5598#define NUM2INTEGRAL scm_num2ulong
5599#define INTEGRAL2NUM scm_ulong2num
5600#define INTEGRAL2BIG scm_i_ulong2big
ca46fb90 5601#define UNSIGNED 1
1be6b49c 5602#define ITYPE unsigned long
ca46fb90 5603#define SIZEOF_ITYPE SIZEOF_UNSIGNED_LONG
1be6b49c
ML		 5604#include "libguile/num2integral.i.c"
5605
1be6b49c
ML		 5606#define NUM2INTEGRAL scm_num2ptrdiff
5607#define INTEGRAL2NUM scm_ptrdiff2num
5608#define INTEGRAL2BIG scm_i_ptrdiff2big
ca46fb90 5609#define UNSIGNED 0
ee33d62a 5610#define ITYPE scm_t_ptrdiff
3d2e8ceb 5611#define UNSIGNED_ITYPE size_t
ee33d62a 5612#define SIZEOF_ITYPE SCM_SIZEOF_SCM_T_PTRDIFF
1be6b49c
ML		 5613#include "libguile/num2integral.i.c"
5614
5615#define NUM2INTEGRAL scm_num2size
5616#define INTEGRAL2NUM scm_size2num
5617#define INTEGRAL2BIG scm_i_size2big
ca46fb90 5618#define UNSIGNED 1
1be6b49c 5619#define ITYPE size_t
3d2e8ceb 5620#define SIZEOF_ITYPE SIZEOF_SIZE_T
1be6b49c 5621#include "libguile/num2integral.i.c"
0f2d19dd 5622
ee33d62a 5623#if SCM_SIZEOF_LONG_LONG != 0
1cc91f1b 5624
caf08e65
MV		 5625#ifndef ULONG_LONG_MAX
5626#define ULONG_LONG_MAX (~0ULL)
5627#endif
5628
1be6b49c
ML		 5629#define NUM2INTEGRAL scm_num2long_long
5630#define INTEGRAL2NUM scm_long_long2num
5631#define INTEGRAL2BIG scm_i_long_long2big
ca46fb90 5632#define UNSIGNED 0
1be6b49c 5633#define ITYPE long long
3d2e8ceb 5634#define SIZEOF_ITYPE SIZEOF_LONG_LONG
1be6b49c
ML		 5635#include "libguile/num2integral.i.c"
5636
5637#define NUM2INTEGRAL scm_num2ulong_long
5638#define INTEGRAL2NUM scm_ulong_long2num
5639#define INTEGRAL2BIG scm_i_ulong_long2big
ca46fb90 5640#define UNSIGNED 1
1be6b49c 5641#define ITYPE unsigned long long
ca46fb90 5642#define SIZEOF_ITYPE SIZEOF_UNSIGNED_LONG_LONG
1be6b49c 5643#include "libguile/num2integral.i.c"
0f2d19dd 5644
ee33d62a 5645#endif /* SCM_SIZEOF_LONG_LONG != 0 */
caf08e65 5646
5437598b
MD		 5647#define NUM2FLOAT scm_num2float
5648#define FLOAT2NUM scm_float2num
5649#define FTYPE float
5650#include "libguile/num2float.i.c"
5651
5652#define NUM2FLOAT scm_num2double
5653#define FLOAT2NUM scm_double2num
5654#define FTYPE double
5655#include "libguile/num2float.i.c"
5656
1be6b49c 5657#ifdef GUILE_DEBUG
caf08e65 5658
6063dc1d
SJ		 5659#ifndef SIZE_MAX
5660#define SIZE_MAX ((size_t) (-1))
5661#endif
5662#ifndef PTRDIFF_MIN
5663#define PTRDIFF_MIN \
b4fb7de8
RB		 5664 ((scm_t_ptrdiff) ((scm_t_ptrdiff) 1 \
5665 << ((sizeof (scm_t_ptrdiff) * SCM_CHAR_BIT) - 1)))
6063dc1d
SJ		 5666#endif
5667#ifndef PTRDIFF_MAX
5668#define PTRDIFF_MAX (~ PTRDIFF_MIN)
5669#endif
5670
0aacf84e
MD		 5671#define CHECK(type, v) \
5672 do \
5673 { \
5674 if ((v) != scm_num2##type (scm_##type##2num (v), 1, "check_sanity")) \
5675 abort (); \
5676 } \
5677 while (0)
caf08e65 5678
1be6b49c
ML		 5679static void
5680check_sanity ()
5681{
5682 CHECK (short, 0);
5683 CHECK (ushort, 0U);
5684 CHECK (int, 0);
5685 CHECK (uint, 0U);
5686 CHECK (long, 0L);
5687 CHECK (ulong, 0UL);
5688 CHECK (size, 0);
5689 CHECK (ptrdiff, 0);
5690
5691 CHECK (short, -1);
5692 CHECK (int, -1);
5693 CHECK (long, -1L);
5694 CHECK (ptrdiff, -1);
5695
5696 CHECK (short, SHRT_MAX);
5697 CHECK (short, SHRT_MIN);
5698 CHECK (ushort, USHRT_MAX);
5699 CHECK (int, INT_MAX);
5700 CHECK (int, INT_MIN);
5701 CHECK (uint, UINT_MAX);
5702 CHECK (long, LONG_MAX);
5703 CHECK (long, LONG_MIN);
5704 CHECK (ulong, ULONG_MAX);
5705 CHECK (size, SIZE_MAX);
5706 CHECK (ptrdiff, PTRDIFF_MAX);
5707 CHECK (ptrdiff, PTRDIFF_MIN);
0f2d19dd 5708
ee33d62a 5709#if SCM_SIZEOF_LONG_LONG != 0
1be6b49c
ML		 5710 CHECK (long_long, 0LL);
5711 CHECK (ulong_long, 0ULL);
1be6b49c 5712 CHECK (long_long, -1LL);
1be6b49c
ML		 5713 CHECK (long_long, LLONG_MAX);
5714 CHECK (long_long, LLONG_MIN);
5715 CHECK (ulong_long, ULLONG_MAX);
5716#endif
0f2d19dd
JB		 5717}
5718
b10586f0
ML		 5719#undef CHECK
5720
5721#define CHECK \
5722 scm_internal_catch (SCM_BOOL_T, check_body, &data, check_handler, &data); \
5723 if (!SCM_FALSEP (data)) abort();
5724
5725static SCM
5726check_body (void *data)
5727{
5728 SCM num = *(SCM *) data;
5729 scm_num2ulong (num, 1, NULL);
5730
5731 return SCM_UNSPECIFIED;
5732}
5733
5734static SCM
5735check_handler (void *data, SCM tag, SCM throw_args)
5736{
5737 SCM *num = (SCM *) data;
5738 *num = SCM_BOOL_F;
5739
5740 return SCM_UNSPECIFIED;
5741}
5742
5743SCM_DEFINE (scm_sys_check_number_conversions, "%check-number-conversions", 0, 0, 0,
b4e15479 5744 (void),
b10586f0
ML		 5745 "Number conversion sanity checking.")
5746#define FUNC_NAME s_scm_sys_check_number_conversions
5747{
5748 SCM data = SCM_MAKINUM (-1);
5749 CHECK;
5750 data = scm_int2num (INT_MIN);
5751 CHECK;
5752 data = scm_ulong2num (ULONG_MAX);
5753 data = scm_difference (SCM_INUM0, data);
5754 CHECK;
5755 data = scm_ulong2num (ULONG_MAX);
5756 data = scm_sum (SCM_MAKINUM (1), data); data = scm_difference (SCM_INUM0, data);
5757 CHECK;
5758 data = scm_int2num (-10000); data = scm_product (data, data); data = scm_product (data, data);
5759 CHECK;
5760
5761 return SCM_UNSPECIFIED;
5762}
5763#undef FUNC_NAME
5764
1be6b49c 5765#endif
0f2d19dd 5766
0f2d19dd
JB		 5767void
5768scm_init_numbers ()
0f2d19dd 5769{
0b799eea
MV		 5770 int i;
5771
713a4259
KR		 5772 mpz_init_set_si (z_negative_one, -1);
5773
a261c0e9
DH		 5774 /* It may be possible to tune the performance of some algorithms by using
5775 * the following constants to avoid the creation of bignums. Please, before
5776 * using these values, remember the two rules of program optimization:
5777 * 1st Rule: Don't do it. 2nd Rule (experts only): Don't do it yet. */
86d31dfe
MV		 5778 scm_c_define ("most-positive-fixnum",
5779 SCM_MAKINUM (SCM_MOST_POSITIVE_FIXNUM));
5780 scm_c_define ("most-negative-fixnum",
5781 SCM_MAKINUM (SCM_MOST_NEGATIVE_FIXNUM));
a261c0e9 5782
f3ae5d60
MD		 5783 scm_add_feature ("complex");
5784 scm_add_feature ("inexact");
5986c47d 5785 scm_flo0 = scm_make_real (0.0);
0b799eea
MV		 5786
5787 /* determine floating point precision */
5788 for(i=2; i <= SCM_MAX_DBL_RADIX; ++i)
5789 {
5790 init_dblprec(&scm_dblprec[i-2],i);
5791 init_fx_radix(fx_per_radix[i-2],i);
5792 }
f872b822 5793#ifdef DBL_DIG
0b799eea
MV		 5794 /* hard code precision for base 10 if the preprocessor tells us to... */
5795 scm_dblprec[10-2] = (DBL_DIG > 20) ? 20 : DBL_DIG;
5796#endif
1be6b49c
ML		 5797
5798#ifdef GUILE_DEBUG
5799 check_sanity ();
5800#endif
f92e85f7
MV		 5801
5802 exactly_one_half = scm_permanent_object (scm_divide (SCM_MAKINUM (1),
5803 SCM_MAKINUM (2)));
a0599745 5804#include "libguile/numbers.x"
0f2d19dd 5805}
89e00824
ML		 5806
5807/*
5808 Local Variables:
5809 c-file-style: "gnu"
5810 End:
5811*/
Unnamed repository; edit this file 'description' to name the repository.