number.c 84 KB
Newer Older
1
/*                                                      -*- coding: utf-8 -*-
eg's avatar
eg committed
2
 *
3
 * n u m b e r . c      -- Numbers management
eg's avatar
eg committed
4
 *
5
 * Copyright © 1993-2018 Erick Gallesio - I3S-CNRS/ESSI <eg@essi.fr>
eg's avatar
eg committed
6 7 8 9 10
 *
 * This program is free software; you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation; either version 2 of the License, or
 * (at your option) any later version.
Erick's avatar
.  
Erick committed
11
 *
eg's avatar
eg committed
12 13 14 15
 * This program is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
Erick's avatar
.  
Erick committed
16
 *
eg's avatar
eg committed
17 18
 * You should have received a copy of the GNU General Public License
 * along with this program; if not, write to the Free Software
Erick's avatar
.  
Erick committed
19
 * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307,
eg's avatar
eg committed
20 21 22 23 24
 * USA.
 *
 *
 *           Author: Erick Gallesio [eg@kaolin.unice.fr]
 *    Creation date: 12-May-1993 10:34
25
 * Last file update: 13-Jun-2018 17:37 (eg)
eg's avatar
eg committed
26 27 28
 */


Erick's avatar
.  
Erick committed
29 30
/* workaround for bad optimisations done in glibc2.1. Thanks to  Andreas Jaeger
 * <aj@suse.de> for it
eg's avatar
eg committed
31
 */
Erick's avatar
.  
Erick committed
32
#define __NO_MATH_INLINES
eg's avatar
eg committed
33 34 35

#include <math.h>
#include <ctype.h>
36
#include <locale.h>
eg's avatar
eg committed
37 38 39 40 41 42 43 44 45 46
#include "stklos.h"

#undef sinc
#if defined(__linux__) && defined(__alpha__)
#  include <signal.h>
#endif


/* Real precision */
static int real_precision = REAL_FORMAT_SIZE;
Erick's avatar
.  
Erick committed
47
static int log10_maxint;
eg's avatar
eg committed
48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65

#define FINITE_REALP(n) ((REAL_VAL(n) != minus_inf) && (REAL_VAL(n) != plus_inf))



/* Declaration of bignums. This is done here instead of stklos.h to avoid
 * to expose the file "gmp.h" in "stklos.h" which is the interface users
 * see to access all the sytem (note that we can also use our version which
 * can be different of the one which is sytem installed, and resolve conflits
 * could be hard).
 */
#include <gmp.h>

struct bignum_obj {
  stk_header header;
  mpz_t val;
};

66
#define BIGNUM_VAL(p)   (((struct bignum_obj *) (p))->val)
eg's avatar
eg committed
67 68 69 70


/*==============================================================================*/

71
#define MY_PI           3.1415926535897932384626433832795029L  /* pi */
eg's avatar
eg committed
72

73 74
#define BIGNUM_FITS_INTEGER(_bn) (mpz_cmp_si((_bn), INT_MIN_VAL) >= 0 &&        \
                                  mpz_cmp_si((_bn), INT_MAX_VAL) <= 0)
eg's avatar
eg committed
75
#define LONG_FITS_INTEGER(_l)    (INT_MIN_VAL <= (_l) && (_l) <= INT_MAX_VAL)
76
#define TYPEOF(n)                (INTP(n)? tc_integer: STYPE(n))
eg's avatar
eg committed
77 78 79 80


#define MINUS_INF "-inf.0"
#define PLUS_INF  "+inf.0"
81
#define NaN       "nan.0"
eg's avatar
eg committed
82 83 84 85 86 87 88 89 90 91 92 93 94 95 96


/* Special IEEE values */
static double plus_inf, minus_inf;
double STk_NaN;

/**** forward declarations ****/
static type_cell convert(SCM *px, SCM *py);

static int zerop(SCM n);
static int negativep(SCM n);
static int positivep(SCM n);
static int isexactp(SCM z);
static SCM gcd2(SCM o1, SCM o2);

97 98
EXTERN_PRIMITIVE("make-rectangular", make_rectangular, subr2, (SCM r, SCM i));
EXTERN_PRIMITIVE("real-part", real_part, subr1, (SCM z));
eg's avatar
eg committed
99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120
EXTERN_PRIMITIVE("magnitude", magnitude, subr1, (SCM z));
EXTERN_PRIMITIVE("angle", angle, subr1, (SCM z));
EXTERN_PRIMITIVE("sqrt", sqrt, subr1, (SCM z));
EXTERN_PRIMITIVE("exact->inexact", ex2inex, subr1, (SCM z));
EXTERN_PRIMITIVE("inexact->exact", inex2ex, subr1, (SCM z));

#define add2 STk_add2
#define mul2 STk_mul2
#define div2 STk_div2
#define sub2 STk_sub2
#define absolute STk_abs
#define exact2inexact STk_ex2inex
#define inexact2exact STk_inex2ex


static SCM int_quotient(SCM o1, SCM o2);
static SCM my_cos(SCM z);
static SCM my_sin(SCM z);



/******************************************************************************
Erick's avatar
.  
Erick committed
121 122
 *
 * Utilities
eg's avatar
eg committed
123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160
 *
 ******************************************************************************/
static void error_bad_number(SCM n)
{
  STk_error("~S is a bad number", n);
}

static void error_at_least_1(void)
{
  STk_error("expects at least one argument");
}

static void error_not_on_a_complex(SCM n)
{
  STk_error("cannot be computed on the complex number ~S", n);
}

static void error_cannot_operate(char *operation, SCM o1, SCM o2)
{
  STk_error("cannot perform %s on ~S and ~S", operation, o1, o2);
}

static void error_divide_by_0(SCM n)
{
  STk_error("cannot divide ~S by 0", n);
}

static void error_incorrect_radix(SCM r)
{
  STk_error("base must be 2, 8, 10 or 16. It was ~S", r);
}


/*
<doc EXT real-precision
 * (real-precision)
 * (real-precision value)
 *
Erick's avatar
.  
Erick committed
161 162
 * This parameter object permits to change the default precision used
 * to print real numbers.
eg's avatar
eg committed
163 164 165 166 167 168 169 170 171 172 173 174
 * @lisp
 * (real-precision)        => 15
 * (define f 0.123456789)
 * (display f)             @print{} 0.123456789
 * (real-precision 3)
 * (display f)             @print{} 0.123
 * @end lisp
doc>
*/
static SCM real_precision_conv(SCM value)
{
  long precision = STk_integer_value(value);
Erick's avatar
.  
Erick committed
175

eg's avatar
eg committed
176 177
  if (precision <= 0 || precision > 50)
    STk_error("real precision must be an integer in ]0 50]. It was ~S",
178
              value);
eg's avatar
eg committed
179 180 181 182 183 184
  real_precision = (int) precision;
  return value;
}


/******************************************************************************
Erick's avatar
.  
Erick committed
185 186
 *
 * Constructor Functions
eg's avatar
eg committed
187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228
 *
 ******************************************************************************/
static Inline SCM Cmake_complex(SCM r, SCM i)
{
  SCM z;

  NEWCELL(z, complex);
  COMPLEX_REAL(z) = r;
  COMPLEX_IMAG(z) = i;
  return z;
}

static Inline SCM make_complex(SCM r, SCM i)
{
  return (zerop(i)) ? r : Cmake_complex(r, i);
}

static Inline SCM make_polar(SCM a, SCM m)
{
  return make_complex(mul2(a, my_cos(m)), mul2(a, my_sin(m)));
}




static Inline SCM Cmake_rational(SCM n, SCM d)
{
  SCM z;

  NEWCELL(z, rational);
  RATIONAL_NUM(z) = n;
  RATIONAL_DEN(z) = d;
  return z;
}


static SCM make_rational(SCM n, SCM d)
{
  SCM gcd;

  if (zerop(d))
    STk_error("cannot make rational with null denominator");
229

eg's avatar
eg committed
230 231 232 233 234
  /* Always keep sign in the denominator */
  if (negativep(d)) {
    n = mul2(n, MAKE_INT(-1));
    d = mul2(d, MAKE_INT(-1));
  }
Erick's avatar
.  
Erick committed
235

eg's avatar
eg committed
236 237 238 239 240 241 242 243 244 245 246
  /* Simplify rational */
  gcd = gcd2(n, d);
  if (gcd != MAKE_INT(1)) {
    if (d == gcd) return int_quotient(n, gcd);
    n = int_quotient(n, gcd);
    d = int_quotient(d, gcd);
  }

  /* Make rational if denominator is not 1 (if n is small, it is not a bignum) */
  if (d ==  MAKE_INT(1))
    return n;
Erick's avatar
.  
Erick committed
247
  else
eg's avatar
eg committed
248 249 250 251 252
    return Cmake_rational(n, d);
}


/******************************************************************************
Erick's avatar
.  
Erick committed
253 254
 *
 * Types declaration
eg's avatar
eg committed
255 256 257 258 259 260 261 262 263 264 265
 *
 ******************************************************************************/
static void print_bignum(SCM n, SCM port, int mode)
{
  char *s;

  s = STk_must_malloc_atomic(mpz_sizeinbase(BIGNUM_VAL(n), 10) + 2);
  mpz_get_str(s, 10, BIGNUM_VAL(n));
  STk_puts(s, port);
  STk_free(s);
}
Erick's avatar
.  
Erick committed
266

eg's avatar
eg committed
267 268 269 270 271 272 273

static void print_rational(SCM n, SCM port, int mode)
{
  STk_print(RATIONAL_NUM(n), port, mode);
  STk_putc('/', port);
  STk_print(RATIONAL_DEN(n), port, mode);
}
Erick's avatar
.  
Erick committed
274

