cgit.sxemacs.org Git - sxemacs/blob - src/ent/floatfns.c

   1 /* Primitive operations on floating point for SXEmacs Lisp interpreter.
   2    Copyright (C) 1988, 1993, 1994 Free Software Foundation, Inc.
   3
   4 This file is part of SXEmacs
   5
   6 SXEmacs is free software: you can redistribute it and/or modify
   7 it under the terms of the GNU General Public License as published by
   8 the Free Software Foundation, either version 3 of the License, or
   9 (at your option) any later version.
  10
  11 SXEmacs is distributed in the hope that it will be useful,
  12 but WITHOUT ANY WARRANTY; without even the implied warranty of
  13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  14 GNU General Public License for more details.
  15
  16 You should have received a copy of the GNU General Public License
  17 along with this program.  If not, see <http://www.gnu.org/licenses/>. */
  18
  19
  20 /* Synched up with: FSF 19.30. */
  21
  22 /* ANSI C requires only these float functions:
  23    acos, asin, atan, atan2, ceil, cos, cosh, exp, fabs, floor, fmod,
  24    frexp, ldexp, log, log10, modf, pow, sin, sinh, sqrt, tan, tanh.
  25
  26    Define HAVE_INVERSE_HYPERBOLIC if you have acosh, asinh, and atanh.
  27    Define HAVE_CBRT if you have cbrt().
  28    Define HAVE_RINT if you have rint().
  29    If you don't define these, then the appropriate routines will be simulated.
  30
  31    Define HAVE_MATHERR if on a system supporting the SysV matherr() callback.
  32    (This should happen automatically.)
  33
  34    Define FLOAT_CHECK_ERRNO if the float library routines set errno.
  35    This has no effect if HAVE_MATHERR is defined.
  36
  37    Define FLOAT_CATCH_SIGILL if the float library routines signal SIGILL.
  38    (What systems actually do this?  Let me know. -jwz)
  39
  40    Define FLOAT_CHECK_DOMAIN if the float library doesn't handle errors by
  41    either setting errno, or signalling SIGFPE/SIGILL.  Otherwise, domain and
  42    range checking will happen before calling the float routines.  This has
  43    no effect if HAVE_MATHERR is defined (since matherr will be called when
  44    a domain error occurs).
  45  */
  46
  47 #include <config.h>
  48 #include "lisp.h"
  49 #include "syssignal.h"
  50
  51 #include "ent.h"
  52
  53 #ifdef HAVE_FPFLOAT
  54
  55 /* The code uses emacs_rint, so that it works to undefine HAVE_RINT
  56    if `rint' exists but does not work right.  */
  57 #ifdef HAVE_RINT
  58 #define emacs_rint rint
  59 #else
  60 static fpfloat emacs_rint(fpfloat x)
  61 {
  62         fpfloat r = floor(x + 0.5);
  63         fpfloat diff = fabs(r - x);
  64         /* Round to even and correct for any roundoff errors.  */
  65         if (diff >= 0.5 && (diff > 0.5 || r != 2.0 * floor(r / 2.0)))
  66                 r += r < x ? 1.0 : -1.0;
  67         return r;
  68 }
  69 #endif
  70
  71 /* Nonzero while executing in floating point.
  72    This tells float_error what to do.  */
  73 static int in_float;
  74
  75 /* If an argument is out of range for a mathematical function,
  76    here is the actual argument value to use in the error message.  */
  77 static Lisp_Object float_error_arg, float_error_arg2;
  78 static const char *float_error_fn_name;
  79
  80 /* Evaluate the floating point expression D, recording NUM
  81    as the original argument for error messages.
  82    D is normally an assignment expression.
  83    Handle errors which may result in signals or may set errno.
  84
  85    Note that float_error may be declared to return void, so you can't
  86    just cast the zero after the colon to (SIGTYPE) to make the types
  87    check properly.  */
  88 #ifdef FLOAT_CHECK_ERRNO
  89 #define IN_FLOAT(d, name, num)                          \
  90   do {                                                  \
  91     float_error_arg = num;                              \
  92     float_error_fn_name = name;                         \
  93     in_float = 1; errno = 0; (d); in_float = 0;         \
  94     if (errno != 0) in_float_error ();                  \
  95   } while (0)
  96 #define IN_FLOAT2(d, name, num, num2)                   \
  97   do {                                                  \
  98     float_error_arg = num;                              \
  99     float_error_arg2 = num2;                            \
 100     float_error_fn_name = name;                         \
 101     in_float = 2; errno = 0; (d); in_float = 0;         \
 102     if (errno != 0) in_float_error ();                  \
 103   } while (0)
 104 #else
 105 #define IN_FLOAT(d, name, num) (in_float = 1, (d), in_float = 0)
 106 #define IN_FLOAT2(d, name, num, num2) (in_float = 2, (d), in_float = 0)
 107 #endif
 108
 109 #define arith_error(op,arg) \
 110   Fsignal (Qarith_error, list2 (build_string (op), arg))
 111 #define range_error(op,arg) \
 112   Fsignal (Qrange_error, list2 (build_string (op), arg))
 113 #define range_error2(op,a1,a2) \
 114   Fsignal (Qrange_error, list3 (build_string (op), a1, a2))
 115 #define domain_error(op,arg) \
 116   Fsignal (Qdomain_error, list2 (build_string (op), arg))
 117 #define domain_error2(op,a1,a2) \
 118   Fsignal (Qdomain_error, list3 (build_string (op), a1, a2))
 119
 120 /* Convert float to Lisp Integer if it fits, else signal a range
 121    error using the given arguments.
 122    If numbers from multi-prec libraries are available, range errors
 123    are never signaled.
 124 */
 125 static Lisp_Object
 126 float_to_int(fpfloat x, const char *name, Lisp_Object num, Lisp_Object num2)
 127 {
 128 #if defined HAVE_MPZ && (defined WITH_GMP || defined WITH_MP)
 129         bigz_set_fpfloat(ent_scratch_bigz, x);
 130         return ent_mpz_downgrade_maybe(ent_scratch_bigz);
 131
 132         if (name || NILP(num) || NILP(num2));
 133         return Qnil;
 134 #else  /* !HAVE_MPZ */
 135         REGISTER EMACS_INT result = (EMACS_INT) x;
 136
 137         if (result > EMACS_INT_MAX || result < EMACS_INT_MIN) {
 138                 if (!UNBOUNDP(num2))
 139                         range_error2(name, num, num2);
 140                 else
 141                         range_error(name, num);
 142         }
 143         return make_int(result);
 144 #endif  /* HAVE_MPZ */
 145 }
 146
 147 static void in_float_error(void)
 148 {
 149         switch (errno) {
 150         case 0:
 151                 break;
 152         case EDOM:
 153                 if (in_float == 2)
 154                         domain_error2(float_error_fn_name, float_error_arg,
 155                                       float_error_arg2);
 156                 else
 157                         domain_error(float_error_fn_name, float_error_arg);
 158                 break;
 159         case ERANGE:
 160                 range_error(float_error_fn_name, float_error_arg);
 161                 break;
 162         default:
 163                 arith_error(float_error_fn_name, float_error_arg);
 164                 break;
 165         }
 166 }
 167
 168 \f
 169 #endif  /* HAVE_FPFLOAT */
 170 \f
 171 /* Trig functions.  */
 172
 173 #if defined HAVE_MPFR && defined WITH_MPFR
 174 #define MPFR_TRIG_FUN(op) do                                            \
 175 {                                                                       \
 176         Lisp_Object bfrnumber;                                          \
 177                                                                         \
 178         if (INDEFP(number))                                             \
 179                 return make_indef(NOT_A_NUMBER);                        \
 180                                                                         \
 181         bigfr_set_prec(ent_scratch_bigfr,                               \
 182                        internal_get_precision(precision));              \
 183                                                                         \
 184         bfrnumber = Fcoerce_number(number, Qbigfr, Qnil);               \
 185         bigfr_##op(ent_scratch_bigfr, XBIGFR_DATA(bfrnumber));          \
 186         return make_bigfr_bfr(ent_scratch_bigfr);                       \
 187 } while (0)
 188 #endif
 189
 190 #if defined(HAVE_MPFR) && defined WITH_MPFR || defined(HAVE_FPFLOAT)
 191
 192 DEFUN("acos", Facos, 1, 2, 0,   /*
 193 Return the inverse cosine of NUMBER.
 194 If optional argument PRECISION is non-nil, its value
 195 (an integer) is used as precision.
 196 */
 197       (number, precision))
 198 {
 199 #if defined HAVE_MPFR && defined WITH_MPFR
 200
 201         MPFR_TRIG_FUN(acos);
 202
 203 #else  /* !HAVE_MPFR */
 204         if (INDEFP(number)) {
 205                 return make_indef(NOT_A_NUMBER);
 206         }
 207
 208         number = ent_lift(number, FLOAT_T, NULL);
 209
 210         if (FLOATP(number)) {
 211                 fpfloat d;
 212                 d = acos(XFLOAT_DATA(number));
 213                 return make_float(d);
 214         } else if (INDEFP(number)) {
 215                 return make_indef(NOT_A_NUMBER);
 216         }
 217
 218         Fsignal(Qarith_error, list1(number));
 219         return Qnil;
 220 #endif  /* HAVE_MPFR */
 221 }
 222
 223 DEFUN("asin", Fasin, 1, 2, 0,   /*
 224 Return the inverse sine of NUMBER.
 225 If optional argument PRECISION is non-nil, its value
 226 (an integer) is used as precision.
 227 */
 228       (number, precision))
 229 {
 230 #if defined HAVE_MPFR && defined WITH_MPFR
 231
 232         MPFR_TRIG_FUN(asin);
 233
 234 #else  /* !HAVE_MPFR */
 235         if (INDEFP(number)) {
 236                 return make_indef(NOT_A_NUMBER);
 237         }
 238
 239         number = ent_lift(number, FLOAT_T, NULL);
 240
 241         if (FLOATP(number)) {
 242                 fpfloat d;
 243                 d = asin(XFLOAT_DATA(number));
 244                 return make_float(d);
 245         } else if (INDEFP(number)) {
 246                 return make_indef(NOT_A_NUMBER);
 247         }
 248
 249         Fsignal(Qarith_error, list1(number));
 250         return Qnil;
 251
 252         if (NILP(precision));
 253 #endif  /* HAVE_MPFR */
 254 }
 255
 256 DEFUN("atan", Fatan, 1, 3, 0,   /*
 257 Return the inverse tangent of NUMBER.
 258 If optional second argument NUMBER2 is provided,
 259 return atan2 (NUMBER, NUMBER2).
 260 If optional argument PRECISION is non-nil, its value
 261 (an integer) is used as precision.
 262 */
 263       (number, number2, precision))
 264 {
 265 #if defined HAVE_MPFR && defined WITH_MPFR
 266         Lisp_Object result;
 267
 268         if (NILP(number2)) {
 269                 Lisp_Object bfrnumber;
 270
 271                 if (INDEFP(number))
 272                         return make_indef(NOT_A_NUMBER);
 273
 274                 bigfr_set_prec(ent_scratch_bigfr,
 275                                internal_get_precision(precision));
 276                 bfrnumber = Fcoerce_number(number, Qbigfr, Qnil);
 277                 bigfr_atan(ent_scratch_bigfr, XBIGFR_DATA(bfrnumber));
 278                 result = make_bigfr_bfr(ent_scratch_bigfr);
 279         } else {
 280                 Lisp_Object bfrn1;
 281                 Lisp_Object bfrn2;
 282
 283                 if (INDEFP(number))
 284                         return make_indef(NOT_A_NUMBER);
 285                 if (INFINITYP(number2))
 286                         return Qzero;
 287                 if (INDEFP(number2))
 288                         return make_indef(NOT_A_NUMBER);
 289
 290                 bigfr_set_prec(ent_scratch_bigfr,
 291                                internal_get_precision(precision));
 292                 bfrn1 = Fcoerce_number(number, Qbigfr, Qnil);
 293                 bfrn2 = Fcoerce_number(number2, Qbigfr, Qnil);
 294                 bigfr_atan2(ent_scratch_bigfr,
 295                             XBIGFR_DATA(bfrn1),
 296                             XBIGFR_DATA(bfrn2));
 297                 result = make_bigfr_bfr(ent_scratch_bigfr);
 298         }
 299
 300         return result;
 301
 302 #else  /* !HAVE_MPFR */
 303         if (INDEFP(number)) {
 304                 return make_indef(NOT_A_NUMBER);
 305         }
 306
 307         number = ent_lift(number, FLOAT_T, NULL);
 308
 309         if (INDEFP(number))
 310                 return make_indef(NOT_A_NUMBER);
 311
 312         if (NILP(number2) && FLOATP(number)) {
 313                 fpfloat d;
 314                 d = atan(XFLOAT_DATA(number));
 315                 return make_float(d);
 316         } else if (FLOATP(number)) {
 317                 number = ent_lift(number2, FLOAT_T, NULL);
 318
 319                 if (FLOATP(number2)) {
 320                         fpfloat d;
 321                         d = atan2(XFLOAT_DATA(number), XFLOAT_DATA(number2));
 322                         return make_float(d);
 323                 } else if (INFINITYP(number2)) {
 324                         return Qzero;
 325                 } else if (INDEFP(number2)) {
 326                         return make_indef(NOT_A_NUMBER);
 327                 }
 328         }
 329
 330         /* Just signal here, I'm not in the mood to distinguish cases here */
 331         Fsignal(Qarith_error, list1(number));
 332         return Qnil;
 333
 334         if (NILP(precision));
 335 #endif  /* HAVE_MPFR */
 336 }
 337
 338 DEFUN("cos", Fcos, 1, 2, 0,     /*
 339 Return the cosine of NUMBER.
 340 If optional argument PRECISION is non-nil, its value
 341 (an integer) is used as precision.
 342 */
 343       (number, precision))
 344 {
 345 #if defined HAVE_MPFR && defined WITH_MPFR
 346
 347         MPFR_TRIG_FUN(cos);
 348
 349 #else  /* !HAVE_MPFR */
 350         if (INDEFP(number)) {
 351                 return make_indef(NOT_A_NUMBER);
 352         }
 353
 354         number = ent_lift(number, FLOAT_T, NULL);
 355
 356         if (FLOATP(number)) {
 357                 fpfloat d;
 358                 d = cos(XFLOAT_DATA(number));
 359                 return make_float(d);
 360         } else if (INDEFP(number)) {
 361                 return make_indef(NOT_A_NUMBER);
 362         }
 363
 364         Fsignal(Qarith_error, list1(number));
 365         return Qnil;
 366
 367         if (NILP(precision));
 368 #endif  /* HAVE_MPFR */
 369 }
 370
 371 DEFUN("sin", Fsin, 1, 2, 0,     /*
 372 Return the sine of NUMBER.
 373 If optional argument PRECISION is non-nil, its value
 374 (an integer) is used as precision.
 375 */
 376       (number, precision))
 377 {
 378 #if defined HAVE_MPFR && defined WITH_MPFR
 379
 380         MPFR_TRIG_FUN(sin);
 381
 382 #else  /* !HAVE_MPFR */
 383         if (INDEFP(number)) {
 384                 return make_indef(NOT_A_NUMBER);
 385         }
 386
 387         number = ent_lift(number, FLOAT_T, NULL);
 388
 389         if (FLOATP(number)) {
 390                 fpfloat d;
 391                 d = sin(XFLOAT_DATA(number));
 392                 return make_float(d);
 393         } else if (INDEFP(number)) {
 394                 return make_indef(NOT_A_NUMBER);
 395         }
 396
 397         Fsignal(Qarith_error, list1(number));
 398         return Qnil;
 399
 400         if (NILP(precision));
 401 #endif  /* HAVE_MPFR */
 402 }
 403
 404 DEFUN("tan", Ftan, 1, 2, 0,     /*
 405 Return the tangent of NUMBER.
 406 If optional argument PRECISION is non-nil, its value
 407 (an integer) is used as precision.
 408 */
 409       (number, precision))
 410 {
 411 #if defined HAVE_MPFR && defined WITH_MPFR
 412
 413         MPFR_TRIG_FUN(tan);
 414
 415 #else  /* !HAVE_MPFR */
 416         if (INDEFP(number)) {
 417                 return make_indef(NOT_A_NUMBER);
 418         }
 419
 420         number = ent_lift(number, FLOAT_T, NULL);
 421
 422         if (FLOATP(number)) {
 423                 fpfloat d;
 424                 d = XFLOAT_DATA(number);
 425                 d = sin(d) / cos(d);
 426                 return make_float(d);
 427         } else if (INDEFP(number)) {
 428                 return make_indef(NOT_A_NUMBER);
 429         }
 430
 431         Fsignal(Qarith_error, list1(number));
 432         return Qnil;
 433
 434         if (NILP(precision));
 435 #endif  /* HAVE_MPFR */
 436 }
 437
 438 #if defined HAVE_MPFR && defined WITH_MPFR
 439 DEFUN("sec", Fsec, 1, 2, 0,     /*
 440 Return the secant of NUMBER.
 441 If optional argument PRECISION is non-nil, its value
 442 (an integer) is used as precision.
 443 */
 444       (number, precision))
 445 {
 446         MPFR_TRIG_FUN(sec);
 447 }
 448
 449 DEFUN("csc", Fcsc, 1, 2, 0,     /*
 450 Return the cosecant of NUMBER.
 451 If optional argument PRECISION is non-nil, its value
 452 (an integer) is used as precision.
 453 */
 454       (number, precision))
 455 {
 456         MPFR_TRIG_FUN(csc);
 457 }
 458
 459 DEFUN("cot", Fcot, 1, 2, 0,     /*
 460 Return the cotangent of NUMBER.
 461 If optional argument PRECISION is non-nil, its value
 462 (an integer) is used as precision.
 463 */
 464       (number, precision))
 465 {
 466         MPFR_TRIG_FUN(cot);
 467 }
 468 #endif  /* HAVE_MPFR */
 469
 470 #endif  /* HAVE_MPFR || HAVE_FPFLOAT (trig functions) */
 471 \f
 472 /* Bessel functions */
 473 #if 0                           /* Leave these out unless we find there's a reason for them.  */
 474 /* #ifdef HAVE_FPFLOAT */
 475
 476 DEFUN("bessel-j0", Fbessel_j0, 1, 1, 0, /*
 477 Return the bessel function j0 of NUMBER.
 478 */
 479       (number))
 480 {
 481         fpfloat d = extract_float(number);
 482         IN_FLOAT(d = j0(d), "bessel-j0", number);
 483         return make_float(d);
 484 }
 485
 486 DEFUN("bessel-j1", Fbessel_j1, 1, 1, 0, /*
 487 Return the bessel function j1 of NUMBER.
 488 */
 489       (number))
 490 {
 491         fpfloat d = extract_float(number);
 492         IN_FLOAT(d = j1(d), "bessel-j1", number);
 493         return make_float(d);
 494 }
 495
 496 DEFUN("bessel-jn", Fbessel_jn, 2, 2, 0, /*
 497 Return the order N bessel function output jn of NUMBER.
 498 The first number (the order) is truncated to an integer.
 499 */
 500       (number1, number2))
 501 {
 502         int i1 = extract_float(number1);
 503         fpfloat f2 = extract_float(number2);
 504
 505         IN_FLOAT(f2 = jn(i1, f2), "bessel-jn", number1);
 506         return make_float(f2);
 507 }
 508
 509 DEFUN("bessel-y0", Fbessel_y0, 1, 1, 0, /*
 510 Return the bessel function y0 of NUMBER.
 511 */
 512       (number))
 513 {
 514         fpfloat d = extract_float(number);
 515         IN_FLOAT(d = y0(d), "bessel-y0", number);
 516         return make_float(d);
 517 }
 518
 519 DEFUN("bessel-y1", Fbessel_y1, 1, 1, 0, /*
 520 Return the bessel function y1 of NUMBER.
 521 */
 522       (number))
 523 {
 524         fpfloat d = extract_float(number);
 525         IN_FLOAT(d = y1(d), "bessel-y0", number);
 526         return make_float(d);
 527 }
 528
 529 DEFUN("bessel-yn", Fbessel_yn, 2, 2, 0, /*
 530 Return the order N bessel function output yn of NUMBER.
 531 The first number (the order) is truncated to an integer.
 532 */
 533       (number1, number2))
 534 {
 535         int i1 = extract_float(number1);
 536         fpfloat f2 = extract_float(number2);
 537
 538         IN_FLOAT(f2 = yn(i1, f2), "bessel-yn", number1);
 539         return make_float(f2);
 540 }
 541
 542 #endif                          /* 0 (bessel functions) */
 543
 544 \f
 545 /* Error functions. */
 546 #if defined(HAVE_MPFR) && defined WITH_MPFR || defined(HAVE_FPFLOAT)
 547 DEFUN("erf", Ferf, 1, 2, 0,     /*
 548 Return the mathematical error function of NUMBER.
 549 */
 550       (number, precision))
 551 {
 552 #if defined HAVE_MPFR && defined WITH_MPFR
 553
 554         MPFR_TRIG_FUN(erf);
 555
 556 #else  /* !HAVE_MPFR */
 557         if (INDEFP(number)) {
 558                 return make_indef(NOT_A_NUMBER);
 559         }
 560
 561         number = ent_lift(number, FLOAT_T, NULL);
 562
 563         if (FLOATP(number)) {
 564                 fpfloat d;
 565                 d = erf(XFLOAT_DATA(number));
 566                 return make_float(d);
 567         } else if (INDEFP(number)) {
 568                 return make_indef(NOT_A_NUMBER);
 569         }
 570
 571         Fsignal(Qarith_error, list1(number));
 572         return Qnil;
 573
 574         if (NILP(precision));
 575 #endif  /* HAVE_MPFR */
 576 }
 577
 578 DEFUN("erfc", Ferfc, 1, 2, 0,   /*
 579 Return the complementary error function of NUMBER.
 580 */
 581       (number, precision))
 582 {
 583 #if defined HAVE_MPFR && defined WITH_MPFR
 584
 585         MPFR_TRIG_FUN(erfc);
 586
 587 #else  /* !HAVE_MPFR */
 588         if (INDEFP(number)) {
 589                 return make_indef(NOT_A_NUMBER);
 590         }
 591
 592         number = ent_lift(number, FLOAT_T, NULL);
 593
 594         if (FLOATP(number)) {
 595                 fpfloat d;
 596                 d = erfc(XFLOAT_DATA(number));
 597                 return make_float(d);
 598         } else if (INDEFP(number)) {
 599                 return make_indef(NOT_A_NUMBER);
 600         }
 601
 602         Fsignal(Qarith_error, list1(number));
 603         return Qnil;
 604
 605         if (NILP(precision));
 606 #endif  /* HAVE_MPFR */
 607 }
 608
 609 DEFUN("log-gamma", Flog_gamma, 1, 2, 0, /*
 610 Return the log gamma of NUMBER.
 611 */
 612       (number, precision))
 613 {
 614 #if defined HAVE_MPFR && defined WITH_MPFR
 615
 616         MPFR_TRIG_FUN(lgamma);
 617
 618 #else  /* !HAVE_MPFR */
 619         if (INDEFP(number)) {
 620                 return make_indef(NOT_A_NUMBER);
 621         }
 622
 623         number = ent_lift(number, FLOAT_T, NULL);
 624
 625         if (FLOATP(number)) {
 626                 fpfloat d;
 627                 d = lgamma(XFLOAT_DATA(number));
 628                 return make_float(d);
 629         } else if (INDEFP(number)) {
 630                 return make_indef(NOT_A_NUMBER);
 631         }
 632
 633         Fsignal(Qarith_error, list1(number));
 634         return Qnil;
 635
 636         if (NILP(precision));
 637 #endif  /* HAVE_MPFR */
 638 }
 639 #endif  /* HAVE_FPFLOAT || HAVE_MPFR */
 640
 641 \f
 642 /* Root and Log functions. */
 643
 644 #if defined(HAVE_FPFLOAT) || defined(HAVE_MPFR) && defined WITH_MPFR
 645 DEFUN("exp", Fexp, 1, 2, 0,     /*
 646 Return the exponential base e of NUMBER.
 647 If optional argument PRECISION is non-nil, its value
 648 (an integer) is used as precision.
 649 */
 650       (number, precision))
 651 {
 652 /* Attention, somehow the precision must be large enough to make the result
 653  * fit, otherwise this is a good memory test :)
 654  */
 655 #if defined(HAVE_MPFR) && defined WITH_MPFR ||  \
 656         defined(HAVE_MPC) && defined WITH_MPC ||        \
 657         defined HAVE_PSEUC && defined WITH_PSEUC
 658
 659         if (INDEFP(number)) {
 660                 if (XINDEF_DATA(number) == POS_INFINITY)
 661                         return number;
 662                 else if (XINDEF_DATA(number) == NEG_INFINITY)
 663                         return Fcoerce_number(Qzero, Qbigfr, precision);
 664                 else
 665                         return number;
 666         }
 667
 668         if (COMPARABLEP(number)) {
 669 #if defined HAVE_MPFR && defined WITH_MPFR
 670                 Lisp_Object bfrnumber;
 671
 672                 bigfr_set_prec(ent_scratch_bigfr,
 673                                internal_get_precision(precision));
 674
 675                 bfrnumber = Fcoerce_number(number, Qbigfr, Qnil);
 676                 bigfr_exp(ent_scratch_bigfr, XBIGFR_DATA(bfrnumber));
 677                 return make_bigfr_bfr(ent_scratch_bigfr);
 678
 679 #endif  /* HAVE_MPFR */
 680 #if defined HAVE_MPC && defined WITH_MPC ||     \
 681         defined HAVE_PSEUC && defined WITH_PSEUC
 682         } else if (BIGCP(number)) {
 683                 bigc_set_prec(ent_scratch_bigc,
 684                               internal_get_precision(precision));
 685
 686                 bigc_exp(ent_scratch_bigc, XBIGC_DATA(number));
 687                 return make_bigc_bc(ent_scratch_bigc);
 688 #endif  /* HAVE_MPC */
 689         }
 690
 691         return wrong_type_argument(Qnumberp, number);
 692 #else  /* !HAVE_MPFR && !HAVE_MPC */
 693         if (INDEFP(number)) {
 694                 goto indefcase;
 695         }
 696
 697         number = ent_lift(number, FLOAT_T, NULL);
 698
 699         if (FLOATP(number)) {
 700                 fpfloat d;
 701                 d = exp(XFLOAT_DATA(number));
 702                 return make_float(d);
 703         } else if (INDEFP(number)) {
 704         indefcase:
 705                 if (XINDEF_DATA(number) == POS_INFINITY)
 706                         return number;
 707                 else if (XINDEF_DATA(number) == NEG_INFINITY)
 708                         return Fcoerce_number(Qzero, Qfloat, precision);
 709                 else
 710                         return number;
 711         }
 712
 713         Fsignal(Qarith_error, list1(number));
 714         return Qnil;
 715
 716         if (NILP(precision));
 717 #endif  /* HAVE_MPFR */
 718 }
 719 #endif  /* HAVE_FPFLOAT || HAVE_MPFR */
 720
 721 DEFUN("2^", Fexp2, 1, 2, 0,     /*
 722 Return the exponential of NUMBER to 2 power.
 723 If optional argument PRECISION is non-nil, its value
 724 \(an integer\) is used as precision in float computations.
 725 */
 726       (number, precision))
 727 {
 728 #if defined HAVE_MPFR && defined WITH_MPFR
 729         Lisp_Object bfrnumber;
 730 #endif
 731 #if defined HAVE_MPZ && (defined WITH_GMP || defined WITH_MP)
 732         if (INTP(number))
 733                 return _ent_binop(ASE_BINARY_OP_POW,
 734                                   INT_T, make_int(2), INT_T, number);
 735 #endif
 736         if (INDEFP(number))
 737                 return _ent_binop(ASE_BINARY_OP_POW,
 738                                   INT_T, make_int(2), INDEF_T, number);
 739
 740 #if defined HAVE_MPFR && defined WITH_MPFR
 741         bfrnumber = Fcoerce_number(number, Qbigfr, Qnil);
 742         bigfr_set_prec(ent_scratch_bigfr,
 743                        internal_get_precision(precision));
 744
 745         bigfr_exp2(ent_scratch_bigfr, XBIGFR_DATA(bfrnumber));
 746         return make_bigfr_bfr(ent_scratch_bigfr);
 747 #endif
 748         /* fallback */
 749         if (NILP(precision));
 750         return Qnil;
 751 }
 752
 753 DEFUN("10^", Fexp10, 1, 2, 0,   /*
 754 Return the exponential of NUMBER to 10 power.
 755 If optional argument PRECISION is non-nil, its value
 756 \(an integer\) is used as precision in float computations.
 757 */
 758       (number, precision))
 759 {
 760 #if defined HAVE_MPFR && defined WITH_MPFR
 761         Lisp_Object bfrnumber;
 762 #endif
 763 #if defined HAVE_MPZ && (defined WITH_GMP || defined WITH_MP)
 764         if (INTP(number))
 765                 return _ent_binop(ASE_BINARY_OP_POW,
 766                                   INT_T, make_int(10), INT_T, number);
 767 #endif
 768         if (INDEFP(number))
 769                 return _ent_binop(ASE_BINARY_OP_POW,
 770                                   INT_T, make_int(10), INDEF_T, number);
 771
 772 #if defined HAVE_MPFR && defined WITH_MPFR
 773         bfrnumber = Fcoerce_number(number, Qbigfr, Qnil);
 774         bigfr_set_prec(ent_scratch_bigfr,
 775                        internal_get_precision(precision));
 776
 777         bigfr_exp10(ent_scratch_bigfr, XBIGFR_DATA(bfrnumber));
 778         return make_bigfr_bfr(ent_scratch_bigfr);
 779 #endif
 780         /* fallback */
 781         if (NILP(precision));
 782         return Qnil;
 783 }
 784
 785 #if defined(HAVE_FPFLOAT) || defined(HAVE_MPFR) && defined WITH_MPFR
 786
 787 #define RETURN_WHEN_INDEF(number)                                   \
 788         if (INDEFP(number)) {                                       \
 789                 if (XINDEF_DATA(number) == POS_INFINITY) {          \
 790                         return number;                              \
 791                 } else if (XINDEF_DATA(number) == NEG_INFINITY) {       \
 792                         return make_indef(NOT_A_NUMBER);            \
 793                 } else {                                            \
 794                         return number;                              \
 795                 }                                                   \
 796         }
 797
 798 DEFUN("log10", Flog10, 1, 2, 0, /*
 799 Return the logarithm base 10 of NUMBER.
 800 If second optional argument PRECISION is given, use its value
 801 (an integer) as precision.
 802 */
 803       (number, precision))
 804 {
 805         RETURN_WHEN_INDEF(number);
 806
 807 #if defined HAVE_MPFR && defined WITH_MPFR
 808         Lisp_Object bfrnumber;
 809
 810         bigfr_set_prec(ent_scratch_bigfr,
 811                        internal_get_precision(precision));
 812
 813         bfrnumber = Fcoerce_number(number, Qbigfr, Qnil);
 814         bigfr_log10(ent_scratch_bigfr, XBIGFR_DATA(bfrnumber));
 815         return make_bigfr_bfr(ent_scratch_bigfr);
 816
 817 #else  /* !HAVE_MPFR */
 818         number = ent_lift(number, FLOAT_T, NULL);
 819
 820         RETURN_WHEN_INDEF(number);
 821
 822         if (FLOATP(number)) {
 823                 fpfloat d;
 824 #if HAVE_LOG10
 825                 d = log10(XFLOAT_DATA(number));
 826                 return make_float(d);
 827 #elif HAVE_LOG2
 828         static const fpflot log2_10 = log2(10);
 829         d = log2(XFLOAT_DATA(number))/log2_10;
 830                 return make_float(d);
 831 #elif HAVE_LOG
 832         static const fpflot log_10 - log(10);
 833         d = log(XFLOAT_DATA(number))/log_10;
 834                 return make_float(d);
 835 #else
 836         return ase_unary_operation_undefined(number);
 837 #endif
 838         }
 839
 840         Fsignal(Qarith_error, list1(number));
 841         return Qnil;
 842
 843         if (NILP(precision));
 844 #endif  /* HAVE_MPFR */
 845 }
 846
 847 DEFUN("log2", Flog2, 1, 2, 0,   /*
 848 Return the logarithm base 2 of NUMBER.
 849 If second optional argument PRECISION is given, use its value
 850 (an integer) as precision.
 851 */
 852       (number, precision))
 853 {
 854         RETURN_WHEN_INDEF(number);
 855
 856 #if defined HAVE_MPFR && defined WITH_MPFR
 857         Lisp_Object bfrnumber;
 858
 859         bigfr_set_prec(ent_scratch_bigfr,
 860                        internal_get_precision(precision));
 861
 862         bfrnumber = Fcoerce_number(number, Qbigfr, Qnil);
 863         bigfr_log2(ent_scratch_bigfr, XBIGFR_DATA(bfrnumber));
 864         return make_bigfr_bfr(ent_scratch_bigfr);
 865 #else
 866         number = ent_lift(number, FLOAT_T, NULL);
 867
 868         RETURN_WHEN_INDEF(number);
 869
 870         if (FLOATP(number)) {
 871                 fpfloat d;
 872 #if HAVE_LOG2
 873                 d = log2(XFLOAT_DATA(number));
 874                 return make_float(d);
 875 #elif HAVE_LOG
 876         static const fpflot log_2 - log(2);
 877         d = log(XFLOAT_DATA(number))/log_2;
 878                 return make_float(d);
 879 #else
 880         return ase_unary_operation_undefined(number);
 881 #endif
 882         }
 883
 884         Fsignal(Qarith_error, list1(number));
 885         return Qnil;
 886
 887         if (NILP(precision));
 888 #endif  /* HAVE_MPFR */
 889 }
 890
 891 DEFUN("log", Flog, 1, 3, 0,     /*
 892 Return the natural logarithm of NUMBER.
 893 If second optional argument BASE is given, return the logarithm of
 894 NUMBER using that base.
 895 If third optional argument PRECISION is given, use its value
 896 (an integer) as precision.
 897 */
 898       (number, base, precision))
 899 {
 900         RETURN_WHEN_INDEF(number);
 901
 902         if (INTEGERP(base)) {
 903                 switch(XINT(base)) {
 904                 case 2 : return Flog2 (number, precision);
 905                 case 10: return Flog10(number, precision);
 906                 default: break; /* Intentional Fall through */
 907                 }
 908         }
 909
 910
 911 #if defined HAVE_MPFR && defined WITH_MPFR
 912         if (!NILP(base)) {
 913                 /* Not all bignumber libs optimize log2, for instance
 914                    MPFR implements log2 in function of log. */
 915                 Lisp_Object _logn, _logb;
 916                 _logn = Flog(number, precision);
 917                 if (UNLIKELY(INDEFP(_logn))) {
 918                         return _logn;
 919                 }
 920                 _logb = Flog2(base, precision);
 921                 return ent_binop(ASE_BINARY_OP_QUO, _logn, _logb);
 922         }
 923
 924         Lisp_Object bfrnumber;
 925
 926         bigfr_set_prec(ent_scratch_bigfr,
 927                        internal_get_precision(precision));
 928
 929         bfrnumber = Fcoerce_number(number, Qbigfr, Qnil);
 930         bigfr_log(ent_scratch_bigfr, XBIGFR_DATA(bfrnumber));
 931         return make_bigfr_bfr(ent_scratch_bigfr);
 932
 933 #else  /* !HAVE_MPFR */
 934         if (!NILP(base)) {
 935         /* Processor implementations tend to give an edge to log2 */
 936                 Lisp_Object _logn, _logb;
 937                 _logn = Flog2(number, precision);
 938                 if (UNLIKELY(INDEFP(_logn))) {
 939                         return _logn;
 940                 }
 941                 _logb = Flog2(base, precision);
 942                 return ent_binop(ASE_BINARY_OP_QUO, _logn, _logb);
 943         }
 944
 945         number = ent_lift(number, FLOAT_T, NULL);
 946
 947         RETURN_WHEN_INDEF(number);
 948
 949         if (FLOATP(number)) {
 950                 fpfloat d;
 951                 d = log(XFLOAT_DATA(number));
 952                 return make_float(d);
 953         }
 954
 955         Fsignal(Qarith_error, list1(number));
 956         return Qnil;
 957
 958         if (NILP(precision));
 959 #endif  /* HAVE_MPFR */
 960 }
 961
 962 #undef RETURN_WHEN_INDEF
 963
 964 DEFUN("sqrt", Fsqrt, 1, 2, 0,   /*
 965 Return the square root of NUMBER.
 966 If second optional argument PRECISION is given, use its value
 967 (an integer) as precision.
 968 */
 969       (number, precision))
 970 {
 971 #if defined(HAVE_MPFR) && defined WITH_MPFR ||  \
 972         defined(HAVE_MPC) && defined WITH_MPC ||        \
 973         defined(HAVE_PSEUC) && defined WITH_PSEUC
 974
 975         if (INDEFP(number)) {
 976                 if (XINDEF_DATA(number) == POS_INFINITY)
 977                         return number;
 978                 else if (XINDEF_DATA(number) == NEG_INFINITY)
 979                         return make_indef(COMPLEX_INFINITY);
 980                 else
 981                         return number;
 982         }
 983
 984         if (COMPARABLEP(number)) {
 985 #if defined HAVE_MPFR && defined WITH_MPFR
 986                 bigfr_set_prec(ent_scratch_bigfr,
 987                                internal_get_precision(precision));
 988
 989                 if (NATNUMP(number))
 990                         bigfr_sqrt_ui(ent_scratch_bigfr,
 991                                       (unsigned long)XUINT(number));
 992                 else if (BIGZP(number) &&
 993                          bigz_fits_ulong_p(XBIGZ_DATA(number)) &&
 994                          bigz_sign(XBIGZ_DATA(number)) >= 0) {
 995                         bigfr_sqrt_ui(ent_scratch_bigfr,
 996                                       (unsigned long)bigz_to_ulong(
 997                                               XBIGZ_DATA(number)));
 998                 } else if (!NILP(Fnonnegativep(number))) {
 999                         Lisp_Object bfrnumber;
1000                         bfrnumber = Fcoerce_number(number, Qbigfr, Qnil);
1001                         bigfr_sqrt(ent_scratch_bigfr, XBIGFR_DATA(bfrnumber));
1002                 } else {
1003 #if defined HAVE_MPC && defined WITH_MPC ||     \
1004         defined HAVE_PSEUC && defined WITH_PSEUC
1005                         Lisp_Object bcnumber;
1006                         bigc_set_prec(ent_scratch_bigc,
1007                                       internal_get_precision(precision));
1008                         bcnumber = Fcoerce_number(number, Qbigc, precision);
1009                         bigc_sqrt(ent_scratch_bigc, XBIGC_DATA(bcnumber));
1010                         return make_bigc_bc(ent_scratch_bigc);
1011 #else  /* !HAVE_MPC */
1012                         Lisp_Object bfrnumber;
1013                         bfrnumber = Fcoerce_number(number, Qbigfr, Qnil);
1014                         bigfr_sqrt(ent_scratch_bigfr, XBIGFR_DATA(bfrnumber));
1015 #endif  /* HAVE_MPC */
1016                 }
1017                 return make_bigfr_bfr(ent_scratch_bigfr);
1018 #endif  /* HAVE_MPFR */
1019 #if defined HAVE_MPC && defined WITH_MPC ||     \
1020         defined HAVE_PSEUC && defined WITH_PSEUC
1021         } else if (BIGCP(number) || BIGGP(number)) {
1022                 Lisp_Object bcnumber;
1023                 bigc_set_prec(ent_scratch_bigc,
1024                               internal_get_precision(precision));
1025
1026                 bcnumber = Fcoerce_number(number, Qbigc, precision);
1027                 bigc_sqrt(ent_scratch_bigc, XBIGC_DATA(bcnumber));
1028                 return make_bigc_bc(ent_scratch_bigc);
1029 #endif  /* HAVE_MPC */
1030         }
1031
1032         if (NILP(precision));
1033         return wrong_type_argument(Qnumberp, number);
1034
1035 #else  /* !HAVE_MPFR && !HAVE_MPC */
1036         if (INDEFP(number)) {
1037                 goto indefcase;
1038         }
1039
1040         number = ent_lift(number, FLOAT_T, NULL);
1041
1042         if (FLOATP(number)) {
1043                 fpfloat d;
1044                 d = sqrt(XFLOAT_DATA(number));
1045                 return make_float(d);
1046         } else if (INDEFP(number)) {
1047         indefcase:
1048                 if (XINDEF_DATA(number) == POS_INFINITY)
1049                         return number;
1050                 else if (XINDEF_DATA(number) == NEG_INFINITY)
1051                         return make_indef(COMPLEX_INFINITY);
1052                 else
1053                         return number;
1054         }
1055
1056         Fsignal(Qarith_error, list1(number));
1057         return Qnil;
1058
1059         if (NILP(precision));
1060 #endif  /* HAVE_MPFR */
1061 }
1062
1063 DEFUN("cube-root", Fcube_root, 1, 2, 0, /*
1064 Return the cube root of NUMBER.
1065 If second optional argument PRECISION is given, use its value
1066 (an integer) as precision.
1067 */
1068       (number, precision))
1069 {
1070 #if defined HAVE_MPFR && defined WITH_MPFR
1071         Lisp_Object bfrnumber;
1072
1073         if (INDEFP(number))
1074                 return number;
1075
1076         bigfr_set_prec(ent_scratch_bigfr,
1077                        internal_get_precision(precision));
1078
1079         bfrnumber = Fcoerce_number(number, Qbigfr, Qnil);
1080         bigfr_cbrt(ent_scratch_bigfr, XBIGFR_DATA(bfrnumber));
1081         return make_bigfr_bfr(ent_scratch_bigfr);
1082
1083 #else  /* !HAVE_MPFR */
1084         if (INDEFP(number)) {
1085                 goto indefcase;
1086         }
1087
1088         number = ent_lift(number, FLOAT_T, NULL);
1089
1090         if (FLOATP(number)) {
1091                 fpfloat d;
1092 #ifdef HAVE_CBRT
1093                 d = cbrt(XFLOAT_DATA(number));
1094 #else
1095                 d = XFLOAT_DATA(number);
1096                 if (d >= 0.0)
1097                         d = pow(d, 1.0 / 3.0);
1098                 else
1099                         d = -pow(-d, 1.0 / 3.0);
1100 #endif
1101                 return make_float(d);
1102         } else if (INDEFP(number)) {
1103         indefcase:
1104                 return number;
1105         }
1106
1107         Fsignal(Qarith_error, list1(number));
1108         return Qnil;
1109
1110         if (NILP(precision));
1111 #endif  /* HAVE_MPFR */
1112 }
1113 #endif  /* HAVE_FPFLOAT || MPFR */
1114
1115
1116 #if defined HAVE_MPFR && defined WITH_MPFR
1117 DEFUN("root", Froot, 2, 3, 0,   /*
1118 Return the RADIX-th root of NUMBER.
1119 If third optional argument PRECISION is given, use its value
1120 (an integer) as precision.
1121 */
1122       (number, radix, precision))
1123 {
1124         Lisp_Object bfrnumber;
1125
1126         if (!NATNUMP(radix)) {
1127                 dead_wrong_type_argument(Qnatnump, radix);
1128                 return Qnil;
1129         }
1130
1131         if (INDEFP(number)) {
1132                 if (XINDEF_DATA(number) == POS_INFINITY)
1133                         return number;
1134                 else if (XINDEF_DATA(number) == NEG_INFINITY)
1135                         return make_indef(COMPLEX_INFINITY);
1136                 else
1137                         return number;
1138         }
1139
1140         bigfr_set_prec(ent_scratch_bigfr,
1141                        internal_get_precision(precision));
1142
1143         bfrnumber = Fcoerce_number(number, Qbigfr, Qnil);
1144         bigfr_root(ent_scratch_bigfr, XBIGFR_DATA(bfrnumber), XUINT(radix));
1145         return make_bigfr_bfr(ent_scratch_bigfr);
1146 }
1147 #endif  /* HAVE_MPFR */
1148
1149 \f
1150 /* (Inverse) hyperbolic trig functions. */
1151 #if defined(HAVE_FPFLOAT) || defined(HAVE_MPFR) && defined WITH_MPFR
1152
1153 DEFUN("acosh", Facosh, 1, 2, 0, /*
1154 Return the inverse hyperbolic cosine of NUMBER.
1155 If optional argument PRECISION is non-nil, its value
1156 (an integer) is used as precision.
1157 */
1158       (number, precision))
1159 {
1160 #if defined HAVE_MPFR && defined WITH_MPFR
1161
1162         MPFR_TRIG_FUN(acosh);
1163
1164 #else  /* !HAVE_MPFR */
1165         if (INDEFP(number)) {
1166                 return make_indef(NOT_A_NUMBER);
1167         }
1168
1169         number = ent_lift(number, FLOAT_T, NULL);
1170
1171         if (FLOATP(number)) {
1172                 fpfloat d = XFLOAT_DATA(number);
1173 #ifdef HAVE_INVERSE_HYPERBOLIC
1174                 d = acosh(d);
1175 #else
1176                 d = log(d + sqrt(d * d - 1.0));
1177 #endif
1178                 return make_float(d);
1179         } else if (INDEFP(number)) {
1180                 return make_indef(NOT_A_NUMBER);
1181         }
1182
1183         Fsignal(Qarith_error, list1(number));
1184         return Qnil;
1185
1186         if (NILP(precision));
1187 #endif  /* HAVE_MPFR */
1188 }
1189
1190 DEFUN("asinh", Fasinh, 1, 2, 0, /*
1191 Return the inverse hyperbolic sine of NUMBER.
1192 If optional argument PRECISION is non-nil, its value
1193 (an integer) is used as precision.
1194 */
1195       (number, precision))
1196 {
1197 #if defined HAVE_MPFR && defined WITH_MPFR
1198
1199         MPFR_TRIG_FUN(asinh);
1200
1201 #else  /* !HAVE_MPFR */
1202         if (INDEFP(number)) {
1203                 return make_indef(NOT_A_NUMBER);
1204         }
1205
1206         number = ent_lift(number, FLOAT_T, NULL);
1207
1208         if (FLOATP(number)) {
1209                 fpfloat d = XFLOAT_DATA(number);
1210 #ifdef HAVE_INVERSE_HYPERBOLIC
1211                 d = acosh(d);
1212 #else
1213                 d = log(d + sqrt(d * d + 1.0));
1214 #endif
1215                 return make_float(d);
1216         } else if (INDEFP(number)) {
1217                 return make_indef(NOT_A_NUMBER);
1218         }
1219
1220         Fsignal(Qarith_error, list1(number));
1221         return Qnil;
1222
1223         if (NILP(precision));
1224 #endif  /* HAVE_MPFR */
1225 }
1226
1227 DEFUN("atanh", Fatanh, 1, 2, 0, /*
1228 Return the inverse hyperbolic tangent of NUMBER.
1229 If optional argument PRECISION is non-nil, its value
1230 (an integer) is used as precision.
1231 */
1232       (number, precision))
1233 {
1234 #if defined HAVE_MPFR && defined WITH_MPFR
1235
1236         MPFR_TRIG_FUN(atanh);
1237
1238 #else  /* !HAVE_MPFR */
1239         if (INDEFP(number)) {
1240                 return make_indef(NOT_A_NUMBER);
1241         }
1242
1243         number = ent_lift(number, FLOAT_T, NULL);
1244
1245         if (FLOATP(number)) {
1246                 fpfloat d = XFLOAT_DATA(number);
1247 #ifdef HAVE_INVERSE_HYPERBOLIC
1248                 d = atanh(d);
1249 #else
1250                 d = 0.5 * log((1.0 + d) / (1.0 - d));
1251 #endif
1252                 return make_float(d);
1253         } else if (INDEFP(number)) {
1254                 return make_indef(NOT_A_NUMBER);
1255         }
1256
1257         Fsignal(Qarith_error, list1(number));
1258         return Qnil;
1259
1260         if (NILP(precision));
1261 #endif  /* HAVE_MPFR */
1262 }
1263
1264 DEFUN("cosh", Fcosh, 1, 2, 0,   /*
1265 Return the hyperbolic cosine of NUMBER.
1266 If optional argument PRECISION is non-nil, its value
1267 (an integer) is used as precision.
1268 */
1269       (number, precision))
1270 {
1271 #if defined HAVE_MPFR && defined WITH_MPFR
1272
1273         MPFR_TRIG_FUN(cosh);
1274
1275 #else  /* !HAVE_MPFR */
1276         if (INDEFP(number)) {
1277                 return make_indef(NOT_A_NUMBER);
1278         }
1279
1280         number = ent_lift(number, FLOAT_T, NULL);
1281
1282         if (FLOATP(number)) {
1283                 fpfloat d;
1284                 d = cosh(XFLOAT_DATA(number));
1285                 return make_float(d);
1286         } else if (INDEFP(number)) {
1287                 return make_indef(NOT_A_NUMBER);
1288         }
1289
1290         Fsignal(Qarith_error, list1(number));
1291         return Qnil;
1292
1293         if (NILP(precision));
1294 #endif  /* HAVE_MPFR */
1295 }
1296
1297 DEFUN("sinh", Fsinh, 1, 2, 0,   /*
1298 Return the hyperbolic sine of NUMBER.
1299 If optional argument PRECISION is non-nil, its value
1300 (an integer) is used as precision.
1301 */
1302       (number, precision))
1303 {
1304 #if defined HAVE_MPFR && defined WITH_MPFR
1305
1306         MPFR_TRIG_FUN(sinh);
1307
1308 #else  /* !HAVE_MPFR */
1309         if (INDEFP(number)) {
1310                 return make_indef(NOT_A_NUMBER);
1311         }
1312
1313         number = ent_lift(number, FLOAT_T, NULL);
1314
1315         if (FLOATP(number)) {
1316                 fpfloat d;
1317                 d = sinh(XFLOAT_DATA(number));
1318                 return make_float(d);
1319         } else if (INDEFP(number)) {
1320                 return make_indef(NOT_A_NUMBER);
1321         }
1322
1323         Fsignal(Qarith_error, list1(number));
1324         return Qnil;
1325
1326         if (NILP(precision));
1327 #endif  /* HAVE_MFPR */
1328 }
1329
1330 DEFUN("tanh", Ftanh, 1, 2, 0,   /*
1331 Return the hyperbolic tangent of NUMBER.
1332 If optional argument PRECISION is non-nil, its value
1333 (an integer) is used as precision.
1334 */
1335       (number, precision))
1336 {
1337 #if defined HAVE_MPFR && defined WITH_MPFR
1338
1339         MPFR_TRIG_FUN(tanh);
1340
1341 #else  /* !HAVE_MPFR */
1342         if (INDEFP(number)) {
1343                 return make_indef(NOT_A_NUMBER);
1344         }
1345
1346         number = ent_lift(number, FLOAT_T, NULL);
1347
1348         if (FLOATP(number)) {
1349                 fpfloat d = XFLOAT_DATA(number);
1350                 d = tanh(d);
1351                 return make_float(d);
1352         } else if (INDEFP(number)) {
1353                 return make_indef(NOT_A_NUMBER);
1354         }
1355
1356         Fsignal(Qarith_error, list1(number));
1357         return Qnil;
1358
1359         if (NILP(precision));
1360 #endif  /* HAVE_MPFR */
1361 }
1362
1363 #if defined HAVE_MPFR && defined WITH_MPFR
1364
1365 DEFUN("sech", Fsech, 1, 2, 0,   /*
1366 Return the hyperbolic secant of NUMBER.
1367 If optional argument PRECISION is non-nil, its value
1368 (an integer) is used as precision.
1369 */
1370       (number, precision))
1371 {
1372         MPFR_TRIG_FUN(sech);
1373 }
1374
1375 DEFUN("csch", Fcsch, 1, 2, 0,   /*
1376 Return the hyperbolic cosecant of NUMBER.
1377 If optional argument PRECISION is non-nil, its value
1378 (an integer) is used as precision.
1379 */
1380       (number, precision))
1381 {
1382         MPFR_TRIG_FUN(csch);
1383 }
1384
1385 DEFUN("coth", Fcoth, 1, 2, 0,   /*
1386 Return the hyperbolic cotangent of NUMBER.
1387 If optional argument PRECISION is non-nil, its value
1388 (an integer) is used as precision.
1389 */
1390       (number, precision))
1391 {
1392         MPFR_TRIG_FUN(coth);
1393 }
1394 #endif  /* HAVE_MPFR */
1395
1396 #endif  /* HAVE_MPFR || HAVE_FPFLOAT (inverse trig functions) */
1397
1398 \f
1399 /* Rounding functions */
1400
1401 DEFUN("abs", Fabs, 1, 1, 0,     /*
1402 Return the absolute value of NUMBER.
1403 */
1404       (number))
1405 {
1406 #ifdef HAVE_FPFLOAT
1407         if (FLOATP(number)) {
1408                 return make_float(fabs(XFLOAT_DATA(number)));
1409         }
1410 #endif                          /* HAVE_FPFLOAT */
1411
1412         if (INTP(number)) {
1413 #if defined HAVE_MPZ && (defined WITH_GMP || defined WITH_MP)
1414                 /* The most negative Lisp int will overflow */
1415                 return (XINT(number) >= 0)
1416                         ? number : make_integer(-XINT(number));
1417 #else  /* !HAVE_MPZ */
1418                 return (XINT(number) >= 0) ? number : make_int(-XINT(number));
1419 #endif  /* HAVE_MPZ */
1420         }
1421
1422 #if defined HAVE_MPZ && (defined WITH_GMP || defined WITH_MP)
1423         if (BIGZP(number)) {
1424                 if (bigz_sign(XBIGZ_DATA(number)) >= 0)
1425                         return number;
1426
1427                 bigz_abs(ent_scratch_bigz, XBIGZ_DATA(number));
1428                 return make_bigz_bz(ent_scratch_bigz);
1429         }
1430 #endif  /* HAVE_MPZ */
1431
1432 #if defined HAVE_MPQ && defined WITH_GMP
1433         if (BIGQP(number)) {
1434                 if (bigq_sign(XBIGQ_DATA(number)) >= 0)
1435                         return number;
1436
1437                 bigq_abs(ent_scratch_bigq, XBIGQ_DATA(number));
1438                 return make_bigq_bq(ent_scratch_bigq);
1439         }
1440 #endif  /* HAVE_MPQ */
1441
1442 #if defined HAVE_MPF && defined WITH_GMP
1443         if (BIGFP(number)) {
1444                 if (bigf_sign(XBIGF_DATA (number)) >= 0)
1445                         return number;
1446
1447                 bigf_set_prec(ent_scratch_bigf, XBIGF_GET_PREC(number));
1448
1449                 bigf_abs(ent_scratch_bigf, XBIGF_DATA(number));
1450                 return make_bigf_bf(ent_scratch_bigf);
1451         }
1452 #endif  /* HAVE_MPF */
1453
1454 #if defined HAVE_MPFR && defined WITH_MPFR
1455         if (BIGFRP(number)) {
1456                 if (bigfr_sign(XBIGFR_DATA (number)) >= 0)
1457                         return number;
1458
1459                 bigfr_set_prec(ent_scratch_bigfr, XBIGFR_GET_PREC(number));
1460
1461                 bigfr_abs(ent_scratch_bigfr, XBIGFR_DATA(number));
1462                 return make_bigfr_bfr(ent_scratch_bigfr);
1463         }
1464 #endif  /* HAVE_MPFR */
1465
1466 #if defined(HAVE_PSEUG) && defined WITH_PSEUG && defined(HAVE_MPFR)
1467         if (BIGGP(number)) {
1468                 bigfr_set_prec(ent_scratch_bigfr,
1469                                internal_get_precision(Qnil));
1470
1471                 bigg_abs(ent_scratch_bigfr, XBIGG_DATA(number));
1472                 return make_bigfr_bfr(ent_scratch_bigfr);
1473         }
1474 #endif  /* HAVE_PSEUG && HAVE_MPFR */
1475
1476 #if defined HAVE_MPC && defined WITH_MPC ||     \
1477         defined HAVE_PSEUC && defined WITH_PSEUC
1478         if (BIGCP(number)) {
1479                 bigfr_set_prec(ent_scratch_bigfr, XBIGC_GET_PREC(number));
1480
1481                 if (bigc_nan_p(XBIGC_DATA(number)))
1482                         return make_indef(NOT_A_NUMBER);
1483                 else if (bigc_inf_p(XBIGC_DATA(number)))
1484                         return make_indef(POS_INFINITY);
1485                 else
1486                         bigc_abs(ent_scratch_bigfr, XBIGC_DATA(number));
1487
1488                 return make_bigfr_bfr(ent_scratch_bigfr);
1489         }
1490 #endif  /* HAVE_PSEUG */
1491
1492         if (INDEFP(number)) {
1493                 if (XINDEF_DATA(number) == POS_INFINITY)
1494                         return number;
1495                 else if (XINDEF_DATA(number) == NEG_INFINITY)
1496                         return make_indef(POS_INFINITY);
1497                 else
1498                         return number;
1499         }
1500
1501         return Fabs(wrong_type_argument(Qnumberp, number));
1502 }
1503
1504 #if defined(HAVE_FPFLOAT)
1505 /* fuck fuck fuck, I want this in number.el */
1506 DEFUN("float", Ffloat, 1, 1, 0, /*
1507 Return the floating point number numerically equal to NUMBER.
1508 */
1509       (number))
1510 {
1511         /* Just create the float in order of preference */
1512         return Fcoerce_number(number, Qfloat, Qnil);
1513 }
1514 #endif  /* HAVE_FPFLOAT */
1515
1516 #ifdef HAVE_FPFLOAT
1517 DEFUN("logb", Flogb, 1, 1, 0,   /*
1518 Return largest integer <= the base 2 log of the magnitude of NUMBER.
1519 This is the same as the exponent of a float.
1520 */
1521       (number))
1522 {
1523         fpfloat f = extract_float(number);
1524
1525         if (f == 0.0)
1526                 return make_int(EMACS_INT_MIN);
1527 #ifdef HAVE_LOGB
1528         {
1529                 fpfloat _lb = logb(f);
1530                 Lisp_Object val;
1531                 IN_FLOAT(val = make_int((EMACS_INT)_lb), "logb", number);
1532                 return val;
1533         }
1534 #else
1535 #ifdef HAVE_FREXP
1536         {
1537                 int exqp;
1538                 IN_FLOAT(frexp(f, &exqp), "logb", number);
1539                 return make_int(exqp - 1);
1540         }
1541 #else
1542         {
1543                 int i;
1544                 fpfloat d;
1545                 EMACS_INT val;
1546                 if (f < 0.0)
1547                         f = -f;
1548                 val = -1;
1549                 while (f < 0.5) {
1550                         for (i = 1, d = 0.5; d * d >= f; i += i)
1551                                 d *= d;
1552                         f /= d;
1553                         val -= i;
1554                 }
1555                 while (f >= 1.0) {
1556                         for (i = 1, d = 2.0; d * d <= f; i += i)
1557                                 d *= d;
1558                         f /= d;
1559                         val += i;
1560                 }
1561                 return make_int(val);
1562         }
1563 #endif                          /* ! HAVE_FREXP */
1564 #endif                          /* ! HAVE_LOGB */
1565 }
1566 #endif                          /* HAVE_FPFLOAT */
1567
1568 DEFUN("ceiling", Fceiling, 1, 1, 0,     /*
1569 Return the smallest integer no less than NUMBER.  (Round toward +inf.)
1570 */
1571       (number))
1572 {
1573 #ifdef HAVE_FPFLOAT
1574         if (FLOATP(number)) {
1575                 fpfloat d;
1576                 d = ceil(XFLOAT_DATA(number));
1577                 return (float_to_int(d, "ceiling", number, Qunbound));
1578         }
1579 #endif                          /* HAVE_FPFLOAT */
1580
1581 #if defined HAVE_MPZ && (defined WITH_GMP || defined WITH_MP)
1582         if (INTEGERP(number))
1583 #else  /* !HAVE_MPZ */
1584         if (INTP(number))
1585 #endif  /* HAVE_MPZ */
1586                 return number;
1587
1588 #if defined(HAVE_MPQ) && defined(HAVE_MPZ) && defined WITH_GMP
1589         if (BIGQP(number)) {
1590                 bigz_ceil(ent_scratch_bigz,
1591                           XBIGQ_NUMERATOR(number),
1592                           XBIGQ_DENOMINATOR(number));
1593                 return ent_mpz_downgrade_maybe(ent_scratch_bigz);
1594         }
1595 #endif
1596
1597 #if defined HAVE_MPF && defined WITH_GMP
1598         else if (BIGFP(number)) {
1599 #if defined HAVE_MPZ && (defined WITH_GMP || defined WITH_MP)
1600                 bigf_ceil(ent_scratch_bigf, XBIGF_DATA(number));
1601                 bigz_set_bigf(ent_scratch_bigz, ent_scratch_bigf);
1602                 return ent_mpz_downgrade_maybe(ent_scratch_bigz);
1603 #else  /* !HAVE_MPZ */
1604                 bigf_ceil(ent_scratch_bigf, XBIGF_DATA(number));
1605                 return make_int((EMACS_INT)bigf_to_long(ent_scratch_bigf));
1606 #endif  /* HAVE_MPZ */
1607         }
1608 #endif  /* HAVE_MPF */
1609
1610 #if defined HAVE_MPFR && defined WITH_MPFR
1611         else if (BIGFRP(number)) {
1612 #if defined HAVE_MPZ && (defined WITH_GMP || defined WITH_MP)
1613                 bigfr_ceil(ent_scratch_bigfr, XBIGFR_DATA(number));
1614                 bigz_set_bigfr(ent_scratch_bigz, ent_scratch_bigfr);
1615                 return ent_mpz_downgrade_maybe(ent_scratch_bigz);
1616 #else  /* !HAVE_MPZ */
1617                 bigfr_ceil(ent_scratch_bigfr, XBIGFR_DATA(number));
1618                 return make_int((EMACS_INT)bigfr_to_long(ent_scratch_bigfr));
1619 #endif  /* HAVE_MPZ */
1620         }
1621 #endif  /* HAVE_MPFR */
1622
1623         if (INDEFP(number))
1624                 return number;
1625
1626 #if defined HAVE_MPC && defined WITH_MPC ||     \
1627         defined HAVE_PSEUC && defined WITH_PSEUC ||     \
1628         defined HAVE_PSEUG && defined WITH_PSEUG
1629         return Fceiling(wrong_type_argument(Qcomparablep, number));
1630 #else  /* !HAVE_MPC */
1631         return Fceiling(wrong_type_argument(Qnumberp, number));
1632 #endif  /* HAVE_MPC */
1633 }
1634
1635 DEFUN("floor", Ffloor, 1, 2, 0, /*
1636 Return the largest integer no greater than NUMBER.  (Round towards -inf.)
1637 With optional second argument DIVISOR, return the largest integer no
1638 greater than NUMBER/DIVISOR.
1639 */
1640       (number, divisor))
1641 {
1642         ase_object_type_t ntquo;
1643         Lisp_Object quo;
1644
1645         CHECK_COMPARABLE(number);
1646         if (NILP(divisor)) {
1647                 return Ffloor(number, make_int(1L));
1648
1649         }
1650
1651         /* !NILP(divisor) */
1652
1653         CHECK_COMPARABLE(divisor);
1654
1655         if (INTEGERP(number) && INTEGERP(divisor)) {
1656 #if defined HAVE_MPZ && (defined WITH_GMP || defined WITH_MP)
1657                 /* this is the optimised version, since
1658                  * bigz_floor always takes two arguments
1659                  */
1660                 number = ent_lift(number, BIGZ_T, NULL);
1661                 divisor = ent_lift(divisor, BIGZ_T, NULL);
1662
1663                 bigz_floor(ent_scratch_bigz,
1664                            XBIGZ_DATA(number),
1665                            XBIGZ_DATA(divisor));
1666                 return ent_mpz_downgrade_maybe(ent_scratch_bigz);
1667 #else
1668                 number = ent_lift(number, FLOAT_T, NULL);
1669                 divisor = ent_lift(divisor, FLOAT_T, NULL);
1670 #endif
1671         }
1672
1673         quo = ent_binop(ASE_BINARY_OP_QUO, number, divisor);
1674         ntquo = ase_optable_index(quo);
1675
1676         switch (ntquo) {
1677         case INT_T:             /* trivial */
1678         case BIGZ_T:
1679         case INDEF_T:
1680                 return quo;
1681                 break;
1682         case FLOAT_T: {
1683                 fpfloat d;
1684                 IN_FLOAT((d = floor(XFLOAT_DATA(quo))), "floor", quo);
1685                 return (float_to_int(d, "floor", quo, Qunbound));
1686         }
1687         case BIGQ_T:
1688 #if defined(HAVE_MPQ) && defined(HAVE_MPZ) && defined WITH_GMP
1689                 bigz_floor(ent_scratch_bigz,
1690                            XBIGQ_NUMERATOR(quo), XBIGQ_DENOMINATOR(quo));
1691                 return ent_mpz_downgrade_maybe(ent_scratch_bigz);
1692                 break;
1693 #else
1694                 return quo;
1695 #endif
1696         case BIGF_T:
1697 #if defined HAVE_MPF && defined WITH_GMP
1698                 bigf_floor(ent_scratch_bigf, XBIGF_DATA(quo));
1699 #if defined HAVE_MPZ && (defined WITH_GMP || defined WITH_MP)
1700                 bigz_set_bigf(ent_scratch_bigz, ent_scratch_bigf);
1701                 return ent_mpz_downgrade_maybe(ent_scratch_bigz);
1702 #else  /* !HAVE_MPZ */
1703                 return make_int(
1704                         (EMACS_INT)bigf_to_long(ent_scratch_bigf));
1705 #endif  /* HAVE_MPZ */
1706                 break;
1707 #endif  /* HAVE_MPF */
1708
1709         case BIGFR_T:
1710 #if defined HAVE_MPFR && defined WITH_MPFR
1711                 bigfr_floor(ent_scratch_bigfr, XBIGFR_DATA(quo));
1712 #if defined HAVE_MPZ && (defined WITH_GMP || defined WITH_MP)
1713                 bigz_set_bigfr(ent_scratch_bigz, ent_scratch_bigfr);
1714                 return ent_mpz_downgrade_maybe(ent_scratch_bigz);
1715 #else  /* !HAVE_MPZ */
1716                 return make_int(
1717                         (EMACS_INT)bigfr_to_long(ent_scratch_bigfr));
1718 #endif  /* HAVE_MPZ */
1719                 break;
1720 #endif  /* HAVE_MPFR */
1721
1722         default:
1723                 return quo;
1724         }
1725
1726         return Fsignal(Qdomain_error, Qnil);
1727 }
1728
1729 DEFUN("round", Fround, 1, 1, 0, /*
1730 Return the nearest integer to NUMBER.
1731
1732 NUMBER has to have an archimedian valuation, #'round returns the
1733 integer z for which | number - z | is minimal.
1734 */
1735       (number))
1736 {
1737 #ifdef HAVE_FPFLOAT
1738         if (FLOATP(number)) {
1739                 fpfloat d;
1740                 /* Screw the prevailing rounding mode.  */
1741                 d = emacs_rint(XFLOAT_DATA(number));
1742                 return (float_to_int(d, "round", number, Qunbound));
1743         }
1744 #endif                          /* HAVE_FPFLOAT */
1745
1746 #if defined HAVE_MPZ && (defined WITH_GMP || defined WITH_MP)
1747         if (INTEGERP(number))
1748 #else  /* !HAVE_MPZ */
1749         if (INTP(number))
1750 #endif  /* HAVE_MPZ */
1751                 return number;
1752
1753 #if defined(HAVE_MPQ) && defined(HAVE_MPZ) && defined WITH_GMP
1754         else if (BIGQP(number)) {
1755                 /* first off, let's create the division, remainder as well */
1756                 /* fuck ugly? */
1757                 mpz_tdiv_qr(ent_scratch_bigz,
1758                             bigq_numerator(ent_scratch_bigq),
1759                             XBIGQ_NUMERATOR(number),
1760                             XBIGQ_DENOMINATOR(number));
1761
1762                 /* <- denom(number) * 2 */
1763                 mpz_mul_2exp(bigq_numerator(ent_scratch_bigq),
1764                              bigq_numerator(ent_scratch_bigq), 1);
1765
1766                 /* check if we had to add one */
1767                 if (mpz_cmpabs(bigq_numerator(ent_scratch_bigq),
1768                                XBIGQ_DENOMINATOR(number)) >= 0) {
1769                         /* >= ceil(denom(number) / 2) */
1770                         if (mpz_sgn(bigq_numerator(ent_scratch_bigq)) > 0) {
1771                                 mpz_add_ui(ent_scratch_bigz,
1772                                            ent_scratch_bigz, 1UL);
1773                         } else {
1774                                 mpz_sub_ui(ent_scratch_bigz,
1775                                            ent_scratch_bigz, 1UL);
1776                         }
1777                 }
1778                 return ent_mpz_downgrade_maybe(ent_scratch_bigz);
1779         }
1780 #endif  /* HAVE_MPQ && HAVE_MPZ */
1781
1782 #if defined HAVE_MPF && defined WITH_GMP
1783         else if (BIGFP(number)) {
1784                 warn_when_safe(Qbigf, Qnotice,
1785                                "rounding number of type 'bigf (mpf-floats)"
1786                                "not yet implemented");
1787                 return number;
1788         }
1789 #endif  /* HAVE_MPF */
1790
1791 #if defined HAVE_MPFR && defined WITH_MPFR
1792         else if (BIGFRP(number)) {
1793 #if defined HAVE_MPZ && (defined WITH_GMP || defined WITH_MP)
1794                 bigfr_rint(ent_scratch_bigfr, XBIGFR_DATA(number));
1795                 bigz_set_bigfr(ent_scratch_bigz, ent_scratch_bigfr);
1796                 return ent_mpz_downgrade_maybe(ent_scratch_bigz);
1797 #else  /* !HAVE_MPZ */
1798                 bigfr_rint(ent_scratch_bigfr, XBIGFR_DATA(number));
1799                 return make_int((EMACS_INT)bigfr_to_long(ent_scratch_bigfr));
1800 #endif  /* HAVE_MPZ */
1801         }
1802 #endif  /* HAVE_MPFR */
1803
1804         else if (INDEFP(number))
1805                 return number;
1806
1807 #if defined HAVE_MPC && defined WITH_MPC ||     \
1808         defined HAVE_PSEUC && defined WITH_PSEUC ||     \
1809         defined HAVE_PSEUG && defined WITH_PSEUG
1810         return Fround(wrong_type_argument(Qcomparablep, number));
1811 #else  /* !HAVE_MPC */
1812         return Fround(wrong_type_argument(Qnumberp, number));
1813 #endif  /* HAVE_MPC */
1814 }
1815
1816 DEFUN("truncate", Ftruncate, 1, 1, 0,   /*
1817 Truncate a floating point number to an integer.
1818 Rounds the value toward zero.
1819 */
1820       (number))
1821 {
1822 #ifdef HAVE_FPFLOAT
1823         if (FLOATP(number))
1824                 return float_to_int(XFLOAT_DATA(number), "truncate", number,
1825                                     Qunbound);
1826 #endif                          /* HAVE_FPFLOAT */
1827
1828 #if defined HAVE_MPZ && (defined WITH_GMP || defined WITH_MP)
1829         if (INTEGERP(number))
1830 #else  /* !HAVE_MPZ */
1831         if (INTP(number))
1832 #endif  /* HAVE_MPZ */
1833                 return number;
1834
1835 #if defined HAVE_MPQ && defined WITH_GMP
1836         else if (BIGQP(number)) {
1837                 bigz_div(ent_scratch_bigz,
1838                          XBIGQ_NUMERATOR(number),
1839                          XBIGQ_DENOMINATOR(number));
1840                 return ent_mpz_downgrade_maybe(ent_scratch_bigz);
1841         }
1842 #endif
1843
1844 #if defined HAVE_MPF && defined WITH_GMP
1845         else if (BIGFP(number)) {
1846 #if defined HAVE_MPZ && (defined WITH_GMP || defined WITH_MP)
1847                 bigf_trunc(ent_scratch_bigf, XBIGF_DATA(number));
1848                 bigz_set_bigf(ent_scratch_bigz, ent_scratch_bigf);
1849                 return ent_mpz_downgrade_maybe(ent_scratch_bigz);
1850 #else  /* !HAVE_MPZ */
1851                 bigf_trunc(ent_scratch_bigf, XBIGF_DATA(number));
1852                 return make_int((EMACS_INT)bigf_to_long(ent_scratch_bigf));
1853 #endif  /* HAVE_MPZ */
1854         }
1855 #endif  /* HAVE_MPF */
1856
1857 #if defined HAVE_MPFR && defined WITH_MPFR
1858         else if (BIGFRP(number)) {
1859 #if defined HAVE_MPZ && (defined WITH_GMP || defined WITH_MP)
1860                 bigfr_trunc(ent_scratch_bigfr, XBIGFR_DATA(number));
1861                 bigz_set_bigfr(ent_scratch_bigz, ent_scratch_bigfr);
1862                 return ent_mpz_downgrade_maybe(ent_scratch_bigz);
1863 #else  /* !HAVE_MPZ */
1864                 bigfr_trunc(ent_scratch_bigfr, XBIGFR_DATA(number));
1865                 return make_int((EMACS_INT)bigfr_to_long(ent_scratch_bigfr));
1866 #endif  /* HAVE_MPZ */
1867         }
1868 #endif  /* HAVE_MPFR */
1869
1870         else if (INDEFP(number))
1871                 return number;
1872
1873 #if defined HAVE_MPC && defined WITH_MPC ||     \
1874         defined HAVE_PSEUC && defined WITH_PSEUC ||     \
1875         defined HAVE_PSEUG && defined WITH_PSEUG
1876         return Ftruncate(wrong_type_argument(Qcomparablep, number));
1877 #else  /* !HAVE_MPC */
1878         return Ftruncate(wrong_type_argument(Qnumberp, number));
1879 #endif  /* HAVE_MPC */
1880 }
1881
1882 DEFUN("almost=", Falmost_eq, 2, 3, 0,   /*
1883 Return t if NUMBER1 is almost equal to NUMBER2.
1884
1885 Optional argument THRES can be used to specify the threshold,
1886 float-epsilon by default.
1887 */
1888       (number1, number2, thres))
1889 {
1890 #if defined HAVE_FPFLOAT
1891         if (NILP(thres)) {
1892                 thres = Vfloat_epsilon;
1893         }
1894         CHECK_FLOAT(thres);
1895
1896         if (FLOATP(number1) && FLOATP(number2)) {
1897                 fpfloat n1 = XFLOAT_DATA(number1);
1898                 fpfloat n2 = XFLOAT_DATA(number2);
1899                 fpfloat thr = XFLOAT_DATA(thres);
1900                 fpfloat d;
1901                 if (n1 >= n2) {
1902                         d = n1 - n2;
1903                 } else {
1904                         d = n2 - n1;
1905                 }
1906                 return d < thr ? Qt : Qnil;
1907         }
1908 #endif  /* HAVE_FPFLOAT */
1909         return ent_binrel(ASE_BINARY_REL_EQUALP, number1, number2) ? Qt : Qnil;
1910 }
1911
1912 DEFUN("almost/=", Falmost_neq, 2, 3, 0, /*
1913 Return t if NUMBER1 is clearly different from NUMBER2.
1914
1915 Optional argument THRES can be used to specify the threshold,
1916 float-epsilon by default.
1917 */
1918       (number1, number2, thres))
1919 {
1920 #if defined HAVE_FPFLOAT
1921         if (NILP(thres)) {
1922                 thres = Vfloat_epsilon;
1923         }
1924         CHECK_FLOAT(thres);
1925
1926         if (FLOATP(number1) && FLOATP(number2)) {
1927                 fpfloat n1 = XFLOAT_DATA(number1);
1928                 fpfloat n2 = XFLOAT_DATA(number2);
1929                 fpfloat thr = XFLOAT_DATA(thres);
1930                 fpfloat d;
1931                 if (n1 >= n2) {
1932                         d = n1 - n2;
1933                 } else {
1934                         d = n2 - n1;
1935                 }
1936                 return d < thr ? Qnil : Qt;
1937         }
1938 #endif  /* HAVE_FPFLOAT */
1939         return ent_binrel(ASE_BINARY_REL_NEQP, number1, number2) ? Qt : Qnil;
1940 }
1941
1942 \f
1943 /* misc complex functions */
1944 DEFUN("conjugate", Fconjugate, 1, 1, 0, /*
1945 Return the \(canonical\) conjugate of NUMBER.
1946 If NUMBER is a comparable, just return NUMBER.
1947 */
1948       (number))
1949 {
1950         if (COMPARABLEP(number)) {
1951                 return number;
1952 #if defined HAVE_PSEUG && defined WITH_PSEUG
1953         } else if (BIGGP(number)) {
1954                 bigg_conj(ent_scratch_bigg, XBIGG_DATA(number));
1955                 return make_bigg_bg(ent_scratch_bigg);
1956 #endif
1957 #if defined HAVE_MPC && defined WITH_MPC ||     \
1958         defined HAVE_PSEUC && defined WITH_PSEUC
1959         } else if (BIGCP(number)) {
1960                 bigc_set_prec(ent_scratch_bigc, XBIGC_GET_PREC(number));
1961                 bigc_conj(ent_scratch_bigc, XBIGC_DATA(number));
1962                 return make_bigc_bc(ent_scratch_bigc);
1963 #endif
1964 #if defined HAVE_QUATERN && defined WITH_QUATERN
1965         } else if (QUATERNP(number)) {
1966                 quatern_conj(ent_scratch_quatern, XQUATERN_DATA(number));
1967                 return make_quatern_qu(ent_scratch_quatern);
1968 #endif
1969         } else if (INDEFP(number)) {
1970                 return number;
1971         }
1972
1973         /* what should the rest do? */
1974         return Fconjugate(wrong_type_argument(Qnumberp, number));
1975 }
1976
1977 DEFUN("canonical-norm", Fcanonical_norm, 1, 1, 0,       /*
1978 Return the canonical norm of NUMBER.
1979 */
1980       (number))
1981 {
1982         if (INDEFP(number)) {
1983                 if (INFINITYP(number))
1984                         return make_indef(POS_INFINITY);
1985                 else
1986                         return make_indef(NOT_A_NUMBER);
1987         } else if (COMPARABLEP(number)) {
1988                 return Fabs(number);
1989 #if defined HAVE_PSEUG && defined WITH_PSEUG
1990         } else if (BIGGP(number)) {
1991                 bigg_norm(ent_scratch_bigz, XBIGG_DATA(number));
1992                 return make_bigz_bz(ent_scratch_bigz);
1993 #endif
1994 #if defined HAVE_MPC && defined WITH_MPC ||     \
1995         defined HAVE_PSEUC && defined WITH_PSEUC
1996         } else if (BIGCP(number)) {
1997                 bigfr_set_prec(ent_scratch_bigfr, XBIGC_GET_PREC(number));
1998                 bigc_norm(ent_scratch_bigfr, XBIGC_DATA(number));
1999                 return make_bigfr_bfr(ent_scratch_bigfr);
2000 #endif
2001 #if defined HAVE_QUATERN && defined WITH_QUATERN
2002         } else if (QUATERNP(number)) {
2003                 quatern_norm(ent_scratch_bigz, XQUATERN_DATA(number));
2004                 return make_bigz_bz(ent_scratch_bigz);
2005 #endif
2006         }
2007
2008         /* what should the rest do? */
2009         return Fcanonical_norm(wrong_type_argument(Qnumberp, number));
2010 }
2011
2012 DEFUN("real-part", Freal_part, 1, 1, 0, /*
2013 Return the real part of NUMBER.
2014 */
2015       (number))
2016 {
2017         if (INDEFP(number)) {
2018                 if (COMPARABLE_INDEF_P(number))
2019                         return number;
2020                 else if (INFINITYP(number))
2021                         return make_indef(POS_INFINITY);
2022                 else
2023                         return make_indef(NOT_A_NUMBER);
2024         } else if (COMPARABLEP(number)) {
2025                 return number;
2026 #if defined HAVE_PSEUG && defined WITH_PSEUG
2027         } else if (BIGGP(number)) {
2028                 return make_bigz_bz(bigg_re(XBIGG_DATA(number)));
2029 #endif
2030 #if defined HAVE_MPC && defined WITH_MPC ||     \
2031         defined HAVE_PSEUC && defined WITH_PSEUC
2032         } else if (BIGCP(number)) {
2033                 return make_bigfr_bfr(bigc_re(XBIGC_DATA(number)));
2034 #endif
2035         }
2036
2037         /* what should the rest do? */
2038         return Freal_part(wrong_type_argument(Qnumberp, number));
2039 }
2040
2041 DEFUN("imaginary-part", Fimaginary_part, 1, 1, 0,       /*
2042 Return the imaginary part of NUMBER.
2043 If NUMBER is a comparable, 0 is returned.
2044 */
2045       (number))
2046 {
2047         if (INDEFP(number)) {
2048                 if (COMPARABLE_INDEF_P(number))
2049                         return Qzero;
2050                 else if (INFINITYP(number))
2051                         return make_indef(POS_INFINITY);
2052                 else
2053                         return make_indef(NOT_A_NUMBER);
2054         } else if (RATIONALP(number)) {
2055                 return make_int(0);
2056 #if defined HAVE_MPFR && defined WITH_MPFR
2057         } else if (REALP(number)) {
2058                 return make_bigfr(0.0, 0UL);
2059 #endif
2060 #if defined HAVE_PSEUG && defined WITH_PSEUG
2061         } else if (BIGGP(number)) {
2062                 return make_bigz_bz(bigg_im(XBIGG_DATA(number)));
2063 #endif
2064 #if defined HAVE_MPC && defined WITH_MPC ||     \
2065         defined HAVE_PSEUC && defined WITH_PSEUC
2066         } else if (BIGCP(number)) {
2067                 return make_bigfr_bfr(bigc_im(XBIGC_DATA(number)));
2068 #endif
2069         }
2070
2071         /* what should the rest do? */
2072         return Fimaginary_part(wrong_type_argument(Qnumberp, number));
2073 }
2074
2075 \f
2076 /* Float-rounding functions. */
2077 #if defined(HAVE_FPFLOAT) || defined(HAVE_MPFR) && defined WITH_MPFR || \
2078         defined(HAVE_MPF) && defined WITH_GMP
2079
2080 DEFUN("fceiling", Ffceiling, 1, 1, 0,   /*
2081 Return the smallest integer no less than NUMBER, as a float.
2082 \(Round toward +inf.\)
2083 */
2084       (number))
2085 {
2086 #if defined HAVE_MPF && defined WITH_GMP
2087         if (BIGFP(number)) {
2088                 bigf_set_prec(ent_scratch_bigf, XBIGF_GET_PREC(number));
2089
2090                 bigf_ceil(ent_scratch_bigf, XBIGF_DATA(number));
2091                 return make_bigf_bf(ent_scratch_bigf);
2092         }
2093 #endif  /* HAVE_MPF */
2094
2095 #if defined HAVE_MPFR && defined WITH_MPFR
2096         if (BIGFRP(number)) {
2097                 bigfr_set_prec(ent_scratch_bigfr, XBIGFR_GET_PREC(number));
2098
2099                 bigfr_ceil(ent_scratch_bigfr, XBIGFR_DATA(number));
2100                 return make_bigfr_bfr(ent_scratch_bigfr);
2101         }
2102 #endif  /* HAVE_MPFR */
2103
2104         if (INDEFP(number))
2105                 return number;
2106
2107         number = ent_lift(number, FLOAT_T, NULL);
2108
2109         if (FLOATP(number))
2110                 return make_float(ceil(XFLOAT_DATA(number)));
2111         else
2112                 return number;
2113 }
2114
2115 DEFUN("ffloor", Fffloor, 1, 1, 0,       /*
2116 Return the largest integer no greater than NUMBER, as a float.
2117 \(Round towards -inf.\)
2118 */
2119       (number))
2120 {
2121 #if defined HAVE_MPF && defined WITH_GMP
2122         if (BIGFP(number)) {
2123                 bigf_set_prec(ent_scratch_bigf, XBIGF_GET_PREC(number));
2124
2125                 bigf_floor(ent_scratch_bigf, XBIGF_DATA(number));
2126                 return make_bigf_bf(ent_scratch_bigf);
2127         }
2128 #endif  /* HAVE_MPF */
2129
2130 #if defined HAVE_MPFR && defined WITH_MPFR
2131         if (BIGFRP(number)) {
2132                 bigfr_set_prec(ent_scratch_bigfr, XBIGFR_GET_PREC(number));
2133
2134                 bigfr_floor(ent_scratch_bigfr, XBIGFR_DATA(number));
2135                 return make_bigfr_bfr(ent_scratch_bigfr);
2136         }
2137 #endif  /* HAVE_MPFR */
2138
2139         if (INDEFP(number))
2140                 return number;
2141
2142         number = ent_lift(number, FLOAT_T, NULL);
2143
2144         if (FLOATP(number))
2145                 return make_float(floor(XFLOAT_DATA(number)));
2146         else
2147                 return number;
2148 }
2149
2150 DEFUN("fround", Ffround, 1, 1, 0,       /*
2151 Return the nearest integer to NUMBER, as a float.
2152 */
2153       (number))
2154 {
2155 #if defined HAVE_MPF && defined WITH_GMP
2156         if (BIGFP(number)) {
2157                 warn_when_safe(Qbigf, Qnotice,
2158                                "rounding number of type 'bigf (mpf-floats)"
2159                                "not yet implemented");
2160                 return number;
2161         }
2162 #endif  /* HAVE_MPF */
2163
2164 #if defined HAVE_MPFR && defined WITH_MPFR
2165         if (BIGFRP(number)) {
2166                 bigfr_set_prec(ent_scratch_bigfr, XBIGFR_GET_PREC(number));
2167
2168                 bigfr_rint(ent_scratch_bigfr, XBIGFR_DATA(number));
2169                 return make_bigfr_bfr(ent_scratch_bigfr);
2170         }
2171 #endif  /* HAVE_MPFR */
2172
2173         if (INDEFP(number))
2174                 return number;
2175
2176         number = ent_lift(number, FLOAT_T, NULL);
2177
2178         if (FLOATP(number))
2179                 return make_float(emacs_rint(XFLOAT_DATA(number)));
2180         else
2181                 return number;
2182 }
2183
2184 DEFUN("ftruncate", Fftruncate, 1, 1, 0, /*
2185 Truncate a floating point number to an integral float value.
2186 Rounds the value toward zero.
2187 */
2188       (number))
2189 {
2190         fpfloat d;
2191 #if defined HAVE_MPF && defined WITH_GMP
2192         if (BIGFP(number)) {
2193                 bigf_set_prec(ent_scratch_bigf, XBIGF_GET_PREC(number));
2194
2195                 bigf_trunc(ent_scratch_bigf, XBIGF_DATA(number));
2196                 return make_bigf_bf(ent_scratch_bigf);
2197         }
2198 #endif  /* HAVE_MPF */
2199
2200 #if defined HAVE_MPFR && defined WITH_MPFR
2201         if (BIGFRP(number)) {
2202                 bigfr_set_prec(ent_scratch_bigfr, XBIGFR_GET_PREC(number));
2203
2204                 bigfr_trunc(ent_scratch_bigfr, XBIGFR_DATA(number));
2205                 return make_bigfr_bfr(ent_scratch_bigfr);
2206         }
2207 #endif  /* HAVE_MPFR */
2208
2209         if (INDEFP(number))
2210                 return number;
2211
2212         number = ent_lift(number, FLOAT_T, NULL);
2213
2214         if (FLOATP(number)) {
2215                 d = XFLOAT_DATA(number);
2216                 if (d >= 0.0)
2217                         d = floor(d);
2218                 else
2219                         d = ceil(d);
2220                 return make_float(d);
2221         } else {
2222                 return number;
2223         }
2224 }
2225 #endif  /* HAVE_MPF(R) || HAVE_FPFLOAT (float-rounding functions) */
2226
2227 \f
2228 #ifdef HAVE_FPFLOAT
2229 #ifdef FLOAT_CATCH_SIGILL
2230 static SIGTYPE float_error(int signo)
2231 {
2232         if (!in_float)
2233                 fatal_error_signal(signo);
2234
2235         EMACS_REESTABLISH_SIGNAL(signo, arith_error);
2236         EMACS_UNBLOCK_SIGNAL(signo);
2237
2238         in_float = 0;
2239
2240         /* Was Fsignal(), but it just doesn't make sense for an error
2241            occurring inside a signal handler to be restartable, considering
2242            that anything could happen when the error is signaled and trapped
2243            and considering the asynchronous nature of signal handlers. */
2244         signal_error(Qarith_error, list1(float_error_arg));
2245 }
2246
2247 /* Another idea was to replace the library function `infnan'
2248    where SIGILL is signaled.  */
2249
2250 #endif                          /* FLOAT_CATCH_SIGILL */
2251
2252 /* In C++, it is impossible to determine what type matherr expects
2253    without some more configure magic.
2254    We shouldn't be using matherr anyways - it's a non-standard SYSVism. */
2255 #if defined (HAVE_MATHERR) && !defined(__cplusplus)
2256 int matherr(struct exception *x)
2257 {
2258         Lisp_Object args;
2259         if (!in_float)
2260                 /* Not called from emacs-lisp float routines; do the default thing. */
2261                 return 0;
2262
2263         /* if (!strcmp (x->name, "pow")) x->name = "expt"; */
2264
2265         args = Fcons(build_string(x->name),
2266                      Fcons(make_float(x->arg1), ((in_float == 2)
2267                                                  ? Fcons(make_float(x->arg2),
2268                                                          Qnil)
2269                                                  : Qnil)));
2270         switch (x->type) {
2271         case DOMAIN:
2272                 Fsignal(Qdomain_error, args);
2273                 break;
2274         case SING:
2275                 Fsignal(Qsingularity_error, args);
2276                 break;
2277         case OVERFLOW:
2278                 Fsignal(Qoverflow_error, args);
2279                 break;
2280         case UNDERFLOW:
2281                 Fsignal(Qunderflow_error, args);
2282                 break;
2283         default:
2284                 Fsignal(Qarith_error, args);
2285                 break;
2286         }
2287         return 1;               /* don't set errno or print a message */
2288 }
2289 #endif                          /* HAVE_MATHERR */
2290 #endif                          /* HAVE_FPFLOAT */
2291 \f
2292 void init_floatfns_very_early(void)
2293 {
2294 #ifdef HAVE_FPFLOAT
2295 # ifdef FLOAT_CATCH_SIGILL
2296         signal(SIGILL, float_error);
2297 # endif
2298         in_float = 0;
2299 #endif                          /* HAVE_FPFLOAT */
2300 }
2301
2302 void syms_of_floatfns(void)
2303 {
2304
2305         /* Trig functions.  */
2306
2307 #if defined(HAVE_FPFLOAT) || defined HAVE_MPFR && defined WITH_MPFR
2308         DEFSUBR(Facos);
2309         DEFSUBR(Fasin);
2310         DEFSUBR(Fatan);
2311         DEFSUBR(Fcos);
2312         DEFSUBR(Fsin);
2313         DEFSUBR(Ftan);
2314 #endif  /* HAVE_FPFLOAT || HAVE_MPFR*/
2315 #if defined HAVE_MPFR && defined WITH_MPFR
2316         DEFSUBR(Fsec);
2317         DEFSUBR(Fcsc);
2318         DEFSUBR(Fcot);
2319 #endif
2320
2321         /* Bessel functions */
2322
2323 #if 0
2324         DEFSUBR(Fbessel_y0);
2325         DEFSUBR(Fbessel_y1);
2326         DEFSUBR(Fbessel_yn);
2327         DEFSUBR(Fbessel_j0);
2328         DEFSUBR(Fbessel_j1);
2329         DEFSUBR(Fbessel_jn);
2330 #endif                          /* 0 */
2331
2332         /* Error functions. */
2333
2334 #if 1
2335 #if defined(HAVE_FPFLOAT) || defined(HAVE_MPFR) && defined WITH_MPFR
2336         DEFSUBR(Ferf);
2337         DEFSUBR(Ferfc);
2338         DEFSUBR(Flog_gamma);
2339 #endif
2340 #endif                          /* 0 */
2341
2342         /* Root and Log functions. */
2343
2344 #if defined(HAVE_FPFLOAT) || defined(HAVE_MPFR) && defined WITH_MPFR
2345         DEFSUBR(Fexp);
2346 #endif  /* HAVE_FPFLOAT || HAVE_MPFR */
2347         DEFSUBR(Fexp2);
2348         DEFSUBR(Fexp10);
2349 #if 0
2350         DEFSUBR(Fexpt);
2351 #endif
2352 #if defined(HAVE_FPFLOAT) || defined(HAVE_MPFR) && defined WITH_MPFR
2353         DEFSUBR(Flog);
2354         DEFSUBR(Flog2);
2355         DEFSUBR(Flog10);
2356         DEFSUBR(Fsqrt);
2357         DEFSUBR(Fcube_root);
2358 #if defined HAVE_MPFR && defined WITH_MPFR
2359         DEFSUBR(Froot);
2360 #endif
2361 #endif  /* HAVE_FPFLOAT || HAVE_MPFR*/
2362
2363         /* Inverse trig functions. */
2364
2365 #if defined(HAVE_FPFLOAT) || defined(HAVE_MPFR) && defined WITH_MPFR
2366         DEFSUBR(Facosh);
2367         DEFSUBR(Fasinh);
2368         DEFSUBR(Fatanh);
2369         DEFSUBR(Fcosh);
2370         DEFSUBR(Fsinh);
2371         DEFSUBR(Ftanh);
2372 #endif  /* HAVE_FPFLOAT || HAVE_MPFR */
2373 #if defined HAVE_MPFR && defined WITH_MPFR
2374         DEFSUBR(Fsech);
2375         DEFSUBR(Fcsch);
2376         DEFSUBR(Fcoth);
2377 #endif  /* HAVE_MPFR */
2378
2379         /* Rounding functions */
2380
2381         DEFSUBR(Fabs);
2382 #ifdef HAVE_FPFLOAT
2383         DEFSUBR(Ffloat);
2384         DEFSUBR(Flogb);
2385 #endif                          /* HAVE_FPFLOAT */
2386         DEFSUBR(Fceiling);
2387         DEFSUBR(Ffloor);
2388         DEFSUBR(Fround);
2389         DEFSUBR(Ftruncate);
2390         DEFSUBR(Falmost_eq);
2391         DEFSUBR(Falmost_neq);
2392
2393         /* misc complex functions */
2394         DEFSUBR(Fconjugate);
2395         DEFSUBR(Fcanonical_norm);
2396         DEFSUBR(Freal_part);
2397         DEFSUBR(Fimaginary_part);
2398
2399         /* Float-rounding functions. */
2400
2401 #if defined(HAVE_FPFLOAT) || defined(HAVE_MPF) && defined WITH_GMP ||   \
2402         defined(HAVE_MPFR) && defined WITH_MPFR
2403         DEFSUBR(Ffceiling);
2404         DEFSUBR(Fffloor);
2405         DEFSUBR(Ffround);
2406         DEFSUBR(Fftruncate);
2407 #endif  /* HAVE_FPFLOAT || HAVE_MPF(R) */
2408 }
2409
2410 void vars_of_floatfns(void)
2411 {
2412 }