XEmacs 21.2.30 "Hygeia".
[chise/xemacs-chise.git.1] / src / floatfns.c
1 /* Primitive operations on floating point for XEmacs Lisp interpreter.
2    Copyright (C) 1988, 1993, 1994 Free Software Foundation, Inc.
3
4 This file is part of XEmacs.
5
6 XEmacs is free software; you can redistribute it and/or modify it
7 under the terms of the GNU General Public License as published by the
8 Free Software Foundation; either version 2, or (at your option) any
9 later version.
10
11 XEmacs is distributed in the hope that it will be useful, but WITHOUT
12 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
13 FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
14 for more details.
15
16 You should have received a copy of the GNU General Public License
17 along with XEmacs; 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 /* Synched up with: FSF 19.30. */
22
23 /* ANSI C requires only these float functions:
24    acos, asin, atan, atan2, ceil, cos, cosh, exp, fabs, floor, fmod,
25    frexp, ldexp, log, log10, modf, pow, sin, sinh, sqrt, tan, tanh.
26
27    Define HAVE_INVERSE_HYPERBOLIC if you have acosh, asinh, and atanh.
28    Define HAVE_CBRT if you have cbrt().
29    Define HAVE_RINT if you have rint().
30    If you don't define these, then the appropriate routines will be simulated.
31
32    Define HAVE_MATHERR if on a system supporting the SysV matherr() callback.
33    (This should happen automatically.)
34
35    Define FLOAT_CHECK_ERRNO if the float library routines set errno.
36    This has no effect if HAVE_MATHERR is defined.
37
38    Define FLOAT_CATCH_SIGILL if the float library routines signal SIGILL.
39    (What systems actually do this?  Let me know. -jwz)
40
41    Define FLOAT_CHECK_DOMAIN if the float library doesn't handle errors by
42    either setting errno, or signalling SIGFPE/SIGILL.  Otherwise, domain and
43    range checking will happen before calling the float routines.  This has
44    no effect if HAVE_MATHERR is defined (since matherr will be called when
45    a domain error occurs).
46  */
47
48 #include <config.h>
49 #include "lisp.h"
50 #include "syssignal.h"
51
52 #ifdef LISP_FLOAT_TYPE
53
54 /* Need to define a differentiating symbol -- see sysfloat.h */
55 #define THIS_FILENAME floatfns
56 #include "sysfloat.h"
57
58 /* The code uses emacs_rint, so that it works to undefine HAVE_RINT
59    if `rint' exists but does not work right.  */
60 #ifdef HAVE_RINT
61 #define emacs_rint rint
62 #else
63 static double
64 emacs_rint (double x)
65 {
66   double r = floor (x + 0.5);
67   double diff = fabs (r - x);
68   /* Round to even and correct for any roundoff errors.  */
69   if (diff >= 0.5 && (diff > 0.5 || r != 2.0 * floor (r / 2.0)))
70     r += r < x ? 1.0 : -1.0;
71   return r;
72 }
73 #endif
74
75 /* Nonzero while executing in floating point.
76    This tells float_error what to do.  */
77 static int in_float;
78
79 /* If an argument is out of range for a mathematical function,
80    here is the actual argument value to use in the error message.  */
81 static Lisp_Object float_error_arg, float_error_arg2;
82 static const char *float_error_fn_name;
83
84 /* Evaluate the floating point expression D, recording NUM
85    as the original argument for error messages.
86    D is normally an assignment expression.
87    Handle errors which may result in signals or may set errno.
88
89    Note that float_error may be declared to return void, so you can't
90    just cast the zero after the colon to (SIGTYPE) to make the types
91    check properly.  */
92 #ifdef FLOAT_CHECK_ERRNO
93 #define IN_FLOAT(d, name, num)                          \
94   do {                                                  \
95     float_error_arg = num;                              \
96     float_error_fn_name = name;                         \
97     in_float = 1; errno = 0; (d); in_float = 0;         \
98     if (errno != 0) in_float_error ();                  \
99   } while (0)
100 #define IN_FLOAT2(d, name, num, num2)                   \
101   do {                                                  \
102     float_error_arg = num;                              \
103     float_error_arg2 = num2;                            \
104     float_error_fn_name = name;                         \
105     in_float = 2; errno = 0; (d); in_float = 0;         \
106     if (errno != 0) in_float_error ();                  \
107   } while (0)
108 #else
109 #define IN_FLOAT(d, name, num) (in_float = 1, (d), in_float = 0)
110 #define IN_FLOAT2(d, name, num, num2) (in_float = 2, (d), in_float = 0)
111 #endif
112
113
114 #define arith_error(op,arg) \
115   Fsignal (Qarith_error, list2 (build_string (op), arg))
116 #define range_error(op,arg) \
117   Fsignal (Qrange_error, list2 (build_string (op), arg))
118 #define range_error2(op,a1,a2) \
119   Fsignal (Qrange_error, list3 (build_string (op), a1, a2))
120 #define domain_error(op,arg) \
121   Fsignal (Qdomain_error, list2 (build_string (op), arg))
122 #define domain_error2(op,a1,a2) \
123   Fsignal (Qdomain_error, list3 (build_string (op), a1, a2))
124
125
126 /* Convert float to Lisp Integer if it fits, else signal a range
127    error using the given arguments.  */
128 static Lisp_Object
129 float_to_int (double x, const char *name, Lisp_Object num, Lisp_Object num2)
130 {
131   if (x >= ((EMACS_INT) 1 << (VALBITS-1))
132       || x <= - ((EMACS_INT) 1 << (VALBITS-1)) - (EMACS_INT) 1)
133   {
134     if (!UNBOUNDP (num2))
135       range_error2 (name, num, num2);
136     else
137       range_error (name, num);
138   }
139   return (make_int ((EMACS_INT) x));
140 }
141
142
143 static void
144 in_float_error (void)
145 {
146   switch (errno)
147   {
148   case 0:
149     break;
150   case EDOM:
151     if (in_float == 2)
152       domain_error2 (float_error_fn_name, float_error_arg, float_error_arg2);
153     else
154       domain_error (float_error_fn_name, float_error_arg);
155     break;
156   case ERANGE:
157     range_error (float_error_fn_name, float_error_arg);
158     break;
159   default:
160     arith_error (float_error_fn_name, float_error_arg);
161     break;
162   }
163 }
164
165 \f
166 static Lisp_Object
167 mark_float (Lisp_Object obj)
168 {
169   return Qnil;
170 }
171
172 static int
173 float_equal (Lisp_Object obj1, Lisp_Object obj2, int depth)
174 {
175   return (extract_float (obj1) == extract_float (obj2));
176 }
177
178 static unsigned long
179 float_hash (Lisp_Object obj, int depth)
180 {
181   /* mod the value down to 32-bit range */
182   /* #### change for 64-bit machines */
183   return (unsigned long) fmod (extract_float (obj), 4e9);
184 }
185
186 static const struct lrecord_description float_description[] = {
187   { XD_END }
188 };
189
190 DEFINE_BASIC_LRECORD_IMPLEMENTATION ("float", float,
191                                      mark_float, print_float, 0, float_equal,
192                                      float_hash, float_description,
193                                      Lisp_Float);
194 \f
195 /* Extract a Lisp number as a `double', or signal an error.  */
196
197 double
198 extract_float (Lisp_Object num)
199 {
200   if (FLOATP (num))
201     return XFLOAT_DATA (num);
202
203   if (INTP (num))
204     return (double) XINT (num);
205
206   return extract_float (wrong_type_argument (Qnumberp, num));
207 }
208 #endif /* LISP_FLOAT_TYPE */
209
210 \f
211 /* Trig functions.  */
212 #ifdef LISP_FLOAT_TYPE
213
214 DEFUN ("acos", Facos, 1, 1, 0, /*
215 Return the inverse cosine of ARG.
216 */
217        (arg))
218 {
219   double d = extract_float (arg);
220 #ifdef FLOAT_CHECK_DOMAIN
221   if (d > 1.0 || d < -1.0)
222     domain_error ("acos", arg);
223 #endif
224   IN_FLOAT (d = acos (d), "acos", arg);
225   return make_float (d);
226 }
227
228 DEFUN ("asin", Fasin, 1, 1, 0, /*
229 Return the inverse sine of ARG.
230 */
231        (arg))
232 {
233   double d = extract_float (arg);
234 #ifdef FLOAT_CHECK_DOMAIN
235   if (d > 1.0 || d < -1.0)
236     domain_error ("asin", arg);
237 #endif
238   IN_FLOAT (d = asin (d), "asin", arg);
239   return make_float (d);
240 }
241
242 DEFUN ("atan", Fatan, 1, 2, 0, /*
243 Return the inverse tangent of ARG.
244 */
245        (arg1, arg2))
246 {
247   double d = extract_float (arg1);
248
249   if (NILP (arg2))
250     IN_FLOAT (d = atan (d), "atan", arg1);
251   else
252     {
253       double d2 = extract_float (arg2);
254 #ifdef FLOAT_CHECK_DOMAIN
255       if (d == 0.0 && d2 == 0.0)
256         domain_error2 ("atan", arg1, arg2);
257 #endif
258       IN_FLOAT2 (d = atan2 (d, d2), "atan", arg1, arg2);
259     }
260   return make_float (d);
261 }
262
263 DEFUN ("cos", Fcos, 1, 1, 0, /*
264 Return the cosine of ARG.
265 */
266        (arg))
267 {
268   double d = extract_float (arg);
269   IN_FLOAT (d = cos (d), "cos", arg);
270   return make_float (d);
271 }
272
273 DEFUN ("sin", Fsin, 1, 1, 0, /*
274 Return the sine of ARG.
275 */
276        (arg))
277 {
278   double d = extract_float (arg);
279   IN_FLOAT (d = sin (d), "sin", arg);
280   return make_float (d);
281 }
282
283 DEFUN ("tan", Ftan, 1, 1, 0, /*
284 Return the tangent of ARG.
285 */
286        (arg))
287 {
288   double d = extract_float (arg);
289   double c = cos (d);
290 #ifdef FLOAT_CHECK_DOMAIN
291   if (c == 0.0)
292     domain_error ("tan", arg);
293 #endif
294   IN_FLOAT (d = (sin (d) / c), "tan", arg);
295   return make_float (d);
296 }
297 #endif /* LISP_FLOAT_TYPE (trig functions) */
298
299 \f
300 /* Bessel functions */
301 #if 0 /* Leave these out unless we find there's a reason for them.  */
302 /* #ifdef LISP_FLOAT_TYPE */
303
304 DEFUN ("bessel-j0", Fbessel_j0, 1, 1, 0, /*
305 Return the bessel function j0 of ARG.
306 */
307        (arg))
308 {
309   double d = extract_float (arg);
310   IN_FLOAT (d = j0 (d), "bessel-j0", arg);
311   return make_float (d);
312 }
313
314 DEFUN ("bessel-j1", Fbessel_j1, 1, 1, 0, /*
315 Return the bessel function j1 of ARG.
316 */
317        (arg))
318 {
319   double d = extract_float (arg);
320   IN_FLOAT (d = j1 (d), "bessel-j1", arg);
321   return make_float (d);
322 }
323
324 DEFUN ("bessel-jn", Fbessel_jn, 2, 2, 0, /*
325 Return the order N bessel function output jn of ARG.
326 The first arg (the order) is truncated to an integer.
327 */
328        (arg1, arg2))
329 {
330   int i1 = extract_float (arg1);
331   double f2 = extract_float (arg2);
332
333   IN_FLOAT (f2 = jn (i1, f2), "bessel-jn", arg1);
334   return make_float (f2);
335 }
336
337 DEFUN ("bessel-y0", Fbessel_y0, 1, 1, 0, /*
338 Return the bessel function y0 of ARG.
339 */
340        (arg))
341 {
342   double d = extract_float (arg);
343   IN_FLOAT (d = y0 (d), "bessel-y0", arg);
344   return make_float (d);
345 }
346
347 DEFUN ("bessel-y1", Fbessel_y1, 1, 1, 0, /*
348 Return the bessel function y1 of ARG.
349 */
350        (arg))
351 {
352   double d = extract_float (arg);
353   IN_FLOAT (d = y1 (d), "bessel-y0", arg);
354   return make_float (d);
355 }
356
357 DEFUN ("bessel-yn", Fbessel_yn, 2, 2, 0, /*
358 Return the order N bessel function output yn of ARG.
359 The first arg (the order) is truncated to an integer.
360 */
361        (arg1, arg2))
362 {
363   int i1 = extract_float (arg1);
364   double f2 = extract_float (arg2);
365
366   IN_FLOAT (f2 = yn (i1, f2), "bessel-yn", arg1);
367   return make_float (f2);
368 }
369
370 #endif /* 0 (bessel functions) */
371 \f
372 /* Error functions. */
373 #if 0 /* Leave these out unless we see they are worth having.  */
374 /* #ifdef LISP_FLOAT_TYPE */
375
376 DEFUN ("erf", Ferf, 1, 1, 0, /*
377 Return the mathematical error function of ARG.
378 */
379        (arg))
380 {
381   double d = extract_float (arg);
382   IN_FLOAT (d = erf (d), "erf", arg);
383   return make_float (d);
384 }
385
386 DEFUN ("erfc", Ferfc, 1, 1, 0, /*
387 Return the complementary error function of ARG.
388 */
389        (arg))
390 {
391   double d = extract_float (arg);
392   IN_FLOAT (d = erfc (d), "erfc", arg);
393   return make_float (d);
394 }
395
396 DEFUN ("log-gamma", Flog_gamma, 1, 1, 0, /*
397 Return the log gamma of ARG.
398 */
399        (arg))
400 {
401   double d = extract_float (arg);
402   IN_FLOAT (d = lgamma (d), "log-gamma", arg);
403   return make_float (d);
404 }
405
406 #endif /* 0 (error functions) */
407
408 \f
409 /* Root and Log functions. */
410
411 #ifdef LISP_FLOAT_TYPE
412 DEFUN ("exp", Fexp, 1, 1, 0, /*
413 Return the exponential base e of ARG.
414 */
415        (arg))
416 {
417   double d = extract_float (arg);
418 #ifdef FLOAT_CHECK_DOMAIN
419   if (d > 709.7827)   /* Assume IEEE doubles here */
420     range_error ("exp", arg);
421   else if (d < -709.0)
422     return make_float (0.0);
423   else
424 #endif
425     IN_FLOAT (d = exp (d), "exp", arg);
426   return make_float (d);
427 }
428 #endif /* LISP_FLOAT_TYPE */
429
430
431 DEFUN ("expt", Fexpt, 2, 2, 0, /*
432 Return the exponential ARG1 ** ARG2.
433 */
434        (arg1, arg2))
435 {
436   if (INTP (arg1) && /* common lisp spec */
437       INTP (arg2)) /* don't promote, if both are ints */
438     {
439       EMACS_INT retval;
440       EMACS_INT x = XINT (arg1);
441       EMACS_INT y = XINT (arg2);
442
443       if (y < 0)
444         {
445           if (x == 1)
446             retval = 1;
447           else if (x == -1)
448             retval = (y & 1) ? -1 : 1;
449           else
450             retval = 0;
451         }
452       else
453         {
454           retval = 1;
455           while (y > 0)
456             {
457               if (y & 1)
458                 retval *= x;
459               x *= x;
460               y = (EMACS_UINT) y >> 1;
461             }
462         }
463       return make_int (retval);
464     }
465
466 #ifdef LISP_FLOAT_TYPE
467   {
468     double f1 = extract_float (arg1);
469     double f2 = extract_float (arg2);
470     /* Really should check for overflow, too */
471     if (f1 == 0.0 && f2 == 0.0)
472       f1 = 1.0;
473 # ifdef FLOAT_CHECK_DOMAIN
474     else if ((f1 == 0.0 && f2 < 0.0) || (f1 < 0 && f2 != floor(f2)))
475       domain_error2 ("expt", arg1, arg2);
476 # endif /* FLOAT_CHECK_DOMAIN */
477     IN_FLOAT2 (f1 = pow (f1, f2), "expt", arg1, arg2);
478     return make_float (f1);
479   }
480 #else
481   CHECK_INT_OR_FLOAT (arg1);
482   CHECK_INT_OR_FLOAT (arg2);
483   return Fexpt (arg1, arg2);
484 #endif /* LISP_FLOAT_TYPE */
485 }
486
487 #ifdef LISP_FLOAT_TYPE
488 DEFUN ("log", Flog, 1, 2, 0, /*
489 Return the natural logarithm of ARG.
490 If second optional argument BASE is given, return log ARG using that base.
491 */
492        (arg, base))
493 {
494   double d = extract_float (arg);
495 #ifdef FLOAT_CHECK_DOMAIN
496   if (d <= 0.0)
497     domain_error2 ("log", arg, base);
498 #endif
499   if (NILP (base))
500     IN_FLOAT (d = log (d), "log", arg);
501   else
502     {
503       double b = extract_float (base);
504 #ifdef FLOAT_CHECK_DOMAIN
505       if (b <= 0.0 || b == 1.0)
506         domain_error2 ("log", arg, base);
507 #endif
508       if (b == 10.0)
509         IN_FLOAT2 (d = log10 (d), "log", arg, base);
510       else
511         IN_FLOAT2 (d = (log (d) / log (b)), "log", arg, base);
512     }
513   return make_float (d);
514 }
515
516
517 DEFUN ("log10", Flog10, 1, 1, 0, /*
518 Return the logarithm base 10 of ARG.
519 */
520        (arg))
521 {
522   double d = extract_float (arg);
523 #ifdef FLOAT_CHECK_DOMAIN
524   if (d <= 0.0)
525     domain_error ("log10", arg);
526 #endif
527   IN_FLOAT (d = log10 (d), "log10", arg);
528   return make_float (d);
529 }
530
531
532 DEFUN ("sqrt", Fsqrt, 1, 1, 0, /*
533 Return the square root of ARG.
534 */
535        (arg))
536 {
537   double d = extract_float (arg);
538 #ifdef FLOAT_CHECK_DOMAIN
539   if (d < 0.0)
540     domain_error ("sqrt", arg);
541 #endif
542   IN_FLOAT (d = sqrt (d), "sqrt", arg);
543   return make_float (d);
544 }
545
546
547 DEFUN ("cube-root", Fcube_root, 1, 1, 0, /*
548 Return the cube root of ARG.
549 */
550        (arg))
551 {
552   double d = extract_float (arg);
553 #ifdef HAVE_CBRT
554   IN_FLOAT (d = cbrt (d), "cube-root", arg);
555 #else
556   if (d >= 0.0)
557     IN_FLOAT (d = pow (d, 1.0/3.0), "cube-root", arg);
558   else
559     IN_FLOAT (d = -pow (-d, 1.0/3.0), "cube-root", arg);
560 #endif
561   return make_float (d);
562 }
563 #endif /* LISP_FLOAT_TYPE */
564
565 \f
566 /* Inverse trig functions. */
567 #ifdef LISP_FLOAT_TYPE
568 /* #if 0  Not clearly worth adding...  */
569
570 DEFUN ("acosh", Facosh, 1, 1, 0, /*
571 Return the inverse hyperbolic cosine of ARG.
572 */
573        (arg))
574 {
575   double d = extract_float (arg);
576 #ifdef FLOAT_CHECK_DOMAIN
577   if (d < 1.0)
578     domain_error ("acosh", arg);
579 #endif
580 #ifdef HAVE_INVERSE_HYPERBOLIC
581   IN_FLOAT (d = acosh (d), "acosh", arg);
582 #else
583   IN_FLOAT (d = log (d + sqrt (d*d - 1.0)), "acosh", arg);
584 #endif
585   return make_float (d);
586 }
587
588 DEFUN ("asinh", Fasinh, 1, 1, 0, /*
589 Return the inverse hyperbolic sine of ARG.
590 */
591        (arg))
592 {
593   double d = extract_float (arg);
594 #ifdef HAVE_INVERSE_HYPERBOLIC
595   IN_FLOAT (d = asinh (d), "asinh", arg);
596 #else
597   IN_FLOAT (d = log (d + sqrt (d*d + 1.0)), "asinh", arg);
598 #endif
599   return make_float (d);
600 }
601
602 DEFUN ("atanh", Fatanh, 1, 1, 0, /*
603 Return the inverse hyperbolic tangent of ARG.
604 */
605        (arg))
606 {
607   double d = extract_float (arg);
608 #ifdef FLOAT_CHECK_DOMAIN
609   if (d >= 1.0 || d <= -1.0)
610     domain_error ("atanh", arg);
611 #endif
612 #ifdef HAVE_INVERSE_HYPERBOLIC
613   IN_FLOAT (d = atanh (d), "atanh", arg);
614 #else
615   IN_FLOAT (d = 0.5 * log ((1.0 + d) / (1.0 - d)), "atanh", arg);
616 #endif
617   return make_float (d);
618 }
619
620 DEFUN ("cosh", Fcosh, 1, 1, 0, /*
621 Return the hyperbolic cosine of ARG.
622 */
623        (arg))
624 {
625   double d = extract_float (arg);
626 #ifdef FLOAT_CHECK_DOMAIN
627   if (d > 710.0 || d < -710.0)
628     range_error ("cosh", arg);
629 #endif
630   IN_FLOAT (d = cosh (d), "cosh", arg);
631   return make_float (d);
632 }
633
634 DEFUN ("sinh", Fsinh, 1, 1, 0, /*
635 Return the hyperbolic sine of ARG.
636 */
637        (arg))
638 {
639   double d = extract_float (arg);
640 #ifdef FLOAT_CHECK_DOMAIN
641   if (d > 710.0 || d < -710.0)
642     range_error ("sinh", arg);
643 #endif
644   IN_FLOAT (d = sinh (d), "sinh", arg);
645   return make_float (d);
646 }
647
648 DEFUN ("tanh", Ftanh, 1, 1, 0, /*
649 Return the hyperbolic tangent of ARG.
650 */
651        (arg))
652 {
653   double d = extract_float (arg);
654   IN_FLOAT (d = tanh (d), "tanh", arg);
655   return make_float (d);
656 }
657 #endif /* LISP_FLOAT_TYPE (inverse trig functions) */
658 \f
659 /* Rounding functions */
660
661 DEFUN ("abs", Fabs, 1, 1, 0, /*
662 Return the absolute value of ARG.
663 */
664        (arg))
665 {
666 #ifdef LISP_FLOAT_TYPE
667   if (FLOATP (arg))
668     {
669       IN_FLOAT (arg = make_float (fabs (XFLOAT_DATA (arg))),
670                 "abs", arg);
671       return arg;
672     }
673 #endif /* LISP_FLOAT_TYPE */
674
675   if (INTP (arg))
676     return (XINT (arg) >= 0) ? arg : make_int (- XINT (arg));
677
678   return Fabs (wrong_type_argument (Qnumberp, arg));
679 }
680
681 #ifdef LISP_FLOAT_TYPE
682 DEFUN ("float", Ffloat, 1, 1, 0, /*
683 Return the floating point number numerically equal to ARG.
684 */
685        (arg))
686 {
687   if (INTP (arg))
688     return make_float ((double) XINT (arg));
689
690   if (FLOATP (arg))             /* give 'em the same float back */
691     return arg;
692
693   return Ffloat (wrong_type_argument (Qnumberp, arg));
694 }
695 #endif /* LISP_FLOAT_TYPE */
696
697
698 #ifdef LISP_FLOAT_TYPE
699 DEFUN ("logb", Flogb, 1, 1, 0, /*
700 Return largest integer <= the base 2 log of the magnitude of ARG.
701 This is the same as the exponent of a float.
702 */
703        (arg))
704 {
705   double f = extract_float (arg);
706
707   if (f == 0.0)
708     return make_int (- (EMACS_INT)(((EMACS_UINT) 1) << (VALBITS - 1))); /* most-negative-fixnum */
709 #ifdef HAVE_LOGB
710   {
711     Lisp_Object val;
712     IN_FLOAT (val = make_int ((EMACS_INT) logb (f)), "logb", arg);
713     return val;
714   }
715 #else
716 #ifdef HAVE_FREXP
717   {
718     int exqp;
719     IN_FLOAT (frexp (f, &exqp), "logb", arg);
720     return make_int (exqp - 1);
721   }
722 #else
723   {
724     int i;
725     double d;
726     EMACS_INT val;
727     if (f < 0.0)
728       f = -f;
729     val = -1;
730     while (f < 0.5)
731       {
732         for (i = 1, d = 0.5; d * d >= f; i += i)
733           d *= d;
734         f /= d;
735         val -= i;
736       }
737     while (f >= 1.0)
738       {
739         for (i = 1, d = 2.0; d * d <= f; i += i)
740           d *= d;
741         f /= d;
742         val += i;
743       }
744     return make_int (val);
745   }
746 #endif /* ! HAVE_FREXP */
747 #endif /* ! HAVE_LOGB */
748 }
749 #endif /* LISP_FLOAT_TYPE */
750
751
752 DEFUN ("ceiling", Fceiling, 1, 1, 0, /*
753 Return the smallest integer no less than ARG.  (Round toward +inf.)
754 */
755        (arg))
756 {
757 #ifdef LISP_FLOAT_TYPE
758   if (FLOATP (arg))
759     {
760       double d;
761       IN_FLOAT ((d = ceil (XFLOAT_DATA (arg))), "ceiling", arg);
762       return (float_to_int (d, "ceiling", arg, Qunbound));
763     }
764 #endif /* LISP_FLOAT_TYPE */
765
766   if (INTP (arg))
767     return arg;
768
769   return Fceiling (wrong_type_argument (Qnumberp, arg));
770 }
771
772
773 DEFUN ("floor", Ffloor, 1, 2, 0, /*
774 Return the largest integer no greater than ARG.  (Round towards -inf.)
775 With optional DIVISOR, return the largest integer no greater than ARG/DIVISOR.
776 */
777        (arg, divisor))
778 {
779   CHECK_INT_OR_FLOAT (arg);
780
781   if (! NILP (divisor))
782     {
783       EMACS_INT i1, i2;
784
785       CHECK_INT_OR_FLOAT (divisor);
786
787 #ifdef LISP_FLOAT_TYPE
788       if (FLOATP (arg) || FLOATP (divisor))
789         {
790           double f1 = extract_float (arg);
791           double f2 = extract_float (divisor);
792
793           if (f2 == 0)
794             Fsignal (Qarith_error, Qnil);
795
796           IN_FLOAT2 (f1 = floor (f1 / f2), "floor", arg, divisor);
797           return float_to_int (f1, "floor", arg, divisor);
798         }
799 #endif /* LISP_FLOAT_TYPE */
800
801       i1 = XINT (arg);
802       i2 = XINT (divisor);
803
804       if (i2 == 0)
805         Fsignal (Qarith_error, Qnil);
806
807       /* With C's /, the result is implementation-defined if either operand
808          is negative, so use only nonnegative operands.  */
809       i1 = (i2 < 0
810             ? (i1 <= 0  ?  -i1 / -i2  :  -1 - ((i1 - 1) / -i2))
811             : (i1 < 0  ?  -1 - ((-1 - i1) / i2)  :  i1 / i2));
812
813       return (make_int (i1));
814     }
815
816 #ifdef LISP_FLOAT_TYPE
817   if (FLOATP (arg))
818     {
819       double d;
820       IN_FLOAT ((d = floor (XFLOAT_DATA (arg))), "floor", arg);
821       return (float_to_int (d, "floor", arg, Qunbound));
822     }
823 #endif /* LISP_FLOAT_TYPE */
824
825   return arg;
826 }
827
828 DEFUN ("round", Fround, 1, 1, 0, /*
829 Return the nearest integer to ARG.
830 */
831        (arg))
832 {
833 #ifdef LISP_FLOAT_TYPE
834   if (FLOATP (arg))
835     {
836       double d;
837       /* Screw the prevailing rounding mode.  */
838       IN_FLOAT ((d = emacs_rint (XFLOAT_DATA (arg))), "round", arg);
839       return (float_to_int (d, "round", arg, Qunbound));
840     }
841 #endif /* LISP_FLOAT_TYPE */
842
843   if (INTP (arg))
844     return arg;
845
846   return Fround (wrong_type_argument (Qnumberp, arg));
847 }
848
849 DEFUN ("truncate", Ftruncate, 1, 1, 0, /*
850 Truncate a floating point number to an integer.
851 Rounds the value toward zero.
852 */
853        (arg))
854 {
855 #ifdef LISP_FLOAT_TYPE
856   if (FLOATP (arg))
857     return float_to_int (XFLOAT_DATA (arg), "truncate", arg, Qunbound);
858 #endif /* LISP_FLOAT_TYPE */
859
860   if (INTP (arg))
861     return arg;
862
863   return Ftruncate (wrong_type_argument (Qnumberp, arg));
864 }
865 \f
866 /* Float-rounding functions. */
867 #ifdef LISP_FLOAT_TYPE
868 /* #if 1  It's not clear these are worth adding... */
869
870 DEFUN ("fceiling", Ffceiling, 1, 1, 0, /*
871 Return the smallest integer no less than ARG, as a float.
872 \(Round toward +inf.\)
873 */
874        (arg))
875 {
876   double d = extract_float (arg);
877   IN_FLOAT (d = ceil (d), "fceiling", arg);
878   return make_float (d);
879 }
880
881 DEFUN ("ffloor", Fffloor, 1, 1, 0, /*
882 Return the largest integer no greater than ARG, as a float.
883 \(Round towards -inf.\)
884 */
885        (arg))
886 {
887   double d = extract_float (arg);
888   IN_FLOAT (d = floor (d), "ffloor", arg);
889   return make_float (d);
890 }
891
892 DEFUN ("fround", Ffround, 1, 1, 0, /*
893 Return the nearest integer to ARG, as a float.
894 */
895        (arg))
896 {
897   double d = extract_float (arg);
898   IN_FLOAT (d = emacs_rint (d), "fround", arg);
899   return make_float (d);
900 }
901
902 DEFUN ("ftruncate", Fftruncate, 1, 1, 0, /*
903 Truncate a floating point number to an integral float value.
904 Rounds the value toward zero.
905 */
906        (arg))
907 {
908   double d = extract_float (arg);
909   if (d >= 0.0)
910     IN_FLOAT (d = floor (d), "ftruncate", arg);
911   else
912     IN_FLOAT (d = ceil (d), "ftruncate", arg);
913   return make_float (d);
914 }
915
916 #endif /* LISP_FLOAT_TYPE (float-rounding functions) */
917
918 \f
919 #ifdef LISP_FLOAT_TYPE
920 #ifdef FLOAT_CATCH_SIGILL
921 static SIGTYPE
922 float_error (int signo)
923 {
924   if (! in_float)
925     fatal_error_signal (signo);
926
927   EMACS_REESTABLISH_SIGNAL (signo, arith_error);
928   EMACS_UNBLOCK_SIGNAL (signo);
929
930   in_float = 0;
931
932   /* Was Fsignal(), but it just doesn't make sense for an error
933      occurring inside a signal handler to be restartable, considering
934      that anything could happen when the error is signaled and trapped
935      and considering the asynchronous nature of signal handlers. */
936   signal_error (Qarith_error, list1 (float_error_arg));
937 }
938
939 /* Another idea was to replace the library function `infnan'
940    where SIGILL is signaled.  */
941
942 #endif /* FLOAT_CATCH_SIGILL */
943
944 /* In C++, it is impossible to determine what type matherr expects
945    without some more configure magic.
946    We shouldn't be using matherr anyways - it's a non-standard SYSVism. */
947 #if defined (HAVE_MATHERR) && !defined(__cplusplus)
948 int
949 matherr (struct exception *x)
950 {
951   Lisp_Object args;
952   if (! in_float)
953     /* Not called from emacs-lisp float routines; do the default thing. */
954     return 0;
955
956   /* if (!strcmp (x->name, "pow")) x->name = "expt"; */
957
958   args = Fcons (build_string (x->name),
959                 Fcons (make_float (x->arg1),
960                        ((in_float == 2)
961                         ? Fcons (make_float (x->arg2), Qnil)
962                         : Qnil)));
963   switch (x->type)
964     {
965     case DOMAIN:    Fsignal (Qdomain_error,      args); break;
966     case SING:      Fsignal (Qsingularity_error, args); break;
967     case OVERFLOW:  Fsignal (Qoverflow_error,    args); break;
968     case UNDERFLOW: Fsignal (Qunderflow_error,   args); break;
969     default:        Fsignal (Qarith_error,       args); break;
970     }
971   return 1;     /* don't set errno or print a message */
972 }
973 #endif /* HAVE_MATHERR */
974 #endif /* LISP_FLOAT_TYPE */
975
976 \f
977 void
978 init_floatfns_very_early (void)
979 {
980 #ifdef LISP_FLOAT_TYPE
981 # ifdef FLOAT_CATCH_SIGILL
982   signal (SIGILL, float_error);
983 # endif
984   in_float = 0;
985 #endif /* LISP_FLOAT_TYPE */
986 }
987
988 void
989 syms_of_floatfns (void)
990 {
991   INIT_LRECORD_IMPLEMENTATION (float);
992
993   /* Trig functions.  */
994
995 #ifdef LISP_FLOAT_TYPE
996   DEFSUBR (Facos);
997   DEFSUBR (Fasin);
998   DEFSUBR (Fatan);
999   DEFSUBR (Fcos);
1000   DEFSUBR (Fsin);
1001   DEFSUBR (Ftan);
1002 #endif /* LISP_FLOAT_TYPE */
1003
1004   /* Bessel functions */
1005
1006 #if 0
1007   DEFSUBR (Fbessel_y0);
1008   DEFSUBR (Fbessel_y1);
1009   DEFSUBR (Fbessel_yn);
1010   DEFSUBR (Fbessel_j0);
1011   DEFSUBR (Fbessel_j1);
1012   DEFSUBR (Fbessel_jn);
1013 #endif /* 0 */
1014
1015   /* Error functions. */
1016
1017 #if 0
1018   DEFSUBR (Ferf);
1019   DEFSUBR (Ferfc);
1020   DEFSUBR (Flog_gamma);
1021 #endif /* 0 */
1022
1023   /* Root and Log functions. */
1024
1025 #ifdef LISP_FLOAT_TYPE
1026   DEFSUBR (Fexp);
1027 #endif /* LISP_FLOAT_TYPE */
1028   DEFSUBR (Fexpt);
1029 #ifdef LISP_FLOAT_TYPE
1030   DEFSUBR (Flog);
1031   DEFSUBR (Flog10);
1032   DEFSUBR (Fsqrt);
1033   DEFSUBR (Fcube_root);
1034 #endif /* LISP_FLOAT_TYPE */
1035
1036   /* Inverse trig functions. */
1037
1038 #ifdef LISP_FLOAT_TYPE
1039   DEFSUBR (Facosh);
1040   DEFSUBR (Fasinh);
1041   DEFSUBR (Fatanh);
1042   DEFSUBR (Fcosh);
1043   DEFSUBR (Fsinh);
1044   DEFSUBR (Ftanh);
1045 #endif /* LISP_FLOAT_TYPE */
1046
1047   /* Rounding functions */
1048
1049   DEFSUBR (Fabs);
1050 #ifdef LISP_FLOAT_TYPE
1051   DEFSUBR (Ffloat);
1052   DEFSUBR (Flogb);
1053 #endif /* LISP_FLOAT_TYPE */
1054   DEFSUBR (Fceiling);
1055   DEFSUBR (Ffloor);
1056   DEFSUBR (Fround);
1057   DEFSUBR (Ftruncate);
1058
1059   /* Float-rounding functions. */
1060
1061 #ifdef LISP_FLOAT_TYPE
1062   DEFSUBR (Ffceiling);
1063   DEFSUBR (Fffloor);
1064   DEFSUBR (Ffround);
1065   DEFSUBR (Fftruncate);
1066 #endif /* LISP_FLOAT_TYPE */
1067 }
1068
1069 void
1070 vars_of_floatfns (void)
1071 {
1072 #ifdef LISP_FLOAT_TYPE
1073   Fprovide (intern ("lisp-float-type"));
1074 #endif
1075 }