eg's avatar
eg committed
275 276 277
static void print_complex(SCM n, SCM port, int mode)
{
  SCM imag = COMPLEX_IMAG(n);
Erick's avatar
.  
Erick committed
278

eg's avatar
eg committed
279 280 281 282 283 284 285 286 287 288 289 290 291 292 293 294 295 296 297 298 299 300 301 302 303 304
  STk_print(COMPLEX_REAL(n), port, mode);
  if (positivep(imag))
    STk_putc('+', port);
  STk_print(imag, port, mode);
  STk_putc('i', port);
}



static struct extended_type_descr xtype_bignum = {
  "bignum",
  print_bignum
};

static struct extended_type_descr xtype_complex = {
  "complex",
  print_complex
};

static struct extended_type_descr xtype_rational = {
  "rational",
  print_rational
};


/******************************************************************************
Erick's avatar
.  
Erick committed
305 306
 *
 * Conversion Functions
eg's avatar
eg committed
307 308 309 310 311 312
 *
 ******************************************************************************/

static Inline SCM long2scheme_bignum(long x)
{
  SCM z;
Erick's avatar
.  
Erick committed
313

eg's avatar
eg committed
314 315 316 317 318 319 320 321 322 323 324 325 326 327 328
  NEWCELL(z, bignum);
  mpz_init_set_si(BIGNUM_VAL(z), x);
  return z;
}


static Inline SCM long2integer(long x)
{
  return (INT_MIN_VAL<=x && x<=INT_MAX_VAL) ?  MAKE_INT(x): long2scheme_bignum(x);
}


static Inline SCM double2real(double x)
{
  SCM z;
Erick's avatar
.  
Erick committed
329

eg's avatar
eg committed
330 331 332 333 334
  NEWCELL(z, real);
  REAL_VAL(z) = x;
  return z;
}

335
static SCM double2integer(double n)     /* small or big depending of n's size */
eg's avatar
eg committed
336 337 338 339 340 341 342 343 344 345
{
  int i, j;
  size_t size = 20;
  char *tmp = NULL;
  SCM z;

  /* Try first to convert n to a long */
  if (((double) INT_MIN_VAL <= n) && (n <= (double) INT_MAX_VAL))
    return MAKE_INT((long) n);

Erick's avatar
.  
Erick committed
346
  /* n doesn't fit in a long => build a bignum. THIS IS VERY INEFFICIENT */
eg's avatar
eg committed
347 348 349 350 351 352 353 354 355 356 357 358 359 360 361 362 363 364 365 366 367 368 369 370 371 372 373 374 375 376 377 378 379 380 381 382 383 384 385
  tmp = STk_must_malloc(size);
  i = 0;
  if (n < 0.0) { tmp[i++] = '-'; n = -n; }
  do {
    if (i >= size) tmp = STk_must_realloc(tmp, size *= 2);
    tmp[i++] = (int) fmod(n, (double) 10) + '0';
    n = floor(n / 10.0);
  }
  while (n > 0.0);
  tmp[i] = 0;

  /* Reverse the content of string tmp */
  for (i=i-1, j=(tmp[0]=='-'); i > j; i--, j++) {
    char c = tmp[i];
    tmp[i] = tmp[j];
    tmp[j] = c;
  }

  /* tmp contains a textual representation of n. Convert it to a bignum */
  z = STk_Cstr2number(tmp, 10L);
  if (tmp) STk_free(tmp);
  return z;
}

static SCM double2rational(double d)
{
  double fraction, i;
  SCM int_part, num, den, res;
  int negative = 0;

  if (d < 0.0) { negative = 1; d = -d; }
  fraction = modf(d, &i);
  int_part = double2integer(i);

  if (!fraction) {
    res = int_part;
  } else {
    num = MAKE_INT(0);
    den = MAKE_INT(1);
Erick's avatar
.  
Erick committed
386

eg's avatar
eg committed
387 388
    while (fraction) {
      num      = mul2(num, MAKE_INT(2));
Erick's avatar
.  
Erick committed
389
      den      = mul2(den, MAKE_INT(2));
eg's avatar
eg committed
390 391
      fraction = modf(ldexp(fraction, 1), &i);
      if (i)
392
        num = add2(num, MAKE_INT(1));
eg's avatar
eg committed
393 394 395
    }
    res = add2(int_part, div2(num, den));
  }
Erick's avatar
.  
Erick committed
396

eg's avatar
eg committed
397 398 399 400 401 402 403 404 405 406 407 408 409 410 411 412 413 414 415 416 417 418 419 420
  return negative? mul2(res, MAKE_INT(-1)): res;
}

static Inline SCM bignum2scheme_bignum(mpz_t n)
{
  SCM z;

  NEWCELL(z, bignum);
  mpz_init_set(BIGNUM_VAL(z), n);
  return z;
}

static Inline SCM bignum2integer(mpz_t n)
{
  return MAKE_INT(mpz_get_si(n));
}

static Inline SCM bignum2number(mpz_t n)  /* => int or bignum */
{
  return (BIGNUM_FITS_INTEGER(n)) ? bignum2integer(n): bignum2scheme_bignum(n);
}

static Inline double bignum2double(mpz_t n)
{
Erick's avatar
.  
Erick committed
421 422
  /* I do not use the function mpz_get_d since it gives an unspecified value
   * when converting a number which is +inf or -inf
eg's avatar
eg committed
423 424 425 426
   */
 char *s = STk_must_malloc_atomic(mpz_sizeinbase(n, 10) + 2);
 return atof(mpz_get_str(s, 10, n));
}
Erick's avatar
.  
Erick committed
427

eg's avatar
eg committed
428 429 430 431 432 433 434 435 436 437 438 439 440 441 442 443 444 445 446 447 448 449 450 451 452 453 454 455 456 457 458 459 460 461 462 463 464 465 466 467 468 469 470 471 472

static Inline double scheme_bignum2double(SCM b)
{
  return bignum2double(BIGNUM_VAL(b));
}


static Inline SCM scheme_bignum2real(SCM bn)
{
  return double2real(scheme_bignum2double(bn));
}


static double rational2double(SCM r)
{
  SCM num = RATIONAL_NUM(r);
  SCM den = RATIONAL_DEN(r);

  switch (convert(&num, &den)) {
    case tc_integer: return ((double) INT_VAL(num)) / ((double) INT_VAL(den));
    case tc_bignum:  return scheme_bignum2double(num) / scheme_bignum2double(den);
    default:         STk_panic("bad rational ~S", r);
  }
  return 0.0; /* never reached */
}

static Inline SCM rational2real(SCM r)
{
  return double2real(rational2double(r));
}

static Inline SCM real2integer(SCM r)
{
  double v = REAL_VAL(r);

  if (floor(v) != v) {
    /* This is not an inexact integer (weak test) */
    STk_error("bad number (~s) in an integer division", r);
  }
  return double2integer(v);
}

void STk_double2Cstr(char *buffer, double n)
{
  sprintf(buffer, "%.*g", real_precision, n);
Erick's avatar
.  
Erick committed
473
  if (strchr(buffer, '.') == NULL && strchr(buffer, 'e') == NULL)
eg's avatar
eg committed
474 475 476 477 478 479 480 481 482 483 484 485
    strcat(buffer, ".0");
  /* Treat special cases of +nan.0 and +inf.0 */
  if (isalpha(buffer[0])) {
    if (strcmp(buffer, "inf.0") == 0) strcpy(buffer, "+inf.0");
    if (strcmp(buffer, "nan.0") == 0) strcpy(buffer, "+nan.0");
  }
}

/* Convert a number to a C-string. Result must be freed if != from buffer */
static char *number2Cstr(SCM n, long base, char buffer[])
{
  char *s = buffer;
Erick's avatar
.  
Erick committed
486

eg's avatar
eg committed
487 488 489
  switch (TYPEOF(n)) {
    case tc_integer:
      {
490 491 492 493 494 495 496 497 498 499 500 501 502 503 504 505 506 507
        long tmp, val = INT_VAL(n);
        int u;

        if (val < 0) {
          val  = -val;
          *s++ = '-';
        }
        /* Find how much digit we need */
        for (s++, tmp=val; tmp >= base; tmp /= base) s++;

        *s = '\0'; tmp = val;
        do {
          u = tmp % base;
          *(--s) = u + ((u < 10) ? '0' : 'a'-10);
          tmp   /= base;
        }
        while (tmp);
        return buffer;
eg's avatar
eg committed
508 509 510 511 512 513
      }
    case tc_bignum:
      s = STk_must_malloc(mpz_sizeinbase(BIGNUM_VAL(n), base) + 2);
      s = mpz_get_str(s, base, BIGNUM_VAL(n));
      return s;
    case tc_rational:
Erick's avatar
.  
Erick committed
514
      {
515 516 517 518 519 520 521 522
        char *s1, *s2, *s3, tmp[100];

        s1 = number2Cstr(RATIONAL_NUM(n), base, buffer);
        s2 = number2Cstr(RATIONAL_DEN(n), base, tmp);
        s3 = STk_must_malloc(strlen(s1) + strlen(s2) + 2);
        sprintf(s3, "%s/%s", s1, s2);
        if (s2!=tmp) STk_free(s2); /*buffer will event. be deallocated by caller*/
        return s3;
eg's avatar
eg committed
523 524 525
      }
    case tc_complex:
      {
526 527 528 529 530 531 532 533
        char *s1, *s2, *s3, tmp[100];

        s1 = number2Cstr(COMPLEX_REAL(n), base, buffer);
        s2 = number2Cstr(COMPLEX_IMAG(n), base, tmp);
        s3 = STk_must_malloc(strlen(s1) + strlen(s2) + 3);
        sprintf(s3, "%s%s%si", s1, ((*s2 == '-') ? "": "+"), s2);
        if (s2!=tmp) STk_free(s2); /*buffer will event. be deallocated by caller*/
        return s3;
eg's avatar
eg committed
534 535 536 537 538 539 540 541 542 543 544 545 546 547 548 549 550 551 552 553 554
      }
    case tc_real:
      if (base != 10) STk_error("base must be 10 for this number", n);
      STk_double2Cstr(buffer, REAL_VAL(n));
      return buffer;

    default: return NULL; /* never reached */
  }
}


