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