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