/*===== The general conversion routine ==== */

static type_cell convert(SCM *px, SCM *py)
{
  SCM x = *px;
  SCM y = *py;

  if (TYPEOF(x)==TYPEOF(y)) return(TYPEOF(x)); /* avoid testing on current cases */
  switch (TYPEOF(x)) {
    case tc_complex:
555 556 557 558 559 560 561 562 563
            switch (TYPEOF(y)) {
              case tc_complex: /*already done */
              case tc_real:
              case tc_rational:
              case tc_bignum:
              case tc_integer:  *py = Cmake_complex(y, MAKE_INT(0)); break;
              default:          error_bad_number(y);                 break;
            }
            break;
eg's avatar
eg committed
564
    case tc_real:
565 566 567 568 569 570 571 572 573
            switch (TYPEOF(y)) {
              case tc_complex:  *px = Cmake_complex(x, MAKE_INT(0));    break;
              case tc_real:    /*already done */ ;                      break;
              case tc_rational: *py = rational2real(y);                 break;
              case tc_bignum:   *py = scheme_bignum2real(y);            break;
              case tc_integer:  *py = double2real((double) INT_VAL(y)); break;
              default:          error_bad_number(y);                    break;
            }
            break;
eg's avatar
eg committed
574
    case tc_rational:
575 576 577 578 579 580 581 582 583
            switch (TYPEOF(y)) {
              case tc_complex:  *px = Cmake_complex(x, MAKE_INT(0));   break;
              case tc_real:     *px = rational2real(x);                break;
              case tc_rational: /*already done */ ;                    break;
              case tc_bignum:   /* no break */
              case tc_integer:  *py = Cmake_rational(y , MAKE_INT(1)); break;
              default:          error_bad_number(y);                   break;
            }
            break;
eg's avatar
eg committed
584
    case tc_bignum:
585 586 587 588 589 590 591 592 593
            switch (TYPEOF(y)) {
              case tc_complex:  *px = Cmake_complex(x, MAKE_INT(0));    break;
              case tc_real:     *px = scheme_bignum2real(x);            break;
              case tc_rational: *px = Cmake_rational(x , MAKE_INT(1));  break;
              case tc_bignum:   /* already done */                      break;
              case tc_integer:  *py = long2scheme_bignum(INT_VAL(y));   break;
              default:          error_bad_number(y);                    break;
            }
            break;
Erick's avatar
.  
Erick committed
594
    case tc_integer:
595 596 597 598 599 600 601 602 603
            switch (TYPEOF(y)) {
              case tc_complex:  *px = Cmake_complex(x, MAKE_INT(0));    break;
              case tc_real:     *px = double2real((double) INT_VAL(x)); break;
              case tc_rational: *px = Cmake_rational(x,  MAKE_INT(1));  break;
              case tc_bignum:   *px = long2scheme_bignum(INT_VAL(x));   break;
              case tc_integer:  /* already done */                      break;
              default:          error_bad_number(y);                    break;
            }
            break;
eg's avatar
eg committed
604 605 606 607 608 609 610 611 612 613 614 615 616 617 618 619 620 621 622 623 624 625 626 627 628 629 630 631 632 633 634 635 636 637 638 639 640 641 642 643 644 645
    default: error_bad_number(x);
  }
  return TYPEOF(*px);
}


long STk_integer_value(SCM x) /* Returns LONG_MIN if not representable as long */
{
  if (INTP(x)) return INT_VAL(x);
  if (BIGNUMP(x)) {
    mpz_t *v = &BIGNUM_VAL(x);
    if (mpz_cmp_si(*v, LONG_MIN) > 0 && mpz_cmp_si(*v, LONG_MAX) <= 0)
      return  mpz_get_si(*v);
  }
  return LONG_MIN;
}

unsigned long STk_uinteger_value(SCM x) /* Returns ULONG_MAX if not an ulong */
{
  if (INTP(x) && (x > 0)) return INT_VAL(x); /* sign(INTEGER_VAL(x)) == sign(x) */
  if (BIGNUMP(x)) {
    mpz_t *v = &BIGNUM_VAL(x);
    if (mpz_cmp_ui(*v, 0) >= 0 && mpz_cmp_ui(*v, ULONG_MAX) < 0)
      return mpz_get_ui(*v);
  }
  return ULONG_MAX;
}


SCM STk_long2integer(long n)
{
  return long2integer(n);
}


SCM STk_ulong2integer(unsigned long n)
{
  if (n <= INT_MAX_VAL) {  /* n  >= 0 since it is an ulong */
    return MAKE_INT(n);
  }
  else {
    SCM z;
Erick's avatar
.  
Erick committed
646

eg's avatar
eg committed
647 648 649 650 651 652 653 654 655 656 657 658
    NEWCELL(z, bignum);
    mpz_init_set_ui(BIGNUM_VAL(z), n);
    return z;
  }
}


long STk_integer2int32(SCM n, int *overflow)
{
  *overflow = 0;

  if (INTP(n)) {
659
#if (LONG_MAX == INT32_MAX)       /* longs are on 32 bits */
eg's avatar
eg committed
660
    return INT_VAL(n);
661
#else                             /* longs are more than 32 bits (probably 64) */
eg's avatar
eg committed
662
    long val = INT_VAL(n);
Erick's avatar
.  
Erick committed
663

eg's avatar
eg committed
664 665 666 667 668 669 670 671 672 673 674 675 676 677 678 679 680 681 682 683 684
    if ((- INT32_MAX - 1) <= val && val < INT32_MAX)
      return val;
    else {
      *overflow = 1;
      return 0;
    }
#endif
  }
  if (BIGNUMP(n)) {
    mpz_t *v = &BIGNUM_VAL(n);
    if (mpz_cmp_si(*v, (- INT32_MAX - 1)) >= 0 && mpz_cmp_si(*v, INT32_MAX) <= 0)
      return mpz_get_si(*v);
  }
  *overflow = 1;
  return 0;
}


unsigned long STk_integer2uint32(SCM n, int *overflow)
{
  *overflow = 0;
Erick's avatar
.  
Erick committed
685

eg's avatar
eg committed
686 687
  if (INTP(n)) {
    long val = INT_VAL(n);
Erick's avatar
.  
Erick committed
688

eg's avatar
eg committed
689
    if (val >= 0) {
690
#if (ULONG_MAX == UINT32_MAX)   /* unsigned longs are on 32 bits */
eg's avatar
eg committed
691
      return (unsigned long) INT_VAL(n);
692
#else                           /* longs are more than 32 bits (probably 64) */
eg's avatar
eg committed
693
      if (val < UINT32_MAX)
694
        return val;
eg's avatar
eg committed
695
      else {
696 697
        *overflow = 1;
        return 0;
eg's avatar
eg committed
698 699 700 701 702 703 704 705 706 707 708 709 710 711 712 713 714 715 716 717 718
      }
#endif
    }
  }
  if (BIGNUMP(n)) {
    mpz_t *v = &BIGNUM_VAL(n);
    if (mpz_cmp_ui(*v, 0) >= 0 && mpz_cmp_ui(*v, UINT32_MAX) <= 0)
      return mpz_get_ui(*v);
  }
  *overflow = 1;
  return 0;
}



SCM STk_double2real(double d)
{
  return double2real(d); /* use the inlined function */
}


719
double STk_number2double(SCM n) /* returns NaN if not convertible */
eg's avatar
eg committed
720 721 722 723 724 725 726 727 728 729
{
  switch (TYPEOF(n)) {
    case tc_real:     return REAL_VAL(n);
    case tc_rational: return REAL_VAL(rational2real(n));
    case tc_bignum:   return REAL_VAL(scheme_bignum2real(n));
    case tc_integer:  return (double) INT_VAL(n);
    default:          return STk_NaN;
  }
}

Erick's avatar
.  
Erick committed
730

eg's avatar
eg committed
731 732

/******************************************************************************
Erick's avatar
.  
Erick committed
733 734
 *
 * Utilities
eg's avatar
eg committed
735 736 737 738 739 740
 *
 ******************************************************************************/

