src/floatfns.c

   1 /* Primitive operations on floating point for GNU Emacs Lisp interpreter.
   2    Copyright (C) 1988, 1993 Free Software Foundation, Inc.
   3
   4 This file is part of GNU Emacs.
   5
   6 GNU Emacs 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 2, or (at your option)
   9 any later version.
  10
  11 GNU Emacs 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 GNU Emacs; see the file COPYING.  If not, write to
  18 the Free Software Foundation, 675 Mass Ave, Cambridge, MA 02139, USA.  */
  19
  20
  21 /* ANSI C requires only these float functions:
  22    acos, asin, atan, atan2, ceil, cos, cosh, exp, fabs, floor, fmod,
  23    frexp, ldexp, log, log10, modf, pow, sin, sinh, sqrt, tan, tanh.
  24
  25    Define HAVE_INVERSE_HYPERBOLIC if you have acosh, asinh, and atanh.
  26    Define HAVE_CBRT if you have cbrt.
  27    Define HAVE_RINT if you have rint.
  28    If you don't define these, then the appropriate routines will be simulated.
  29
  30    Define HAVE_MATHERR if on a system supporting the SysV matherr callback.
  31    (This should happen automatically.)
  32
  33    Define FLOAT_CHECK_ERRNO if the float library routines set errno.
  34    This has no effect if HAVE_MATHERR is defined.
  35
  36    Define FLOAT_CATCH_SIGILL if the float library routines signal SIGILL.
  37    (What systems actually do this?  Please let us know.)
  38
  39    Define FLOAT_CHECK_DOMAIN if the float library doesn't handle errors by
  40    either setting errno, or signalling SIGFPE/SIGILL.  Otherwise, domain and
  41    range checking will happen before calling the float routines.  This has
  42    no effect if HAVE_MATHERR is defined (since matherr will be called when
  43    a domain error occurs.)
  44  */
  45
  46 #include <signal.h>
  47
  48 #include <config.h>
  49 #include "lisp.h"
  50 #include "syssignal.h"
  51
  52 Lisp_Object Qarith_error;
  53
  54 #ifdef LISP_FLOAT_TYPE
  55
  56 #ifdef MSDOS
  57 /* These are redefined (correctly, but differently) in values.h.  */
  58 #undef INTBITS
  59 #undef LONGBITS
  60 #undef SHORTBITS
  61 #endif
  62
  63 /* Work around a problem that happens because math.h on hpux 7
  64    defines two static variables--which, in Emacs, are not really static,
  65    because `static' is defined as nothing.  The problem is that they are
  66    defined both here and in lread.c.
  67    These macros prevent the name conflict.  */
  68 #if defined (HPUX) && !defined (HPUX8)
  69 #define _MAXLDBL floatfns_maxldbl
  70 #define _NMAXLDBL floatfns_nmaxldbl
  71 #endif
  72
  73 #include <math.h>
  74
  75 /* This declaration is omitted on some systems, like Ultrix.  */
  76 #if !defined (HPUX) && defined (HAVE_LOGB)
  77 extern double logb ();
  78 #endif /* not HPUX and HAVE_LOGB */
  79
  80 #ifndef MSDOS
  81 #if defined(DOMAIN) && defined(SING) && defined(OVERFLOW)
  82     /* If those are defined, then this is probably a `matherr' machine. */
  83 # ifndef HAVE_MATHERR
  84 #  define HAVE_MATHERR
  85 # endif
  86 #endif
  87 #endif
  88
  89 #ifdef NO_MATHERR
  90 #undef HAVE_MATHERR
  91 #endif
  92
  93 #ifdef HAVE_MATHERR
  94 # ifdef FLOAT_CHECK_ERRNO
  95 #  undef FLOAT_CHECK_ERRNO
  96 # endif
  97 # ifdef FLOAT_CHECK_DOMAIN
  98 #  undef FLOAT_CHECK_DOMAIN
  99 # endif
 100 #endif
 101
 102 #ifndef NO_FLOAT_CHECK_ERRNO
 103 #define FLOAT_CHECK_ERRNO
 104 #endif
 105
 106 #ifdef FLOAT_CHECK_ERRNO
 107 # include <errno.h>
 108
 109 extern int errno;
 110 #endif
 111
 112 /* Avoid traps on VMS from sinh and cosh.
 113    All the other functions set errno instead.  */
 114
 115 #ifdef VMS
 116 #undef cosh
 117 #undef sinh
 118 #define cosh(x) ((exp(x)+exp(-x))*0.5)
 119 #define sinh(x) ((exp(x)-exp(-x))*0.5)
 120 #endif /* VMS */
 121
 122 #ifndef HAVE_RINT
 123 #define rint(x) (floor((x)+0.5))
 124 #endif
 125
 126 static SIGTYPE float_error ();
 127
 128 /* Nonzero while executing in floating point.
 129    This tells float_error what to do.  */
 130
 131 static int in_float;
 132
 133 /* If an argument is out of range for a mathematical function,
 134    here is the actual argument value to use in the error message.  */
 135
 136 static Lisp_Object float_error_arg, float_error_arg2;
 137
 138 static char *float_error_fn_name;
 139
 140 /* Evaluate the floating point expression D, recording NUM
 141    as the original argument for error messages.
 142    D is normally an assignment expression.
 143    Handle errors which may result in signals or may set errno.
 144
 145    Note that float_error may be declared to return void, so you can't
 146    just cast the zero after the colon to (SIGTYPE) to make the types
 147    check properly.  */
 148
 149 #ifdef FLOAT_CHECK_ERRNO
 150 #define IN_FLOAT(d, name, num)                          \
 151   do {                                                  \
 152     float_error_arg = num;                              \
 153     float_error_fn_name = name;                         \
 154     in_float = 1; errno = 0; (d); in_float = 0;         \
 155     switch (errno) {                                    \
 156     case 0: break;                                      \
 157     case EDOM:   domain_error (float_error_fn_name, float_error_arg);   \
 158     case ERANGE: range_error (float_error_fn_name, float_error_arg);    \
 159     default:     arith_error (float_error_fn_name, float_error_arg);    \
 160     }                                                   \
 161   } while (0)
 162 #define IN_FLOAT2(d, name, num, num2)                   \
 163   do {                                                  \
 164     float_error_arg = num;                              \
 165     float_error_arg2 = num2;                            \
 166     float_error_fn_name = name;                         \
 167     in_float = 1; errno = 0; (d); in_float = 0;         \
 168     switch (errno) {                                    \
 169     case 0: break;                                      \
 170     case EDOM:   domain_error (float_error_fn_name, float_error_arg);   \
 171     case ERANGE: range_error (float_error_fn_name, float_error_arg);    \
 172     default:     arith_error (float_error_fn_name, float_error_arg);    \
 173     }                                                   \
 174   } while (0)
 175 #else
 176 #define IN_FLOAT(d, name, num) (in_float = 1, (d), in_float = 0)
 177 #define IN_FLOAT2(d, name, num, num2) (in_float = 1, (d), in_float = 0)
 178 #endif
 179
 180 /* Convert float to Lisp_Int if it fits, else signal a range error
 181    using the given arguments.  */
 182 #define FLOAT_TO_INT(x, i, name, num)                                   \
 183   do                                                                    \
 184     {                                                                   \
 185       if ((x) >= (1 << (VALBITS-1)) || (x) <= - (1 << (VALBITS-1)) - 1) \
 186         range_error (name, num);                                        \
 187       XSET (i, Lisp_Int,  (int)(x));                                    \
 188     }                                                                   \
 189   while (0)
 190 #define FLOAT_TO_INT2(x, i, name, num1, num2)                           \
 191   do                                                                    \
 192     {                                                                   \
 193       if ((x) >= (1 << (VALBITS-1)) || (x) <= - (1 << (VALBITS-1)) - 1) \
 194         range_error2 (name, num1, num2);                                \
 195       XSET (i, Lisp_Int,  (int)(x));                                    \
 196     }                                                                   \
 197   while (0)
 198
 199 #define arith_error(op,arg) \
 200   Fsignal (Qarith_error, Fcons (build_string ((op)), Fcons ((arg), Qnil)))
 201 #define range_error(op,arg) \
 202   Fsignal (Qrange_error, Fcons (build_string ((op)), Fcons ((arg), Qnil)))
 203 #define range_error2(op,a1,a2) \
 204   Fsignal (Qrange_error, Fcons (build_string ((op)), \
 205                                 Fcons ((a1), Fcons ((a2), Qnil))))
 206 #define domain_error(op,arg) \
 207   Fsignal (Qdomain_error, Fcons (build_string ((op)), Fcons ((arg), Qnil)))
 208 #define domain_error2(op,a1,a2) \
 209   Fsignal (Qdomain_error, Fcons (build_string ((op)), \
 210                                  Fcons ((a1), Fcons ((a2), Qnil))))
 211
 212 /* Extract a Lisp number as a `double', or signal an error.  */
 213
 214 double
 215 extract_float (num)
 216      Lisp_Object num;
 217 {
 218   CHECK_NUMBER_OR_FLOAT (num, 0);
 219
 220   if (XTYPE (num) == Lisp_Float)
 221     return XFLOAT (num)->data;
 222   return (double) XINT (num);
 223 }
 224 \f
 225 /* Trig functions.  */
 226
 227 DEFUN ("acos", Facos, Sacos, 1, 1, 0,
 228   "Return the inverse cosine of ARG.")
 229   (arg)
 230      register Lisp_Object arg;
 231 {
 232   double d = extract_float (arg);
 233 #ifdef FLOAT_CHECK_DOMAIN
 234   if (d > 1.0 || d < -1.0)
 235     domain_error ("acos", arg);
 236 #endif
 237   IN_FLOAT (d = acos (d), "acos", arg);
 238   return make_float (d);
 239 }
 240
 241 DEFUN ("asin", Fasin, Sasin, 1, 1, 0,
 242   "Return the inverse sine of ARG.")
 243   (arg)
 244      register Lisp_Object arg;
 245 {
 246   double d = extract_float (arg);
 247 #ifdef FLOAT_CHECK_DOMAIN
 248   if (d > 1.0 || d < -1.0)
 249     domain_error ("asin", arg);
 250 #endif
 251   IN_FLOAT (d = asin (d), "asin", arg);
 252   return make_float (d);
 253 }
 254
 255 DEFUN ("atan", Fatan, Satan, 1, 1, 0,
 256   "Return the inverse tangent of ARG.")
 257   (arg)
 258      register Lisp_Object arg;
 259 {
 260   double d = extract_float (arg);
 261   IN_FLOAT (d = atan (d), "atan", arg);
 262   return make_float (d);
 263 }
 264
 265 DEFUN ("cos", Fcos, Scos, 1, 1, 0,
 266   "Return the cosine of ARG.")
 267   (arg)
 268      register Lisp_Object arg;
 269 {
 270   double d = extract_float (arg);
 271   IN_FLOAT (d = cos (d), "cos", arg);
 272   return make_float (d);
 273 }
 274
 275 DEFUN ("sin", Fsin, Ssin, 1, 1, 0,
 276   "Return the sine of ARG.")
 277   (arg)
 278      register Lisp_Object arg;
 279 {
 280   double d = extract_float (arg);
 281   IN_FLOAT (d = sin (d), "sin", arg);
 282   return make_float (d);
 283 }
 284
 285 DEFUN ("tan", Ftan, Stan, 1, 1, 0,
 286   "Return the tangent of ARG.")
 287   (arg)
 288      register Lisp_Object arg;
 289 {
 290   double d = extract_float (arg);
 291   double c = cos (d);
 292 #ifdef FLOAT_CHECK_DOMAIN
 293   if (c == 0.0)
 294     domain_error ("tan", arg);
 295 #endif
 296   IN_FLOAT (d = sin (d) / c, "tan", arg);
 297   return make_float (d);
 298 }
 299 \f
 300 #if 0 /* Leave these out unless we find there's a reason for them.  */
 301
 302 DEFUN ("bessel-j0", Fbessel_j0, Sbessel_j0, 1, 1, 0,
 303   "Return the bessel function j0 of ARG.")
 304   (arg)
 305      register Lisp_Object arg;
 306 {
 307   double d = extract_float (arg);
 308   IN_FLOAT (d = j0 (d), "bessel-j0", arg);
 309   return make_float (d);
 310 }
 311
 312 DEFUN ("bessel-j1", Fbessel_j1, Sbessel_j1, 1, 1, 0,
 313   "Return the bessel function j1 of ARG.")
 314   (arg)
 315      register Lisp_Object arg;
 316 {
 317   double d = extract_float (arg);
 318   IN_FLOAT (d = j1 (d), "bessel-j1", arg);
 319   return make_float (d);
 320 }
 321
 322 DEFUN ("bessel-jn", Fbessel_jn, Sbessel_jn, 2, 2, 0,
 323   "Return the order N bessel function output jn of ARG.\n\
 324 The first arg (the order) is truncated to an integer.")
 325   (arg1, arg2)
 326      register Lisp_Object arg1, arg2;
 327 {
 328   int i1 = extract_float (arg1);
 329   double f2 = extract_float (arg2);
 330
 331   IN_FLOAT (f2 = jn (i1, f2), "bessel-jn", arg1);
 332   return make_float (f2);
 333 }
 334
 335 DEFUN ("bessel-y0", Fbessel_y0, Sbessel_y0, 1, 1, 0,
 336   "Return the bessel function y0 of ARG.")
 337   (arg)
 338      register Lisp_Object arg;
 339 {
 340   double d = extract_float (arg);
 341   IN_FLOAT (d = y0 (d), "bessel-y0", arg);
 342   return make_float (d);
 343 }
 344
 345 DEFUN ("bessel-y1", Fbessel_y1, Sbessel_y1, 1, 1, 0,
 346   "Return the bessel function y1 of ARG.")
 347   (arg)
 348      register Lisp_Object arg;
 349 {
 350   double d = extract_float (arg);
 351   IN_FLOAT (d = y1 (d), "bessel-y0", arg);
 352   return make_float (d);
 353 }
 354
 355 DEFUN ("bessel-yn", Fbessel_yn, Sbessel_yn, 2, 2, 0,
 356   "Return the order N bessel function output yn of ARG.\n\
 357 The first arg (the order) is truncated to an integer.")
 358   (arg1, arg2)
 359      register Lisp_Object arg1, arg2;
 360 {
 361   int i1 = extract_float (arg1);
 362   double f2 = extract_float (arg2);
 363
 364   IN_FLOAT (f2 = yn (i1, f2), "bessel-yn", arg1);
 365   return make_float (f2);
 366 }
 367
 368 #endif
 369 \f
 370 #if 0 /* Leave these out unless we see they are worth having.  */
 371
 372 DEFUN ("erf", Ferf, Serf, 1, 1, 0,
 373   "Return the mathematical error function of ARG.")
 374   (arg)
 375      register Lisp_Object arg;
 376 {
 377   double d = extract_float (arg);
 378   IN_FLOAT (d = erf (d), "erf", arg);
 379   return make_float (d);
 380 }
 381
 382 DEFUN ("erfc", Ferfc, Serfc, 1, 1, 0,
 383   "Return the complementary error function of ARG.")
 384   (arg)
 385      register Lisp_Object arg;
 386 {
 387   double d = extract_float (arg);
 388   IN_FLOAT (d = erfc (d), "erfc", arg);
 389   return make_float (d);
 390 }
 391
 392 DEFUN ("log-gamma", Flog_gamma, Slog_gamma, 1, 1, 0,
 393   "Return the log gamma of ARG.")
 394   (arg)
 395      register Lisp_Object arg;
 396 {
 397   double d = extract_float (arg);
 398   IN_FLOAT (d = lgamma (d), "log-gamma", arg);
 399   return make_float (d);
 400 }
 401
 402 DEFUN ("cube-root", Fcube_root, Scube_root, 1, 1, 0,
 403   "Return the cube root of ARG.")
 404   (arg)
 405      register Lisp_Object arg;
 406 {
 407   double d = extract_float (arg);
 408 #ifdef HAVE_CBRT
 409   IN_FLOAT (d = cbrt (d), "cube-root", arg);
 410 #else
 411   if (d >= 0.0)
 412     IN_FLOAT (d = pow (d, 1.0/3.0), "cube-root", arg);
 413   else
 414     IN_FLOAT (d = -pow (-d, 1.0/3.0), "cube-root", arg);
 415 #endif
 416   return make_float (d);
 417 }
 418
 419 #endif
 420 \f
 421 DEFUN ("exp", Fexp, Sexp, 1, 1, 0,
 422   "Return the exponential base e of ARG.")
 423   (arg)
 424      register Lisp_Object arg;
 425 {
 426   double d = extract_float (arg);
 427 #ifdef FLOAT_CHECK_DOMAIN
 428   if (d > 709.7827)   /* Assume IEEE doubles here */
 429     range_error ("exp", arg);
 430   else if (d < -709.0)
 431     return make_float (0.0);
 432   else
 433 #endif
 434     IN_FLOAT (d = exp (d), "exp", arg);
 435   return make_float (d);
 436 }
 437
 438 DEFUN ("expt", Fexpt, Sexpt, 2, 2, 0,
 439   "Return the exponential X ** Y.")
 440   (arg1, arg2)
 441      register Lisp_Object arg1, arg2;
 442 {
 443   double f1, f2;
 444
 445   CHECK_NUMBER_OR_FLOAT (arg1, 0);
 446   CHECK_NUMBER_OR_FLOAT (arg2, 0);
 447   if (XTYPE (arg1) == Lisp_Int     /* common lisp spec */
 448       && XTYPE (arg2) == Lisp_Int) /* don't promote, if both are ints */
 449     {                           /* this can be improved by pre-calculating */
 450       int acc, x, y;            /* some binary powers of x then accumulating */
 451       Lisp_Object val;
 452
 453       x = XINT (arg1);
 454       y = XINT (arg2);
 455       acc = 1;
 456
 457       if (y < 0)
 458         {
 459           if (x == 1)
 460             acc = 1;
 461           else if (x == -1)
 462             acc = (y & 1) ? -1 : 1;
 463           else
 464             acc = 0;
 465         }
 466       else
 467         {
 468           while (y > 0)
 469             {
 470               if (y & 1)
 471                 acc *= x;
 472               x *= x;
 473               y = (unsigned)y >> 1;
 474             }
 475         }
 476       XSET (val, Lisp_Int, acc);
 477       return val;
 478     }
 479   f1 = (XTYPE (arg1) == Lisp_Float) ? XFLOAT (arg1)->data : XINT (arg1);
 480   f2 = (XTYPE (arg2) == Lisp_Float) ? XFLOAT (arg2)->data : XINT (arg2);
 481   /* Really should check for overflow, too */
 482   if (f1 == 0.0 && f2 == 0.0)
 483     f1 = 1.0;
 484 #ifdef FLOAT_CHECK_DOMAIN
 485   else if ((f1 == 0.0 && f2 < 0.0) || (f1 < 0 && f2 != floor(f2)))
 486     domain_error2 ("expt", arg1, arg2);
 487 #endif
 488   IN_FLOAT2 (f1 = pow (f1, f2), "expt", arg1, arg2);
 489   return make_float (f1);
 490 }
 491
 492 DEFUN ("log", Flog, Slog, 1, 2, 0,
 493   "Return the natural logarithm of ARG.\n\
 494 If second optional argument BASE is given, return log ARG using that base.")
 495   (arg, base)
 496      register Lisp_Object arg, base;
 497 {
 498   double d = extract_float (arg);
 499
 500 #ifdef FLOAT_CHECK_DOMAIN
 501   if (d <= 0.0)
 502     domain_error2 ("log", arg, base);
 503 #endif
 504   if (NILP (base))
 505     IN_FLOAT (d = log (d), "log", arg);
 506   else
 507     {
 508       double b = extract_float (base);
 509
 510 #ifdef FLOAT_CHECK_DOMAIN
 511       if (b <= 0.0 || b == 1.0)
 512         domain_error2 ("log", arg, base);
 513 #endif
 514       if (b == 10.0)
 515         IN_FLOAT2 (d = log10 (d), "log", arg, base);
 516       else
 517         IN_FLOAT2 (d = log (d) / log (b), "log", arg, base);
 518     }
 519   return make_float (d);
 520 }
 521
 522 DEFUN ("log10", Flog10, Slog10, 1, 1, 0,
 523   "Return the logarithm base 10 of ARG.")
 524   (arg)
 525      register Lisp_Object arg;
 526 {
 527   double d = extract_float (arg);
 528 #ifdef FLOAT_CHECK_DOMAIN
 529   if (d <= 0.0)
 530     domain_error ("log10", arg);
 531 #endif
 532   IN_FLOAT (d = log10 (d), "log10", arg);
 533   return make_float (d);
 534 }
 535
 536 DEFUN ("sqrt", Fsqrt, Ssqrt, 1, 1, 0,
 537   "Return the square root of ARG.")
 538   (arg)
 539      register Lisp_Object arg;
 540 {
 541   double d = extract_float (arg);
 542 #ifdef FLOAT_CHECK_DOMAIN
 543   if (d < 0.0)
 544     domain_error ("sqrt", arg);
 545 #endif
 546   IN_FLOAT (d = sqrt (d), "sqrt", arg);
 547   return make_float (d);
 548 }
 549 \f
 550 #if 0 /* Not clearly worth adding.  */
 551
 552 DEFUN ("acosh", Facosh, Sacosh, 1, 1, 0,
 553   "Return the inverse hyperbolic cosine of ARG.")
 554   (arg)
 555      register Lisp_Object arg;
 556 {
 557   double d = extract_float (arg);
 558 #ifdef FLOAT_CHECK_DOMAIN
 559   if (d < 1.0)
 560     domain_error ("acosh", arg);
 561 #endif
 562 #ifdef HAVE_INVERSE_HYPERBOLIC
 563   IN_FLOAT (d = acosh (d), "acosh", arg);
 564 #else
 565   IN_FLOAT (d = log (d + sqrt (d*d - 1.0)), "acosh", arg);
 566 #endif
 567   return make_float (d);
 568 }
 569
 570 DEFUN ("asinh", Fasinh, Sasinh, 1, 1, 0,
 571   "Return the inverse hyperbolic sine of ARG.")
 572   (arg)
 573      register Lisp_Object arg;
 574 {
 575   double d = extract_float (arg);
 576 #ifdef HAVE_INVERSE_HYPERBOLIC
 577   IN_FLOAT (d = asinh (d), "asinh", arg);
 578 #else
 579   IN_FLOAT (d = log (d + sqrt (d*d + 1.0)), "asinh", arg);
 580 #endif
 581   return make_float (d);
 582 }
 583
 584 DEFUN ("atanh", Fatanh, Satanh, 1, 1, 0,
 585   "Return the inverse hyperbolic tangent of ARG.")
 586   (arg)
 587      register Lisp_Object arg;
 588 {
 589   double d = extract_float (arg);
 590 #ifdef FLOAT_CHECK_DOMAIN
 591   if (d >= 1.0 || d <= -1.0)
 592     domain_error ("atanh", arg);
 593 #endif
 594 #ifdef HAVE_INVERSE_HYPERBOLIC
 595   IN_FLOAT (d = atanh (d), "atanh", arg);
 596 #else
 597   IN_FLOAT (d = 0.5 * log ((1.0 + d) / (1.0 - d)), "atanh", arg);
 598 #endif
 599   return make_float (d);
 600 }
 601
 602 DEFUN ("cosh", Fcosh, Scosh, 1, 1, 0,
 603   "Return the hyperbolic cosine of ARG.")
 604   (arg)
 605      register Lisp_Object arg;
 606 {
 607   double d = extract_float (arg);
 608 #ifdef FLOAT_CHECK_DOMAIN
 609   if (d > 710.0 || d < -710.0)
 610     range_error ("cosh", arg);
 611 #endif
 612   IN_FLOAT (d = cosh (d), "cosh", arg);
 613   return make_float (d);
 614 }
 615
 616 DEFUN ("sinh", Fsinh, Ssinh, 1, 1, 0,
 617   "Return the hyperbolic sine of ARG.")
 618   (arg)
 619      register Lisp_Object arg;
 620 {
 621   double d = extract_float (arg);
 622 #ifdef FLOAT_CHECK_DOMAIN
 623   if (d > 710.0 || d < -710.0)
 624     range_error ("sinh", arg);
 625 #endif
 626   IN_FLOAT (d = sinh (d), "sinh", arg);
 627   return make_float (d);
 628 }
 629
 630 DEFUN ("tanh", Ftanh, Stanh, 1, 1, 0,
 631   "Return the hyperbolic tangent of ARG.")
 632   (arg)
 633      register Lisp_Object arg;
 634 {
 635   double d = extract_float (arg);
 636   IN_FLOAT (d = tanh (d), "tanh", arg);
 637   return make_float (d);
 638 }
 639 #endif
 640 \f
 641 DEFUN ("abs", Fabs, Sabs, 1, 1, 0,
 642   "Return the absolute value of ARG.")
 643   (arg)
 644      register Lisp_Object arg;
 645 {
 646   CHECK_NUMBER_OR_FLOAT (arg, 0);
 647
 648   if (XTYPE (arg) == Lisp_Float)
 649     IN_FLOAT (arg = make_float (fabs (XFLOAT (arg)->data)), "abs", arg);
 650   else if (XINT (arg) < 0)
 651     XSETINT (arg, - XFASTINT (arg));
 652
 653   return arg;
 654 }
 655
 656 DEFUN ("float", Ffloat, Sfloat, 1, 1, 0,
 657   "Return the floating point number equal to ARG.")
 658   (arg)
 659      register Lisp_Object arg;
 660 {
 661   CHECK_NUMBER_OR_FLOAT (arg, 0);
 662
 663   if (XTYPE (arg) == Lisp_Int)
 664     return make_float ((double) XINT (arg));
 665   else                          /* give 'em the same float back */
 666     return arg;
 667 }
 668
 669 DEFUN ("logb", Flogb, Slogb, 1, 1, 0,
 670   "Returns largest integer <= the base 2 log of the magnitude of ARG.\n\
 671 This is the same as the exponent of a float.")
 672      (arg)
 673      Lisp_Object arg;
 674 {
 675   Lisp_Object val;
 676   int value;
 677   double f = extract_float (arg);
 678
 679   if (f == 0.0)
 680     value = -(VALMASK >> 1);
 681   else
 682     {
 683 #ifdef HAVE_LOGB
 684       IN_FLOAT (value = logb (f), "logb", arg);
 685 #else
 686 #ifdef HAVE_FREXP
 687       IN_FLOAT (frexp (f, &value), "logb", arg);
 688       value--;
 689 #else
 690       int i;
 691       double d;
 692       if (f < 0.0)
 693         f = -f;
 694       value = -1;
 695       while (f < 0.5)
 696         {
 697           for (i = 1, d = 0.5; d * d >= f; i += i)
 698             d *= d;
 699           f /= d;
 700           value -= i;
 701         }
 702       while (f >= 1.0)
 703         {
 704           for (i = 1, d = 2.0; d * d <= f; i += i)
 705             d *= d;
 706           f /= d;
 707           value += i;
 708         }
 709 #endif
 710 #endif
 711     }
 712   XSET (val, Lisp_Int, value);
 713   return val;
 714 }
 715
 716 /* the rounding functions  */
 717
 718 DEFUN ("ceiling", Fceiling, Sceiling, 1, 1, 0,
 719   "Return the smallest integer no less than ARG.  (Round toward +inf.)")
 720   (arg)
 721      register Lisp_Object arg;
 722 {
 723   CHECK_NUMBER_OR_FLOAT (arg, 0);
 724
 725   if (XTYPE (arg) == Lisp_Float)
 726     {
 727       double d;
 728
 729       IN_FLOAT (d = ceil (XFLOAT (arg)->data), "ceiling", arg);
 730       FLOAT_TO_INT (d, arg, "ceiling", arg);
 731     }
 732
 733   return arg;
 734 }
 735
 736 #endif /* LISP_FLOAT_TYPE */
 737
 738
 739 DEFUN ("floor", Ffloor, Sfloor, 1, 2, 0,
 740   "Return the largest integer no greater than ARG.  (Round towards -inf.)\n\
 741 With optional DIVISOR, return the largest integer no greater than ARG/DIVISOR.")
 742   (arg, divisor)
 743      register Lisp_Object arg, divisor;
 744 {
 745   CHECK_NUMBER_OR_FLOAT (arg, 0);
 746
 747   if (! NILP (divisor))
 748     {
 749       int i1, i2;
 750
 751       CHECK_NUMBER_OR_FLOAT (divisor, 1);
 752
 753 #ifdef LISP_FLOAT_TYPE
 754       if (XTYPE (arg) == Lisp_Float || XTYPE (divisor) == Lisp_Float)
 755         {
 756           double f1, f2;
 757
 758           f1 = XTYPE (arg) == Lisp_Float ? XFLOAT (arg)->data : XINT (arg);
 759           f2 = (XTYPE (divisor) == Lisp_Float
 760                 ? XFLOAT (divisor)->data : XINT (divisor));
 761           if (f2 == 0)
 762             Fsignal (Qarith_error, Qnil);
 763
 764           IN_FLOAT2 (f1 = floor (f1 / f2), "floor", arg, divisor);
 765           FLOAT_TO_INT2 (f1, arg, "floor", arg, divisor);
 766           return arg;
 767         }
 768 #endif
 769
 770       i1 = XINT (arg);
 771       i2 = XINT (divisor);
 772
 773       if (i2 == 0)
 774         Fsignal (Qarith_error, Qnil);
 775
 776       /* With C's /, the result is implementation-defined if either operand
 777          is negative, so use only nonnegative operands.  */
 778       i1 = (i2 < 0
 779             ? (i1 <= 0  ?  -i1 / -i2  :  -1 - ((i1 - 1) / -i2))
 780             : (i1 < 0  ?  -1 - ((-1 - i1) / i2)  :  i1 / i2));
 781
 782       XSET (arg, Lisp_Int, i1);
 783       return arg;
 784     }
 785
 786 #ifdef LISP_FLOAT_TYPE
 787   if (XTYPE (arg) == Lisp_Float)
 788     {
 789       double d;
 790       IN_FLOAT (d = floor (XFLOAT (arg)->data), "floor", arg);
 791       FLOAT_TO_INT (d, arg, "floor", arg);
 792     }
 793 #endif
 794
 795   return arg;
 796 }
 797
 798 #ifdef LISP_FLOAT_TYPE
 799
 800 DEFUN ("round", Fround, Sround, 1, 1, 0,
 801   "Return the nearest integer to ARG.")
 802   (arg)
 803      register Lisp_Object arg;
 804 {
 805   CHECK_NUMBER_OR_FLOAT (arg, 0);
 806
 807   if (XTYPE (arg) == Lisp_Float)
 808     {
 809       double d;
 810
 811       /* Screw the prevailing rounding mode.  */
 812       IN_FLOAT (d = rint (XFLOAT (arg)->data), "round", arg);
 813       FLOAT_TO_INT (d, arg, "round", arg);
 814     }
 815
 816   return arg;
 817 }
 818
 819 DEFUN ("truncate", Ftruncate, Struncate, 1, 1, 0,
 820        "Truncate a floating point number to an int.\n\
 821 Rounds the value toward zero.")
 822   (arg)
 823      register Lisp_Object arg;
 824 {
 825   CHECK_NUMBER_OR_FLOAT (arg, 0);
 826
 827   if (XTYPE (arg) == Lisp_Float)
 828     {
 829       double d;
 830
 831       d = XFLOAT (arg)->data;
 832       FLOAT_TO_INT (d, arg, "truncate", arg);
 833     }
 834
 835   return arg;
 836 }
 837 \f
 838 /* It's not clear these are worth adding.  */
 839
 840 DEFUN ("fceiling", Ffceiling, Sfceiling, 1, 1, 0,
 841   "Return the smallest integer no less than ARG, as a float.\n\
 842 \(Round toward +inf.\)")
 843   (arg)
 844      register Lisp_Object arg;
 845 {
 846   double d = extract_float (arg);
 847   IN_FLOAT (d = ceil (d), "fceiling", arg);
 848   return make_float (d);
 849 }
 850
 851 DEFUN ("ffloor", Fffloor, Sffloor, 1, 1, 0,
 852   "Return the largest integer no greater than ARG, as a float.\n\
 853 \(Round towards -inf.\)")
 854   (arg)
 855      register Lisp_Object arg;
 856 {
 857   double d = extract_float (arg);
 858   IN_FLOAT (d = floor (d), "ffloor", arg);
 859   return make_float (d);
 860 }
 861
 862 DEFUN ("fround", Ffround, Sfround, 1, 1, 0,
 863   "Return the nearest integer to ARG, as a float.")
 864   (arg)
 865      register Lisp_Object arg;
 866 {
 867   double d = extract_float (arg);
 868   IN_FLOAT (d = rint (d), "fround", arg);
 869   return make_float (d);
 870 }
 871
 872 DEFUN ("ftruncate", Fftruncate, Sftruncate, 1, 1, 0,
 873        "Truncate a floating point number to an integral float value.\n\
 874 Rounds the value toward zero.")
 875   (arg)
 876      register Lisp_Object arg;
 877 {
 878   double d = extract_float (arg);
 879   if (d >= 0.0)
 880     IN_FLOAT (d = floor (d), "ftruncate", arg);
 881   else
 882     IN_FLOAT (d = ceil (d), "ftruncate", arg);
 883   return make_float (d);
 884 }
 885 \f
 886 #ifdef FLOAT_CATCH_SIGILL
 887 static SIGTYPE
 888 float_error (signo)
 889      int signo;
 890 {
 891   if (! in_float)
 892     fatal_error_signal (signo);
 893
 894 #ifdef BSD
 895 #ifdef BSD4_1
 896   sigrelse (SIGILL);
 897 #else /* not BSD4_1 */
 898   sigsetmask (SIGEMPTYMASK);
 899 #endif /* not BSD4_1 */
 900 #else
 901   /* Must reestablish handler each time it is called.  */
 902   signal (SIGILL, float_error);
 903 #endif /* BSD */
 904
 905   in_float = 0;
 906
 907   Fsignal (Qarith_error, Fcons (float_error_arg, Qnil));
 908 }
 909
 910 /* Another idea was to replace the library function `infnan'
 911    where SIGILL is signaled.  */
 912
 913 #endif /* FLOAT_CATCH_SIGILL */
 914
 915 #ifdef HAVE_MATHERR
 916 int
 917 matherr (x)
 918      struct exception *x;
 919 {
 920   Lisp_Object args;
 921   if (! in_float)
 922     /* Not called from emacs-lisp float routines; do the default thing. */
 923     return 0;
 924   if (!strcmp (x->name, "pow"))
 925     x->name = "expt";
 926
 927   args
 928     = Fcons (build_string (x->name),
 929              Fcons (make_float (x->arg1),
 930                     ((!strcmp (x->name, "log") || !strcmp (x->name, "pow"))
 931                      ? Fcons (make_float (x->arg2), Qnil)
 932                      : Qnil)));
 933   switch (x->type)
 934     {
 935     case DOMAIN:        Fsignal (Qdomain_error, args);          break;
 936     case SING:          Fsignal (Qsingularity_error, args);     break;
 937     case OVERFLOW:      Fsignal (Qoverflow_error, args);        break;
 938     case UNDERFLOW:     Fsignal (Qunderflow_error, args);       break;
 939     default:            Fsignal (Qarith_error, args);           break;
 940     }
 941   return (1);   /* don't set errno or print a message */
 942 }
 943 #endif /* HAVE_MATHERR */
 944
 945 init_floatfns ()
 946 {
 947 #ifdef FLOAT_CATCH_SIGILL
 948   signal (SIGILL, float_error);
 949 #endif
 950   in_float = 0;
 951 }
 952
 953 #else /* not LISP_FLOAT_TYPE */
 954
 955 init_floatfns ()
 956 {}
 957
 958 #endif /* not LISP_FLOAT_TYPE */
 959
 960 syms_of_floatfns ()
 961 {
 962 #ifdef LISP_FLOAT_TYPE
 963   defsubr (&Sacos);
 964   defsubr (&Sasin);
 965   defsubr (&Satan);
 966   defsubr (&Scos);
 967   defsubr (&Ssin);
 968   defsubr (&Stan);
 969 #if 0
 970   defsubr (&Sacosh);
 971   defsubr (&Sasinh);
 972   defsubr (&Satanh);
 973   defsubr (&Scosh);
 974   defsubr (&Ssinh);
 975   defsubr (&Stanh);
 976   defsubr (&Sbessel_y0);
 977   defsubr (&Sbessel_y1);
 978   defsubr (&Sbessel_yn);
 979   defsubr (&Sbessel_j0);
 980   defsubr (&Sbessel_j1);
 981   defsubr (&Sbessel_jn);
 982   defsubr (&Serf);
 983   defsubr (&Serfc);
 984   defsubr (&Slog_gamma);
 985   defsubr (&Scube_root);
 986 #endif
 987   defsubr (&Sfceiling);
 988   defsubr (&Sffloor);
 989   defsubr (&Sfround);
 990   defsubr (&Sftruncate);
 991   defsubr (&Sexp);
 992   defsubr (&Sexpt);
 993   defsubr (&Slog);
 994   defsubr (&Slog10);
 995   defsubr (&Ssqrt);
 996
 997   defsubr (&Sabs);
 998   defsubr (&Sfloat);
 999   defsubr (&Slogb);
1000   defsubr (&Sceiling);
1001   defsubr (&Sround);
1002   defsubr (&Struncate);
1003 #endif /* LISP_FLOAT_TYPE */
1004   defsubr (&Sfloor);
1005 }