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