static long do_compare(SCM x, SCM y)
{
  double d1=0, d2=0;
Erick's avatar
.  
Erick committed
741

eg's avatar
eg committed
742
  switch (TYPEOF(x)) {
Erick's avatar
.  
Erick committed
743
    case tc_real:
744 745 746 747 748 749 750 751 752 753
            switch (TYPEOF(y)) {
              case tc_complex:  goto general_diff;
              case tc_real:     d1 = REAL_VAL(x); d2 = REAL_VAL(y);
                                goto double_diff;
              case tc_rational:
              case tc_bignum:   goto general_diff;
              case tc_integer:  d1 = REAL_VAL(x); d2 =  INT_VAL(y);
                                goto double_diff;
              default:          break;
            }
eg's avatar
eg committed
754
    case tc_integer:
755 756 757 758 759 760 761 762 763
            switch (TYPEOF(y)) {
              case tc_complex:  goto general_diff;
              case tc_real:     d1 = INT_VAL(x); d2 = REAL_VAL(y);
                                goto double_diff;
              case tc_rational:
              case tc_bignum:   goto general_diff;
              case tc_integer:  return (INT_VAL(x) - INT_VAL(y));
              default:          break;
            }
eg's avatar
eg committed
764 765 766
    case tc_complex:
    case tc_rational:
    case tc_bignum:
767 768 769 770 771 772 773 774
            switch (TYPEOF(y)) {
              case tc_complex:
              case tc_real:
              case tc_rational:
              case tc_bignum:
              case tc_integer:  goto general_diff;
              default:          break;
            }
eg's avatar
eg committed
775 776 777 778 779 780 781
    default: break;
  }
  /* if we are here, it s that x and y cannot be compared */
  STk_error("comparison between ~S and ~S impossible", x,  y);
double_diff:
  if (isnan(d1) && isnan(d2))
    return 0;
Erick's avatar
.  
Erick committed
782
  return (d1 == d2) ? 0 : ((d1 < d2)?  -1 : 1);
eg's avatar
eg committed
783 784 785
general_diff:
  {
    SCM d = sub2(x, y);
Erick's avatar
.  
Erick committed
786

eg's avatar
eg committed
787 788 789 790 791 792 793
    if (zerop(d)) return 0;
    /* complex numbers cannot be compared => return always 1 */
    return COMPLEXP(d) ? 1 : (negativep(d) ? -1: 1);
  }
}


Erick's avatar
.  
Erick committed
794
static SCM int_quotient(SCM x, SCM y)
eg's avatar
eg committed
795 796 797 798 799 800 801
/* Specialized version for rationals. Accepts only integer or bignums as params */
{
  mpz_t q, r;

  if (INTP(x)) {
    if (INTP(y))
      return MAKE_INT(INT_VAL(x)/INT_VAL(y));
Erick's avatar
.  
Erick committed
802
    else
803
      x = long2scheme_bignum(INT_VAL(x));
eg's avatar
eg committed
804 805 806 807 808 809
  } else {
    if (INTP(y))
      y = long2scheme_bignum(INT_VAL(y));
  }
  /* Here x and y are both bignum */
  mpz_init(q); mpz_init(r);
Erick's avatar
.  
Erick committed
810
  mpz_tdiv_qr(q, r, BIGNUM_VAL(x), BIGNUM_VAL(y));
eg's avatar
eg committed
811 812 813 814 815 816 817 818 819 820 821 822 823 824 825
  return bignum2number(q);
}

static int digitp(char c, long base)
{
  c = ('0' <= c && c <= '9') ? c - '0':
      ('a' <= c && c <= 'f') ? c - 'a' + 10:
      ('A' <= c && c <= 'F') ? c - 'A' + 10:
      (c == '#')             ? 0           :
      100;
  return (c < base);
}


/******************************************************************************
Erick's avatar
.  
Erick committed
826 827
 *
 * Number parser
eg's avatar
eg committed
828 829 830 831 832 833 834 835 836 837 838 839 840 841
 *
 ******************************************************************************/

static SCM compute_exact_real(char *s, char *p1, char *p2, char *p3, char *p4)
{
  SCM int_part, fract_part, exp_part;
  mpz_t tmp;

  mpz_init(tmp);
  int_part   = MAKE_INT(0);
  fract_part = MAKE_INT(0);
  exp_part   = MAKE_INT(1);

  /* Representation of the given number (number is '\0' terminated)
Erick's avatar
.  
Erick committed
842
   *
eg's avatar
eg committed
843 844
   *        +xxxxxxxxxxxxxxxxx.yyyyyyyyyyyyyE+zzzzz
   *        ^                 ^^            ^^
845
   *        |                 ||            ||
eg's avatar
eg committed
846 847 848
   *        +-str          p1-++-p2      p3-++-p4
   */

Erick's avatar
.  
Erick committed
849 850
  /* patch the given string so that splitting the various parts of the number
   * is easy
eg's avatar
eg committed
851 852 853 854
   */
  if (p1) *p1 = '\0';
  if (p3) *p3 = '\0';

855
  if (p1) {             /* compute integer part */
eg's avatar
eg committed
856 857 858
    if (mpz_init_set_str(tmp, s, 10L) < 0) return STk_false;
    int_part = bignum2number(tmp);
  }
Erick's avatar
.  
Erick committed
859

860
  if (p3 > p2) {        /* compute decimal part as a rational 0.12 => 6/5 */
eg's avatar
eg committed
861 862 863 864 865 866 867 868 869 870
    SCM num, den;

    if (mpz_init_set_str(tmp, p2, 10L) < 0) return STk_false;
    num = bignum2number(tmp);

    mpz_ui_pow_ui(tmp, 10UL, strlen(p2));
    den = bignum2number(tmp);

    fract_part = make_rational(num, den);
  }
Erick's avatar
.  
Erick committed
871

872
  if (p4) {             /* compute exposant as a rational 3 => 1000, -3 => 1/1000 */
eg's avatar
eg committed
873
    long expo;
Erick's avatar
.  
Erick committed
874

eg's avatar
eg committed
875 876 877 878 879 880 881 882 883
    expo = atoi(p4);
    if (expo > 0) {
      mpz_ui_pow_ui(tmp, 10UL, expo);
      exp_part = bignum2number(tmp);
    } else {
      mpz_ui_pow_ui(tmp, 10UL, -expo);
      exp_part = div2(MAKE_INT(1), bignum2number(tmp));
    }
  }
Erick's avatar
.  
Erick committed
884

eg's avatar
eg committed
885 886 887 888 889 890 891 892 893 894 895 896 897 898 899 900 901
  /* now return (int_part + fract_part) * exp_part */
  return mul2(add2(int_part, fract_part), exp_part);
}

static SCM read_integer_or_real(char *str, long base, char exact_flag, char **end)
{
  int adigit=0, isint=1;
  char saved_char = '\0', *p = str, *p1, *p2, *p3, *p4;
  SCM res;

  /* see function compute_exact_real for the meaning of these pointers */
  p1 = p2 = p3 = p4 = NULL;

  if (*p == '-' || *p == '+') p+=1;
  if (*p == '#') return STk_false;
  while(digitp(*p, base)) { p+=1; adigit=1; if (*p == '#') isint = 0; }

902
  if (adigit) p1 = p;           /* p1 = end of integral part */
eg's avatar
eg committed
903 904 905 906 907 908 909 910 911 912 913 914 915 916 917 918 919 920 921 922 923 924 925

  if (*p=='.') {
    isint = 0; p += 1;
    p2 = p;
    while(digitp(*p, base)) { p+=1; adigit=1; }
    p3 = p;
  }

  if (!adigit) return STk_false;

  if (*p && strchr("eEsSfFdDlL", *p)) {
    isint = 0;
    p += 1;
    p4 = p;
    if (*p == '-' || *p == '+') p+=1;
    if (!digitp(*p, base)) return STk_false;
    p+=1;
    while (digitp(*p, base)) p+=1;
  }
  if (*p) {
    /* Patch the end of the number with a '\0' (will be restored on exit) */
    saved_char = *p;
    *p = '\0';
Erick's avatar
.  
Erick committed
926
  }
eg's avatar
eg committed
927 928

  if (isint) {
Erick's avatar
.  
Erick committed
929 930 931 932 933
    /* We are sure to have an integer. Read it as a bignum and see if we can
     * convert it in smallnum after that. Small numbers (those with few
     * digits expressed in base 10) are not read as bignums.
     * This optimisation is easily missed (e.g. 000000000000000001 will be
     * read as a bignum), but it avoids allocation for current numbers
934
     * represented in an usual form.
eg's avatar
eg committed
935 936
     */
    mpz_t n;
Erick's avatar
.  
Erick committed
937

eg's avatar
eg committed
938
    if (*str == '+') str+=1; /* mpz_init_set_str doesn't recognize +xyz !!! */
939
    if (strlen(str) <= log10_maxint && base == 10) {
Erick's avatar
.  
Erick committed
940 941
      long num = atol(str);

eg's avatar
eg committed
942 943
      res = (exact_flag == 'i') ? double2real((double) num): MAKE_INT(num);
    }
944 945
    else {
      if (mpz_init_set_str(n, str, base) < 0) {
946 947
        /* Bad syntax for a bignum */
        res = STk_false;
948
      } else if (BIGNUM_FITS_INTEGER(n)) {
949 950
        /* Can be represented as a short integer */
        long num = mpz_get_si(n);
Erick's avatar
.  
Erick committed
951

952
        res = (exact_flag == 'i') ? double2real((double) num): MAKE_INT(num);
953
      } else {
954 955 956
        /* It's a bignum */
        res = bignum2scheme_bignum(n);
        if (exact_flag == 'i') res = scheme_bignum2real(res);
957 958 959
      }
      mpz_clear(n);
    }
eg's avatar
eg committed
960 961 962 963
  } else {
    /* Number is a float */
    if (base == 10) {
      /* Replace sharp signs by 0 */
Erick's avatar
.  
Erick committed
964
      for(p=str; *p; p++)
965 966 967 968 969
        switch (*p) {
          case '#': *p = '0'; break;
          case 's': case 'S': case 'f': case 'F':
          case 'd': case 'D': case 'l': case 'L': *p = 'e';
        }
eg's avatar
eg committed
970
      if (exact_flag == 'e') {
971
        res = compute_exact_real(str, p1, p2, p3, p4);
eg's avatar
eg committed
972
      } else {
973
        res = double2real(strtod(str, &p));
eg's avatar
eg committed
974 975 976 977 978 979 980
      }
    }
    else
      res = STk_false;
  }

  if (saved_char) *p = saved_char;  /* character which ended the number */
981
  *end = p;                         /* position of last analysed character */
eg's avatar
eg committed
982 983 984 985 986 987 988
  return res;
}


static SCM read_rational(SCM num, char *str, long base, char exact_flag, char **end)
{
  SCM den;
Erick's avatar
.  
Erick committed
989

eg's avatar
eg committed
990 991
  den = read_integer_or_real(str, base, exact_flag, end);
  if (den == STk_false) return STk_false;
Erick's avatar
.  
Erick committed
992

eg's avatar
eg committed
993 994 995 996 997 998
  if ((TYPEOF(num) == tc_integer || TYPEOF(num) == tc_bignum) &&
      (TYPEOF(den) == tc_integer || TYPEOF(den) == tc_bignum))
    return make_rational(num, den);
  else
    STk_error("cannot make rational with ~S and ~S", num, den);

999
  return STk_false;             /* never reached */
eg's avatar
eg committed
1000 1001 1002 1003 1004 1005 1006 1007 1008 1009 1010 1011 1012 1013 1014 1015 1016 1017 1018 1019 1020 1021
}

SCM STk_Cstr2number(char *str, long base)
{
  int i, exact, radix, polar, is_signed;
  char *p = str;
  SCM num1, num2;

  is_signed = 0;

  if ((*str == '-' || *str == '+')) {
    is_signed = 1;

    if (isalpha(str[1])) {
      /* Treat special values "+inf.0" -inf.0 and "NaN" as well as +i and -i */
      if (strcmp(str, MINUS_INF)==0) return double2real(minus_inf);
      if (strcmp(str, PLUS_INF)==0)  return double2real(plus_inf);
      if (strcmp(str+1, NaN)==0)     return double2real(STk_NaN);
      if (strcmp(str, "+i")==0)      return make_complex(MAKE_INT(0), MAKE_INT(+1));
      if (strcmp(str, "-i")==0)      return make_complex(MAKE_INT(0), MAKE_INT(-1));
    }
  }
Erick's avatar
.  
Erick committed
1022

eg's avatar
eg committed
1023 1024 1025 1026 1027 1028
  exact = ' ', radix = 0;
  for (i = 0; i < 2; i++) {
    if (*p == '#') {
      p += 1;
      switch (*p++) {
        case 'e': if (exact == ' ') { exact = 'e'; break; }  else return STk_false;
1029 1030 1031 1032 1033 1034
        case 'i': if (exact == ' ') { exact = 'i'; break; }  else return STk_false;
        case 'b': if (!radix) {base = 2;  radix = 1; break;} else return STk_false;
        case 'o': if (!radix) {base = 8;  radix = 1; break;} else return STk_false;
        case 'd': if (!radix) {base = 10; radix = 1; break;} else return STk_false;
        case 'x': if (!radix) {base = 16; radix = 1; break;} else return STk_false;
        default:  return STk_false;
eg's avatar
eg committed
1035 1036 1037 1038 1039
      }
      str += 2;
    }
    if (*p != '#') break;
  }
Erick's avatar
.  
Erick committed
1040

eg's avatar
eg committed
1041 1042
  num1 = read_integer_or_real(p, base, exact, &p);
  if (num1 == STk_false) return STk_false;
Erick's avatar
.  
Erick committed
1043

eg's avatar
eg committed
1044 1045 1046 1047 1048 1049 1050
  if (*p == '/')
    num1 = read_rational(num1, p+1, base, exact, &p);

  if ((*p == '+') || (*p == '-') || (*p == '@')) {
    /* Start to read a complex number */
    if (*p == '+' && p[1] == 'i') {
      p   += 2;
1051
      num1 = make_complex(num1, MAKE_INT(1));   /* special case ...+i */
Erick's avatar
.  
Erick committed
1052
    }
eg's avatar
eg committed
1053 1054
    else if (*p == '-' && p[1] == 'i') {
      p    += 2;
1055
      num1  = make_complex(num1, MAKE_INT(-1)); /* special case ...-i */
eg's avatar
eg committed
1056
    }
1057
    else {                                      /* general case ....[+-@]... */
eg's avatar
eg committed
1058
      polar = (*p == '@') ? (p++,1) : 0;
Erick's avatar
.  
Erick committed
1059

eg's avatar
eg committed
1060 1061
      num2 = read_integer_or_real(p, base, exact, &p);
      if (num2 == STk_false) return STk_false;
Erick's avatar
.  
Erick committed
1062

eg's avatar
eg committed
1063
      if (*p == '/') {
1064 1065 1066
        /* Second member of complex number is a rational */
        num2 = read_rational(num2, p+1, base, exact, &p);
        if (num2 == STk_false) return STk_false;
eg's avatar
eg committed
1067
      }
Erick's avatar
.  
Erick committed
1068

eg's avatar
eg committed
1069
      if (polar) {
1070
        num1 = make_polar(num1, num2);
eg's avatar
eg committed
1071
      } else {
1072 1073 1074 1075
        if (*p == 'i') {
          num1 = make_complex(num1, num2);
          p += 1;
        }
eg's avatar
eg committed
1076 1077 1078 1079 1080 1081 1082
      }
    }
  } else if (*p == 'i' && is_signed) {
    /* We had a number of the form '{+|-}...i' */
    p   += 1;
    num1 = make_complex(MAKE_INT(0), num1);
  }
Erick's avatar
.  
Erick committed
1083

eg's avatar
eg committed
1084 1085 1086 1087 1088 1089
  return (*p) ? STk_false : num1;
}



/******************************************************************************
Erick's avatar
.  
Erick committed
1090
 *
1091
 *                      Scheme primitives and utilities
eg's avatar
eg committed
1092 1093 1094 1095 1096 1097 1098 1099 1100 1101 1102 1103 1104 1105 1106 1107 1108 1109 1110 1111 1112 1113
 *
 ******************************************************************************/


/*
<doc number? complex? real? rational? integer?
 * (number? obj)
 * (complex? obj)
 * (real? obj)
 * (rational? obj)
 * (integer? obj)
 *
 * These numerical type predicates can be applied to any kind of
 * argument, including non-numbers. They return |#t| if the object is of
 * the named type, and otherwise they return |#f|. In general, if a type
 * predicate is true of a number then all higher type predicates are
 * also true of that number. Consequently, if a type predicate is
 * false of a number, then all lower type predicates are also false of
 * that number.
 *
 * If |z| is an inexact complex number, then |(real? z)| is true if and
 * only if |(zero? (imag-part z))| is true. If |x| is an inexact real
Erick's avatar
.  
Erick committed
1114
 * number, then |(integer? x)| is true if and only if
eg's avatar
eg committed
1115
 * |(and (finite? x) (= x (round x)))|
Erick's avatar
.  
Erick committed
1116 1117
 *
 *
eg's avatar
eg committed
1118 1119 1120 1121 1122 1123 1124 1125 1126 1127 1128 1129 1130 1131 1132 1133 1134 1135
 * @lisp
 *   (complex? 3+4i)         =>  #t
 *   (complex? 3)            =>  #t
 *   (real? 3)               =>  #t
 *   (real? -2.5+0.0i)       =>  #t
 *   (real? #e1e10)          =>  #t
 *   (rational? 6/10)        =>  #t
 *   (rational? 6/3)         =>  #t
 *   (integer? 3+0i)         =>  #t
 *   (integer? 3.0)          =>  #t
 *   (integer? 3.2)          =>  #f
 *   (integer? 8/4)          =>  #t
 *   (integer? "no")         =>  #f
 *   (complex? +inf.0)       =>  #t
 *   (real? -inf.0)          =>  #t
 *   (rational? +inf.0)      =>  #f
 *   (integer? -inf.0)       =>  #f
 * @end lisp
Erick's avatar
.  
Erick committed
1136
 *
eg's avatar
eg committed
1137 1138 1139 1140 1141 1142 1143 1144 1145 1146
doc>
 */
DEFINE_PRIMITIVE("number?", numberp, subr1, (SCM x))
{
  switch (TYPEOF (x)) {
    case tc_complex:
    case tc_real:
    case tc_rational:
    case tc_bignum:
    case tc_integer: return STk_true;
1147
    default:         return STk_false;
eg's avatar
eg committed
1148 1149 1150 1151 1152 1153 1154 1155 1156 1157 1158 1159 1160 1161 1162 1163 1164 1165
  }
}


DEFINE_PRIMITIVE("complex?", complexp, subr1, (SCM x))
{
  return STk_numberp(x);
}


DEFINE_PRIMITIVE("real?", realp, subr1, (SCM x))
{
  switch (TYPEOF(x)) {
    case tc_complex: return MAKE_BOOLEAN(zerop(COMPLEX_IMAG(x)));
    case tc_real:
    case tc_rational:
    case tc_bignum:
    case tc_integer: return STk_true;
1166
    default:         return STk_false;
Erick's avatar
.  
Erick committed
1167
  }
eg's avatar
eg committed
1168 1169 1170 1171 1172 1173 1174 1175 1176 1177
}


DEFINE_PRIMITIVE("rational?", rationalp, subr1, (SCM x))
{
  switch (TYPEOF(x)) {
    case tc_real:    return MAKE_BOOLEAN(FINITE_REALP(x));
    case tc_rational:
    case tc_bignum:
    case tc_integer: return STk_true;
1178
    default:         return STk_false;
Erick's avatar
.  
Erick committed
1179
  }
eg's avatar
eg committed
1180 1181 1182 1183 1184 1185
}


/*
<doc EXT bignum?
 * (bignum? x)
Erick's avatar
.  
Erick committed
1186 1187 1188
 *
 * This predicates returns |#t| if |x| is an integer number too large to be
 * represented with a native integer.
eg's avatar
eg committed
1189 1190 1191 1192 1193 1194 1195 1196 1197 1198 1199 1200 1201 1202 1203 1204 1205
 * @lisp
 * (bignum? (expt 2 300))     => |#t|   (very likely)
 * (bignum? 12)               => |#f|
 * (bignum? "no")             => |#f|
 * @end lisp
doc>
a*/
DEFINE_PRIMITIVE("bignum?", bignump, subr1, (SCM x))
{
  return MAKE_BOOLEAN(BIGNUMP(x));
}


DEFINE_PRIMITIVE("integer?", integerp, subr1, (SCM x))
{
  switch (TYPEOF(x)){
    case tc_real:    {
1206 1207 1208 1209 1210
                       double val = REAL_VAL(x);
                       return ((val == minus_inf) || (val == plus_inf)) ?
                                 STk_false:
                                 MAKE_BOOLEAN(floor(val) == val);
                     }
eg's avatar
eg committed
1211 1212
    case tc_bignum:
    case tc_integer: return STk_true;
1213
    default:         return STk_false;
Erick's avatar
.  
Erick committed
1214
  }
eg's avatar
eg committed
1215 1216 1217 1218 1219 1220 1221 1222 1223 1224
}


/*
<doc  exact? inexact?
 * (exact? z)
 * (inexact? z)
 *
 * These numerical predicates provide tests for the exactness of a
 * quantity. For any Scheme number, precisely one of these predicates
Erick's avatar
.  
Erick committed
1225
 * is true.
eg's avatar
eg committed
1226 1227 1228 1229 1230 1231 1232 1233
doc>
 */

static int isexactp(SCM z)
{
  switch (TYPEOF(z)) {
    case tc_complex:  return isexactp(COMPLEX_REAL(z)) && isexactp(COMPLEX_IMAG(z));
    case tc_real:     return FALSE;
Erick's avatar
.  
Erick committed
1234
    case tc_rational:
eg's avatar
eg committed
1235 1236
    case tc_bignum:
    case tc_integer:  return TRUE;
1237
    default:          error_bad_number(z);
eg's avatar
eg committed
1238 1239 1240 1241 1242 1243 1244 1245 1246 1247 1248 1249 1250 1251 1252 1253 1254 1255 1256 1257 1258 1259 1260 1261 1262 1263 1264 1265 1266 1267 1268 1269 1270
  }
  return FALSE; /* never reached */
}


DEFINE_PRIMITIVE("exact?", exactp, subr1, (SCM z))
{
  return MAKE_BOOLEAN(isexactp(z));
}


DEFINE_PRIMITIVE("inexact?", inexactp, subr1, (SCM z))
{
  return MAKE_BOOLEAN(!isexactp(z));
}



/*
<doc  = < > <= >=
 * (= z1 z2 z3 ...)
 * (< x1 x2 x3 ...)
 * (> x1 x2 x3 ...)
 * (<= x1 x2 x3 ...)
 * (>= x1 x2 x3 ...)
 *
 * These procedures return |#t| if their arguments are (respectively):
 * equal, monotonically increasing, monotonically decreasing,
 * monotonically nondecreasing, or monotonically nonincreasing.
 * @lisp
 * (= +inf.0 +inf.0)           =>  #t
 * (= -inf.0 +inf.0)           =>  #f
 * (= -inf.0 -inf.0)           =>  #t
1271
 * @l
eg's avatar
eg committed
1272
 * For any finite real number x:
1273
 * @l
eg's avatar
eg committed
1274 1275 1276 1277 1278
 * (< -inf.0 x +inf.0)         =>  #t
 * (> +inf.0 x -inf.0)         =>  #t
 * @end lisp
doc>
 */
1279 1280 1281 1282 1283 1284 1285 1286 1287 1288 1289 1290 1291
#define COMPARE_NUM(_sname_, _cname_, _max_type_, _operator_)               \
    DEFINE_PRIMITIVE(_sname_, _cname_, vsubr, (int argc, SCM *argv))        \
    {                                                                       \
      SCM previous;                                                         \
                                                                            \
      if (argc == 0) error_at_least_1();                                    \
      if (_max_type_(*argv) == STk_false) error_bad_number(*argv);          \
                                                                            \
      for (previous = *argv--; --argc; previous = *argv--) {                \
        if (_max_type_(*argv) == STk_false) error_bad_number(*argv);        \
        if (do_compare(previous, *argv) _operator_ 0) return STk_false;     \
      }                                                                     \
      return STk_true;                                                      \
eg's avatar
eg committed
1292 1293 1294
    }


1295 1296 1297 1298 1299 1300
#define COMPARE_NUM2(_prim_, _max_type_, _operator_)                        \
    long STk_##_prim_##2(SCM o1, SCM o2)                                    \
    {                                                                       \
      if (_max_type_(o1) == STk_false) error_bad_number(o1);                \
      if (_max_type_(o2) == STk_false) error_bad_number(o2);                \
      return do_compare(o1, o2) _operator_ 0;                               \
eg's avatar
eg committed
1301
    }
Erick's avatar
.  
Erick committed
1302

eg's avatar
eg committed
1303 1304 1305 1306 1307 1308 1309 1310 1311

COMPARE_NUM("=",  numeq, STk_complexp, !=)
COMPARE_NUM("<",  numlt, STk_realp,    >=)
COMPARE_NUM(">",  numgt, STk_realp,    <=)
COMPARE_NUM("<=", numle, STk_realp,    >)
COMPARE_NUM(">=", numge, STk_realp,    <)


/* Version with only two parameters (used by runtime) STk_numeq2, STk_numgt2 ... */
Erick's avatar
.  
Erick committed
1312
COMPARE_NUM2(numeq,   STk_complexp, ==)
eg's avatar
eg committed
1313 1314 1315 1316 1317 1318 1319 1320
COMPARE_NUM2(numlt,   STk_realp,    <)
COMPARE_NUM2(numgt,   STk_realp,    >)
COMPARE_NUM2(numle,   STk_realp,    <=)
COMPARE_NUM2(numge,   STk_realp,    >=)



/*
Erick's avatar
.  
Erick committed
1321
<doc finite? infinite?  zero? positive? negative? odd? even?
eg's avatar
eg committed
1322 1323 1324 1325 1326 1327 1328 1329 1330 1331 1332 1333 1334 1335 1336 1337 1338 1339 1340 1341 1342 1343
 * (finite? z)
 * (infinite? z)
 * (zero? z)
 * (positive? x)
 * (negative? x)
 * (odd? n)
 * (even? n)
 *
 * These numerical predicates test a number for a particular property,
 * returning |#t| or |#f|.
 * @lisp
 * (positive? +inf.0)          ==>  #t
 * (negative? -inf.0)          ==>  #t
 * (finite? -inf.0)            ==>  #f
 * (infinite? +inf.0)          ==>  #t
 * @end lisp
doc>
 */
static Inline int is_oddp(SCM n)
{
  switch (TYPEOF(n)) {
    case tc_integer: return INT_VAL(n) & 1;
1344
    case tc_bignum:  return mpz_odd_p(BIGNUM_VAL(n));
1345
    default:         error_bad_number(n);
eg's avatar
eg committed
1346
  }
1347
  return FALSE; /* never reached */
eg's avatar
eg committed
1348 1349 1350 1351 1352
}

static int zerop(SCM n)
{
  switch (TYPEOF(n)) {
1353
    case tc_integer:  return (INT_VAL(n) == 0);
eg's avatar
eg committed
1354 1355
    case tc_real:     return (REAL_VAL(n) == 0.0);
    case tc_bignum:   return (mpz_cmp_si(BIGNUM_VAL(n), 0L) == 0);
1356 1357
    case tc_complex:  return zerop(COMPLEX_REAL(n)) && zerop(COMPLEX_IMAG(n));
    case tc_rational: return zerop(RATIONAL_NUM(n));
1358
    default:          error_bad_number(n);
eg's avatar
eg committed
1359
  }
1360
  return FALSE; /* never reached */
eg's avatar
eg committed
1361 1362 1363 1364 1365
}

static int positivep(SCM n)
{
  switch (TYPEOF(n)) {
1366
    case tc_integer:  return (INT_VAL(n) > 0);
eg's avatar
eg committed
1367 1368
    case tc_real:     return (REAL_VAL(n) > 0.0);
    case tc_bignum:   return (mpz_cmp_si(BIGNUM_VAL(n), 0L) > 0);
1369
    case tc_rational: return positivep(RATIONAL_NUM(n));
eg's avatar
eg committed
1370
    case tc_complex:  error_not_on_a_complex(n); break;
1371
    default:          error_bad_number(n);
eg's avatar
eg committed
1372
  }
1373
  return FALSE; /* never reached */
eg's avatar
eg committed
1374 1375 1376 1377 1378 1379
}


static int negativep(SCM n)
{
  switch (TYPEOF(n)) {
1380
    case tc_integer:  return (INT_VAL(n) < 0);
eg's avatar
eg committed
1381 1382
    case tc_real:     return (REAL_VAL(n) < 0.0);
    case tc_bignum:   return (mpz_cmp_si(BIGNUM_VAL(n), 0L) < 0);
1383
    case tc_rational: return negativep(RATIONAL_NUM(n));
eg's avatar
eg committed
1384
    case tc_complex:  error_not_on_a_complex(n); break;
1385
    default:          error_bad_number(n);
eg's avatar
eg committed
1386
  }
1387
  return FALSE; /* never reached */
eg's avatar
eg committed
1388 1389 1390 1391 1392 1393 1394
}


static int finitep(SCM n)
{
  switch (TYPEOF(n)) {
    case tc_real:     return (FINITE_REALP(n));
Erick's avatar
.  
Erick committed
1395 1396
    case tc_rational:
    case tc_bignum:
eg's avatar
eg committed
1397
    case tc_integer:  return TRUE;
Erick's avatar
.  
Erick committed
1398
    case tc_complex:  return (finitep(COMPLEX_REAL(n)) &&
1399 1400
                              finitep(COMPLEX_IMAG(n)));
    default:          error_bad_number(n);
eg's avatar
eg committed
1401
  }
1402
  return FALSE; /* never reached */
eg's avatar
eg committed
1403 1404 1405 1406 1407 1408 1409 1410 1411 1412 1413 1414 1415 1416 1417 1418 1419 1420 1421 1422 1423 1424 1425 1426 1427 1428 1429 1430 1431 1432 1433 1434 1435 1436 1437 1438 1439 1440 1441 1442 1443 1444 1445 1446 1447 1448 1449 1450 1451 1452 1453 1454 1455 1456 1457 1458 1459 1460 1461 1462
}


DEFINE_PRIMITIVE("finite?", finitep, subr1, (SCM n))
{
  return MAKE_BOOLEAN(finitep(n));
}


DEFINE_PRIMITIVE("infinite?", infinitep, subr1, (SCM n))
{
  return MAKE_BOOLEAN(!finitep(n));
}


DEFINE_PRIMITIVE("zero?", zerop, subr1, (SCM n))
{
  return MAKE_BOOLEAN(zerop(n));
}


DEFINE_PRIMITIVE("positive?", positivep, subr1, (SCM n))
{
  return MAKE_BOOLEAN(positivep(n));
}


DEFINE_PRIMITIVE("negative?", negativep, subr1, (SCM n))
{
  return MAKE_BOOLEAN(negativep(n));
}


DEFINE_PRIMITIVE("odd?", oddp, subr1, (SCM n))
{
  return MAKE_BOOLEAN(is_oddp(n));
}


DEFINE_PRIMITIVE("even?", evenp, subr1, (SCM n))
{
  return MAKE_BOOLEAN(!is_oddp(n));
}


/*
<doc  max min
 * (max x1 x2 ...)
 * (min x1 x2 ...)
 *
 * These procedures return the maximum or minimum of their arguments.
 *
 * @lisp
 * (max 3 4)              =>  4    ; exact
 * (max 3.9 4)            =>  4.0  ; inexact
 * @end lisp
 * For any real number x:
 * @lisp
 * (max +inf.0 x)         =>  +inf.0
 * (min -inf.0 x)         =>  -inf.0
Erick's avatar
.  
Erick committed
1463
 * @end lisp
1464
 *
eg's avatar
eg committed
1465 1466 1467 1468 1469 1470
 * ,(bold "Note:") If any argument is inexact, then the result will also be
 * inexact
doc>
 */
DEFINE_PRIMITIVE("max", max, vsubr, (int argc, SCM *argv))
{
Erick's avatar
.  
Erick committed
1471
  SCM res;
eg's avatar
eg committed
1472 1473 1474 1475 1476 1477 1478 1479 1480 1481 1482 1483 1484 1485 1486 1487 1488 1489 1490 1491 1492 1493 1494 1495 1496 1497
  int exactp;

  if (argc == 0) error_at_least_1();
  if (argc == 1) {
    if (STk_realp(*argv) == STk_true) return *argv;
    error_bad_number(*argv);
  }

  exactp = isexactp(*argv);

  for (res = *argv--; --argc; argv--) {
    /* See that the argument is a correct number */
    if (STk_realp(*argv) == STk_false) error_bad_number(*argv);

    /* determine if result should be exact or not */
    if (!isexactp(*argv)) exactp = 0;

    /* compute max */
    if (do_compare(res, *argv) < 0) res = *argv;
  }
  return (!exactp && isexactp(res)) ? exact2inexact(res) : res;
}


DEFINE_PRIMITIVE("min", min, vsubr, (int argc, SCM *argv))
{
Erick's avatar
.  
Erick committed
1498
  SCM res;
eg's avatar
eg committed
1499 1500 1501 1502 1503 1504 1505 1506 1507 1508 1509 1510 1511 1512 1513 1514 1515 1516 1517 1518 1519 1520 1521 1522 1523 1524 1525 1526 1527 1528 1529 1530 1531 1532 1533 1534 1535 1536 1537 1538 1539 1540 1541 1542 1543 1544 1545 1546 1547 1548 1549 1550 1551 1552
  int exactp;

  if (argc == 0) error_at_least_1();
  if (argc == 1) {
    if (STk_realp(*argv) == STk_true) return *argv;
    error_bad_number(*argv);
  }

  exactp = isexactp(*argv);

  for (res = *argv--; --argc; argv--) {
    /* See that the argument is a correct number */
    if (STk_realp(*argv) == STk_false) error_bad_number(*argv);

    /* determine if result should be exact or not */
    if (!isexactp(*argv)) exactp = 0;

    /* compute max */
    if (do_compare(res, *argv) > 0) res = *argv;
  }
  return (!exactp && isexactp(res)) ? exact2inexact(res) : res;
}


/*
<doc    + *
 * (+ z1 ...)
 * (* z1 ...)
 *
 * These procedures return the sum or product of their arguments.
 * @lisp
 * (+ 3 4)                 =>  7
 * (+ 3)                   =>  3
 * (+)                     =>  0
 * (+ +inf.0 +inf.0)       =>  +inf.0
 * (+ +inf.0 -inf.0)       =>  +nan.0
 * (* 4)                   =>  4
 * (*)                     =>  1
 * (* 5 +inf.0)            =>  +inf.0
 * (* -5 +inf.0)           =>  -inf.0
 * (* +inf.0 +inf.0)       =>  +inf.0
 * (* +inf.0 -inf.0)       =>  -inf.0
 * (* 0 +inf.0)            =>  +nan.0
 * @end lisp
 * ,(bold "Note:") For any finite number z:
 * @lisp
 * (+ +inf.0 z)            =>  +inf.0
 * (+ -inf.0 z)            =>  -inf.0
 * @end lisp
doc>
 */
SCM STk_add2(SCM o1, SCM o2)
{
  switch (convert(&o1, &o2)) {
Erick's avatar
.  
Erick committed
1553
    case tc_bignum:
eg's avatar
eg committed
1554
        {
1555
          mpz_t add;
Erick's avatar
.  
Erick committed
1556

1557 1558
          mpz_init(add);
          mpz_add(add, BIGNUM_VAL(o1), BIGNUM_VAL(o2));
Erick's avatar
.  
Erick committed
1559

1560 1561 1562 1563
          o1 = bignum2number(add);
          mpz_clear(add);
          break;
        }
Erick's avatar
.  
Erick committed
1564
      case tc_integer:
1565 1566 1567 1568 1569 1570 1571 1572 1573
        {
          long add =  (long) INT_VAL(o1) + INT_VAL(o2);

          if (LONG_FITS_INTEGER(add))
            o1 = MAKE_INT(add);
          else
            o1 = long2scheme_bignum(add);
          break;
        }
eg's avatar
eg committed
1574
      case tc_real:
1575 1576 1577 1578
        {
          o1 = double2real(REAL_VAL(o1) + REAL_VAL(o2));
          break;
        }
eg's avatar
eg committed
1579
      case tc_complex:
1580 1581 1582 1583 1584
        {
          o1 = make_complex(add2(COMPLEX_REAL(o1), COMPLEX_REAL(o2)),
                            add2(COMPLEX_IMAG(o1), COMPLEX_IMAG(o2)));
          break;
        }
eg's avatar
eg committed
1585
      case tc_rational:
1586 1587
        {
          SCM num1, num2, den;
eg's avatar
eg committed
1588

1589 1590 1591
          den  = mul2(RATIONAL_DEN(o1), RATIONAL_DEN(o2));
          num1 = mul2(RATIONAL_NUM(o1), RATIONAL_DEN(o2));
          num2 = mul2(RATIONAL_NUM(o2), RATIONAL_DEN(o1));
eg's avatar
eg committed
1592

1593 1594 1595
          o1 = make_rational(add2(num1, num2), den);
          break;
        }
eg's avatar
eg committed
1596 1597 1598 1599 1600 1601 1602 1603
      default: error_cannot_operate("addition", o1, o2);
  }
  return o1;
}

DEFINE_PRIMITIVE("+", plus, vsubr, (int argc, SCM *argv))
{
  SCM res;
Erick's avatar
.  
Erick committed
1604

eg's avatar
eg committed
1605 1606 1607 1608 1609 1610 1611 1612 1613 1614 1615 1616 1617 1618 1619 1620
  if (argc == 0) return MAKE_INT(0);
  if (argc == 1) return add2(MAKE_INT(0), *argv);

  for (res = *argv--; --argc; argv--)
    res = add2(res, *argv);

  return res;
}


/***
 *** multiplication
 ***/
SCM STk_mul2(SCM o1, SCM o2)
{
  switch (convert(&o1, &o2)) {
Erick's avatar
.  
Erick committed
1621
    case tc_bignum:
1622
      mult_bignum:
eg's avatar
eg committed
1623
      {
1624
        mpz_t prod;
Erick's avatar
.  
Erick committed
1625

1626 1627
        mpz_init(prod);
        mpz_mul(prod, BIGNUM_VAL(o1), BIGNUM_VAL(o2));
Erick's avatar
.  
Erick committed
1628

1629 1630 1631
        o1 = bignum2number(prod);
        mpz_clear(prod);
        break;
eg's avatar
eg committed
1632
      }
Erick's avatar
.  
Erick committed
1633
    case tc_integer:
eg's avatar
eg committed
1634
      {
1635 1636 1637 1638 1639 1640 1641 1642 1643 1644
        long int i1 = INT_VAL(o1);
        long int i2 = INT_VAL(o2);

        o1 = MAKE_INT(i1*i2);
        if (i1 != 0 && (INT_VAL(o1) / i1) != i2) {
          o1 = long2scheme_bignum(i1);
          o2 = long2scheme_bignum(i2);
          goto mult_bignum;
        }
        break;
eg's avatar
eg committed
1645 1646
      }
      case tc_real:
1647 1648 1649 1650
        {
          o1 = double2real(REAL_VAL(o1) * REAL_VAL(o2));
          break;
        }
eg's avatar
eg committed
1651
      case tc_complex:
1652 1653 1654 1655 1656 1657 1658 1659 1660 1661
        {
          SCM r1 = COMPLEX_REAL(o1);
          SCM i1 = COMPLEX_IMAG(o1);
          SCM r2 = COMPLEX_REAL(o2);
          SCM i2 = COMPLEX_IMAG(o2);

          o1 = make_complex(sub2(mul2(r1,r2), mul2(i1, i2)),
                            add2(mul2(r1,i2), mul2(r2, i1)));
          break;
        }
eg's avatar
eg committed
1662
      case tc_rational:
1663 1664 1665 1666 1667
        {
          o1 = make_rational(mul2(RATIONAL_NUM(o1), RATIONAL_NUM(o2)),
                             mul2(RATIONAL_DEN(o1), RATIONAL_DEN(o2)));
          break;
        }
eg's avatar
eg committed
1668 1669 1670 1671 1672 1673 1674 1675 1676 1677 1678 1679 1680 1681
      default: error_cannot_operate("multiplication", o1, o2);
  }
  return o1;
}

DEFINE_PRIMITIVE("*", multiplication, vsubr, (int argc, SCM *argv))
{
  SCM res;

  if (argc == 0) return MAKE_INT(1);
  if (argc == 1) return mul2(MAKE_INT(1), *argv);

  for (res = *argv--; --argc; argv--)
    res = mul2(res, *argv);
Erick's avatar
.  
Erick committed
1682

eg's avatar
eg committed
1683 1684 1685 1686 1687 1688 1689 1690 1691 1692 1693 1694 1695 1696 1697 1698 1699 1700 1701 1702 1703 1704 1705 1706 1707 1708 1709 1710 1711
  return res;
}

/*
<doc   - /
 * (- z)
 * (- z1 z2)
 * (/ z)
 * (/ z1 z2 ...)
 *
 * With two or more arguments, these procedures return the difference or quotient
 * of their arguments, associating to the left. With one argument, however,
 * they return the additive or multiplicative inverse of their argument.
 *
 * @lisp
 * (- 3 4)                 =>  -1
 * (- 3 4 5)               =>  -6
 * (- 3)                   =>  -3
 * (- +inf.0 +inf.0)       => +nan.0
 * (/ 3 4 5)               =>  3/20
 * (/ 3)                   =>  1/3
 * (/ 0.0)                 => +inf.0
 * (/ 0)                   => error (division by 0)
 * @end lisp
doc>
 */
SCM STk_sub2(SCM o1, SCM o2)
{
  switch (convert(&o1, &o2)) {
Erick's avatar
.  
Erick committed
1712
    case tc_bignum:
eg's avatar
eg committed
1713
      {
1714
        mpz_t sub;
Erick's avatar
.  
Erick committed
1715

1716 1717
        mpz_init(sub);
        mpz_sub(sub, BIGNUM_VAL(o1), BIGNUM_VAL(o2));
eg's avatar
eg committed
1718

1719 1720 1721
        o1 = bignum2number(sub),
        mpz_clear(sub);
        break;
eg's avatar
eg committed
1722
      }
Erick's avatar
.  
Erick committed
1723
    case tc_integer:
eg's avatar
eg committed
1724
      {
1725 1726 1727 1728 1729 1730
        long sub = (long) INT_VAL(o1) - INT_VAL(o2);
        if (LONG_FITS_INTEGER(sub))
          o1 = MAKE_INT(sub);
        else
          o1 = long2scheme_bignum(sub);
        break;
eg's avatar
eg committed
1731 1732
      }
      case tc_real:
1733 1734 1735 1736
        {
          o1 = double2real(REAL_VAL(o1) - REAL_VAL(o2));
          break;
        }
eg's avatar
eg committed
1737
      case tc_complex:
1738 1739 1740 1741 1742
        {
          o1 = make_complex(sub2(COMPLEX_REAL(o1), COMPLEX_REAL(o2)),
                            sub2(COMPLEX_IMAG(o1), COMPLEX_IMAG(o2)));
          break;
        }
eg's avatar
eg committed
1743
      case tc_rational:
1744 1745
        {
          SCM num1, num2, den;
eg's avatar
eg committed
1746

1747 1748 1749
          den  = mul2(RATIONAL_DEN(o1), RATIONAL_DEN(o2));
          num1 = mul2(RATIONAL_NUM(o1), RATIONAL_DEN(o2));
          num2 = mul2(RATIONAL_NUM(o2), RATIONAL_DEN(o1));
eg's avatar
eg committed
1750

1751 1752 1753
          o1 = make_rational(sub2(num1, num2), den);
          break;
        }
eg's avatar
eg committed
1754 1755 1756 1757 1758 1759 1760 1761 1762 1763 1764 1765 1766 1767 1768 1769
      default: error_cannot_operate("substraction", o1, o2);
  }
  return o1;
}


DEFINE_PRIMITIVE("-", difference, vsubr, (int argc, SCM *argv))
{
  SCM res;

  if (argc == 0) error_at_least_1();
  if (argc == 1) return sub2(MAKE_INT(0), *argv);

  for (res = *argv-- ; --argc; argv--)
    res = sub2(res, *argv);
  return res;
Erick's avatar
.  
Erick committed
1770
}
eg's avatar
eg committed
1771 1772 1773 1774 1775 1776 1777 1778


/***
 ***   Division
 ***/
SCM STk_div2(SCM o1, SCM o2)
{
  switch (convert(&o1, &o2)) {
Erick's avatar
.  
Erick committed
1779
    case tc_bignum:
eg's avatar
eg committed
1780 1781 1782 1783 1784
    case tc_integer:
      o1 = make_rational(o1, o2);
      break;
    case tc_real:
      {
1785
        double r2 = REAL_VAL(o2);
Erick's avatar
.  
Erick committed
1786

1787 1788 1789
        if (r2 != 1.0)
          o1 = double2real(REAL_VAL(o1) / r2);
        break;
eg's avatar
eg committed
1790 1791 1792
      }
    case tc_rational:
      o1 =  make_rational(mul2(RATIONAL_NUM(o1), RATIONAL_DEN(o2)),
1793
                          mul2(RATIONAL_DEN(o1), RATIONAL_NUM(o2)));
eg's avatar
eg committed
1794
      break;
Erick's avatar
.  
Erick committed
1795
    case tc_complex:
eg's avatar
eg committed
1796
      {
1797 1798 1799 1800 1801 1802 1803 1804 1805 1806 1807 1808 1809 1810
        SCM tmp, new_r, new_i;

        if (!zerop(o1)) {
          tmp   = add2(mul2(COMPLEX_REAL(o2), COMPLEX_REAL(o2)),
                       mul2(COMPLEX_IMAG(o2), COMPLEX_IMAG(o2)));
          new_r = div2(add2(mul2(COMPLEX_REAL(o1), COMPLEX_REAL(o2)),
                            mul2(COMPLEX_IMAG(o1), COMPLEX_IMAG(o2))),
                       tmp);
          new_i = div2(sub2(mul2(COMPLEX_IMAG(o1), COMPLEX_REAL(o2)),
                            mul2(COMPLEX_REAL(o1), COMPLEX_IMAG(o2))),
                       tmp);
          o1 = make_complex(new_r, new_i);
        }
        break;
eg's avatar
eg committed
1811 1812 1813 1814 1815 1816 1817 1818 1819 1820 1821 1822 1823 1824 1825 1826 1827 1828 1829 1830 1831 1832 1833 1834 1835 1836 1837 1838 1839 1840 1841 1842 1843 1844 1845
      }
    default: error_cannot_operate("division", o1, o2);
  }
  return o1;
}

DEFINE_PRIMITIVE("/", division, vsubr, (int argc, SCM *argv))
{
  SCM res;

  if (argc == 0) error_at_least_1();
  if (argc == 1) return div2(MAKE_INT(1), *argv);

  for (res = *argv--; --argc; argv--)
    res = div2(res, *argv);
  return res;
}


/*
<doc  abs
 * (abs x)
 *
 * |Abs| returns the absolute value of its argument.
 * @lisp
 * (abs -7)                =>  7
 * (abs -inf.0)            => +inf.0
 * @end lisp
doc>
 */
DEFINE_PRIMITIVE("abs", abs, subr1, (SCM x))
{
  switch (TYPEOF(x)) {
    case tc_integer:  return (INT_VAL(x) < 0) ? MAKE_INT(-INT_VAL(x)) : x;
    case tc_bignum:   if (mpz_cmp_ui(BIGNUM_VAL(x), 0L) < 0) {
1846 1847 1848 1849 1850 1851 1852 1853
                        mpz_t tmp;

                        mpz_init(tmp);
                        mpz_neg(tmp, BIGNUM_VAL(x));
                        x = bignum2scheme_bignum(tmp);
                        mpz_clear(tmp);
                      }
                      return x;
eg's avatar
eg committed
1854 1855
    case tc_real:     return