Fix problems caught with --enable-gcc-warnings
[emacs.git] / src / floatfns.c
blob46eb85ff8ca7f6227ad5c4981ec1dfc68a66ae0d
1 /* Primitive operations on floating point for GNU Emacs Lisp interpreter.
3 Copyright (C) 1988, 1993-1994, 1999, 2001-2015 Free Software Foundation,
4 Inc.
6 Author: Wolfgang Rupprecht
7 (according to ack.texi)
9 This file is part of GNU Emacs.
11 GNU Emacs is free software: you can redistribute it and/or modify
12 it under the terms of the GNU General Public License as published by
13 the Free Software Foundation, either version 3 of the License, or
14 (at your option) any later version.
16 GNU Emacs is distributed in the hope that it will be useful,
17 but WITHOUT ANY WARRANTY; without even the implied warranty of
18 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
19 GNU General Public License for more details.
21 You should have received a copy of the GNU General Public License
22 along with GNU Emacs. If not, see <http://www.gnu.org/licenses/>. */
25 /* C89 requires only the following math.h functions, and Emacs omits
26 the starred functions since we haven't found a use for them:
27 acos, asin, atan, atan2, ceil, cos, *cosh, exp, fabs, floor, fmod,
28 frexp, ldexp, log, log10 [via (log X 10)], *modf, pow, sin, *sinh,
29 sqrt, tan, *tanh.
31 C99 and C11 require the following math.h functions in addition to
32 the C89 functions. Of these, Emacs currently exports only the
33 starred ones to Lisp, since we haven't found a use for the others:
34 acosh, atanh, cbrt, *copysign, erf, erfc, exp2, expm1, fdim, fma,
35 fmax, fmin, fpclassify, hypot, ilogb, isfinite, isgreater,
36 isgreaterequal, isinf, isless, islessequal, islessgreater, *isnan,
37 isnormal, isunordered, lgamma, log1p, *log2 [via (log X 2)], *logb
38 (approximately), lrint/llrint, lround/llround, nan, nearbyint,
39 nextafter, nexttoward, remainder, remquo, *rint, round, scalbln,
40 scalbn, signbit, tgamma, trunc.
43 #include <config.h>
45 #include "lisp.h"
47 #include <math.h>
49 /* 'isfinite' and 'isnan' cause build failures on Solaris 10 with the
50 bundled GCC in c99 mode. Work around the bugs with simple
51 implementations that are good enough. */
52 #undef isfinite
53 #define isfinite(x) ((x) - (x) == 0)
54 #undef isnan
55 #define isnan(x) ((x) != (x))
57 /* Check that X is a floating point number. */
59 static void
60 CHECK_FLOAT (Lisp_Object x)
62 CHECK_TYPE (FLOATP (x), Qfloatp, x);
65 /* Extract a Lisp number as a `double', or signal an error. */
67 double
68 extract_float (Lisp_Object num)
70 CHECK_NUMBER_OR_FLOAT (num);
72 if (FLOATP (num))
73 return XFLOAT_DATA (num);
74 return (double) XINT (num);
77 /* Trig functions. */
79 DEFUN ("acos", Facos, Sacos, 1, 1, 0,
80 doc: /* Return the inverse cosine of ARG. */)
81 (Lisp_Object arg)
83 double d = extract_float (arg);
84 d = acos (d);
85 return make_float (d);
88 DEFUN ("asin", Fasin, Sasin, 1, 1, 0,
89 doc: /* Return the inverse sine of ARG. */)
90 (Lisp_Object arg)
92 double d = extract_float (arg);
93 d = asin (d);
94 return make_float (d);
97 DEFUN ("atan", Fatan, Satan, 1, 2, 0,
98 doc: /* Return the inverse tangent of the arguments.
99 If only one argument Y is given, return the inverse tangent of Y.
100 If two arguments Y and X are given, return the inverse tangent of Y
101 divided by X, i.e. the angle in radians between the vector (X, Y)
102 and the x-axis. */)
103 (Lisp_Object y, Lisp_Object x)
105 double d = extract_float (y);
107 if (NILP (x))
108 d = atan (d);
109 else
111 double d2 = extract_float (x);
112 d = atan2 (d, d2);
114 return make_float (d);
117 DEFUN ("cos", Fcos, Scos, 1, 1, 0,
118 doc: /* Return the cosine of ARG. */)
119 (Lisp_Object arg)
121 double d = extract_float (arg);
122 d = cos (d);
123 return make_float (d);
126 DEFUN ("sin", Fsin, Ssin, 1, 1, 0,
127 doc: /* Return the sine of ARG. */)
128 (Lisp_Object arg)
130 double d = extract_float (arg);
131 d = sin (d);
132 return make_float (d);
135 DEFUN ("tan", Ftan, Stan, 1, 1, 0,
136 doc: /* Return the tangent of ARG. */)
137 (Lisp_Object arg)
139 double d = extract_float (arg);
140 d = tan (d);
141 return make_float (d);
144 DEFUN ("isnan", Fisnan, Sisnan, 1, 1, 0,
145 doc: /* Return non nil if argument X is a NaN. */)
146 (Lisp_Object x)
148 CHECK_FLOAT (x);
149 return isnan (XFLOAT_DATA (x)) ? Qt : Qnil;
152 #ifdef HAVE_COPYSIGN
153 DEFUN ("copysign", Fcopysign, Scopysign, 2, 2, 0,
154 doc: /* Copy sign of X2 to value of X1, and return the result.
155 Cause an error if X1 or X2 is not a float. */)
156 (Lisp_Object x1, Lisp_Object x2)
158 double f1, f2;
160 CHECK_FLOAT (x1);
161 CHECK_FLOAT (x2);
163 f1 = XFLOAT_DATA (x1);
164 f2 = XFLOAT_DATA (x2);
166 return make_float (copysign (f1, f2));
168 #endif
170 DEFUN ("frexp", Ffrexp, Sfrexp, 1, 1, 0,
171 doc: /* Get significand and exponent of a floating point number.
172 Breaks the floating point number X into its binary significand SGNFCAND
173 (a floating point value between 0.5 (included) and 1.0 (excluded))
174 and an integral exponent EXP for 2, such that:
176 X = SGNFCAND * 2^EXP
178 The function returns the cons cell (SGNFCAND . EXP).
179 If X is zero, both parts (SGNFCAND and EXP) are zero. */)
180 (Lisp_Object x)
182 double f = XFLOATINT (x);
183 int exponent;
184 double sgnfcand = frexp (f, &exponent);
185 return Fcons (make_float (sgnfcand), make_number (exponent));
188 DEFUN ("ldexp", Fldexp, Sldexp, 2, 2, 0,
189 doc: /* Return X * 2**EXP, as a floating point number.
190 EXP must be an integer. */)
191 (Lisp_Object sgnfcand, Lisp_Object exponent)
193 CHECK_NUMBER (exponent);
194 int e = min (max (INT_MIN, XINT (exponent)), INT_MAX);
195 return make_float (ldexp (XFLOATINT (sgnfcand), e));
198 DEFUN ("exp", Fexp, Sexp, 1, 1, 0,
199 doc: /* Return the exponential base e of ARG. */)
200 (Lisp_Object arg)
202 double d = extract_float (arg);
203 d = exp (d);
204 return make_float (d);
207 DEFUN ("expt", Fexpt, Sexpt, 2, 2, 0,
208 doc: /* Return the exponential ARG1 ** ARG2. */)
209 (Lisp_Object arg1, Lisp_Object arg2)
211 double f1, f2, f3;
213 CHECK_NUMBER_OR_FLOAT (arg1);
214 CHECK_NUMBER_OR_FLOAT (arg2);
215 if (INTEGERP (arg1) /* common lisp spec */
216 && INTEGERP (arg2) /* don't promote, if both are ints, and */
217 && XINT (arg2) >= 0) /* we are sure the result is not fractional */
218 { /* this can be improved by pre-calculating */
219 EMACS_INT y; /* some binary powers of x then accumulating */
220 EMACS_UINT acc, x; /* Unsigned so that overflow is well defined. */
221 Lisp_Object val;
223 x = XINT (arg1);
224 y = XINT (arg2);
225 acc = (y & 1 ? x : 1);
227 while ((y >>= 1) != 0)
229 x *= x;
230 if (y & 1)
231 acc *= x;
233 XSETINT (val, acc);
234 return val;
236 f1 = FLOATP (arg1) ? XFLOAT_DATA (arg1) : XINT (arg1);
237 f2 = FLOATP (arg2) ? XFLOAT_DATA (arg2) : XINT (arg2);
238 f3 = pow (f1, f2);
239 return make_float (f3);
242 DEFUN ("log", Flog, Slog, 1, 2, 0,
243 doc: /* Return the natural logarithm of ARG.
244 If the optional argument BASE is given, return log ARG using that base. */)
245 (Lisp_Object arg, Lisp_Object base)
247 double d = extract_float (arg);
249 if (NILP (base))
250 d = log (d);
251 else
253 double b = extract_float (base);
255 if (b == 10.0)
256 d = log10 (d);
257 #if HAVE_LOG2
258 else if (b == 2.0)
259 d = log2 (d);
260 #endif
261 else
262 d = log (d) / log (b);
264 return make_float (d);
267 DEFUN ("sqrt", Fsqrt, Ssqrt, 1, 1, 0,
268 doc: /* Return the square root of ARG. */)
269 (Lisp_Object arg)
271 double d = extract_float (arg);
272 d = sqrt (d);
273 return make_float (d);
276 DEFUN ("abs", Fabs, Sabs, 1, 1, 0,
277 doc: /* Return the absolute value of ARG. */)
278 (register Lisp_Object arg)
280 CHECK_NUMBER_OR_FLOAT (arg);
282 if (FLOATP (arg))
283 arg = make_float (fabs (XFLOAT_DATA (arg)));
284 else if (XINT (arg) < 0)
285 XSETINT (arg, - XINT (arg));
287 return arg;
290 DEFUN ("float", Ffloat, Sfloat, 1, 1, 0,
291 doc: /* Return the floating point number equal to ARG. */)
292 (register Lisp_Object arg)
294 CHECK_NUMBER_OR_FLOAT (arg);
296 if (INTEGERP (arg))
297 return make_float ((double) XINT (arg));
298 else /* give 'em the same float back */
299 return arg;
302 DEFUN ("logb", Flogb, Slogb, 1, 1, 0,
303 doc: /* Returns largest integer <= the base 2 log of the magnitude of ARG.
304 This is the same as the exponent of a float. */)
305 (Lisp_Object arg)
307 Lisp_Object val;
308 EMACS_INT value;
309 double f = extract_float (arg);
311 if (f == 0.0)
312 value = MOST_NEGATIVE_FIXNUM;
313 else if (isfinite (f))
315 int ivalue;
316 frexp (f, &ivalue);
317 value = ivalue - 1;
319 else
320 value = MOST_POSITIVE_FIXNUM;
322 XSETINT (val, value);
323 return val;
327 /* the rounding functions */
329 static Lisp_Object
330 rounding_driver (Lisp_Object arg, Lisp_Object divisor,
331 double (*double_round) (double),
332 EMACS_INT (*int_round2) (EMACS_INT, EMACS_INT),
333 const char *name)
335 CHECK_NUMBER_OR_FLOAT (arg);
337 if (! NILP (divisor))
339 EMACS_INT i1, i2;
341 CHECK_NUMBER_OR_FLOAT (divisor);
343 if (FLOATP (arg) || FLOATP (divisor))
345 double f1, f2;
347 f1 = FLOATP (arg) ? XFLOAT_DATA (arg) : XINT (arg);
348 f2 = (FLOATP (divisor) ? XFLOAT_DATA (divisor) : XINT (divisor));
349 if (! IEEE_FLOATING_POINT && f2 == 0)
350 xsignal0 (Qarith_error);
352 f1 = (*double_round) (f1 / f2);
353 if (FIXNUM_OVERFLOW_P (f1))
354 xsignal3 (Qrange_error, build_string (name), arg, divisor);
355 arg = make_number (f1);
356 return arg;
359 i1 = XINT (arg);
360 i2 = XINT (divisor);
362 if (i2 == 0)
363 xsignal0 (Qarith_error);
365 XSETINT (arg, (*int_round2) (i1, i2));
366 return arg;
369 if (FLOATP (arg))
371 double d = (*double_round) (XFLOAT_DATA (arg));
372 if (FIXNUM_OVERFLOW_P (d))
373 xsignal2 (Qrange_error, build_string (name), arg);
374 arg = make_number (d);
377 return arg;
380 static EMACS_INT
381 ceiling2 (EMACS_INT i1, EMACS_INT i2)
383 return i1 / i2 + ((i1 % i2 != 0) & ((i1 < 0) == (i2 < 0)));
386 static EMACS_INT
387 floor2 (EMACS_INT i1, EMACS_INT i2)
389 return i1 / i2 - ((i1 % i2 != 0) & ((i1 < 0) != (i2 < 0)));
392 static EMACS_INT
393 truncate2 (EMACS_INT i1, EMACS_INT i2)
395 return i1 / i2;
398 static EMACS_INT
399 round2 (EMACS_INT i1, EMACS_INT i2)
401 /* The C language's division operator gives us one remainder R, but
402 we want the remainder R1 on the other side of 0 if R1 is closer
403 to 0 than R is; because we want to round to even, we also want R1
404 if R and R1 are the same distance from 0 and if C's quotient is
405 odd. */
406 EMACS_INT q = i1 / i2;
407 EMACS_INT r = i1 % i2;
408 EMACS_INT abs_r = eabs (r);
409 EMACS_INT abs_r1 = eabs (i2) - abs_r;
410 return q + (abs_r + (q & 1) <= abs_r1 ? 0 : (i2 ^ r) < 0 ? -1 : 1);
413 /* The code uses emacs_rint, so that it works to undefine HAVE_RINT
414 if `rint' exists but does not work right. */
415 #ifdef HAVE_RINT
416 #define emacs_rint rint
417 #else
418 static double
419 emacs_rint (double d)
421 double d1 = d + 0.5;
422 double r = floor (d1);
423 return r - (r == d1 && fmod (r, 2) != 0);
425 #endif
427 static double
428 double_identity (double d)
430 return d;
433 DEFUN ("ceiling", Fceiling, Sceiling, 1, 2, 0,
434 doc: /* Return the smallest integer no less than ARG.
435 This rounds the value towards +inf.
436 With optional DIVISOR, return the smallest integer no less than ARG/DIVISOR. */)
437 (Lisp_Object arg, Lisp_Object divisor)
439 return rounding_driver (arg, divisor, ceil, ceiling2, "ceiling");
442 DEFUN ("floor", Ffloor, Sfloor, 1, 2, 0,
443 doc: /* Return the largest integer no greater than ARG.
444 This rounds the value towards -inf.
445 With optional DIVISOR, return the largest integer no greater than ARG/DIVISOR. */)
446 (Lisp_Object arg, Lisp_Object divisor)
448 return rounding_driver (arg, divisor, floor, floor2, "floor");
451 DEFUN ("round", Fround, Sround, 1, 2, 0,
452 doc: /* Return the nearest integer to ARG.
453 With optional DIVISOR, return the nearest integer to ARG/DIVISOR.
455 Rounding a value equidistant between two integers may choose the
456 integer closer to zero, or it may prefer an even integer, depending on
457 your machine. For example, (round 2.5) can return 3 on some
458 systems, but 2 on others. */)
459 (Lisp_Object arg, Lisp_Object divisor)
461 return rounding_driver (arg, divisor, emacs_rint, round2, "round");
464 DEFUN ("truncate", Ftruncate, Struncate, 1, 2, 0,
465 doc: /* Truncate a floating point number to an int.
466 Rounds ARG toward zero.
467 With optional DIVISOR, truncate ARG/DIVISOR. */)
468 (Lisp_Object arg, Lisp_Object divisor)
470 return rounding_driver (arg, divisor, double_identity, truncate2,
471 "truncate");
475 Lisp_Object
476 fmod_float (Lisp_Object x, Lisp_Object y)
478 double f1, f2;
480 f1 = FLOATP (x) ? XFLOAT_DATA (x) : XINT (x);
481 f2 = FLOATP (y) ? XFLOAT_DATA (y) : XINT (y);
483 f1 = fmod (f1, f2);
485 /* If the "remainder" comes out with the wrong sign, fix it. */
486 if (f2 < 0 ? f1 > 0 : f1 < 0)
487 f1 += f2;
489 return make_float (f1);
492 DEFUN ("fceiling", Ffceiling, Sfceiling, 1, 1, 0,
493 doc: /* Return the smallest integer no less than ARG, as a float.
494 (Round toward +inf.) */)
495 (Lisp_Object arg)
497 double d = extract_float (arg);
498 d = ceil (d);
499 return make_float (d);
502 DEFUN ("ffloor", Fffloor, Sffloor, 1, 1, 0,
503 doc: /* Return the largest integer no greater than ARG, as a float.
504 (Round towards -inf.) */)
505 (Lisp_Object arg)
507 double d = extract_float (arg);
508 d = floor (d);
509 return make_float (d);
512 DEFUN ("fround", Ffround, Sfround, 1, 1, 0,
513 doc: /* Return the nearest integer to ARG, as a float. */)
514 (Lisp_Object arg)
516 double d = extract_float (arg);
517 d = emacs_rint (d);
518 return make_float (d);
521 DEFUN ("ftruncate", Fftruncate, Sftruncate, 1, 1, 0,
522 doc: /* Truncate a floating point number to an integral float value.
523 Rounds the value toward zero. */)
524 (Lisp_Object arg)
526 double d = extract_float (arg);
527 if (d >= 0.0)
528 d = floor (d);
529 else
530 d = ceil (d);
531 return make_float (d);
534 void
535 syms_of_floatfns (void)
537 defsubr (&Sacos);
538 defsubr (&Sasin);
539 defsubr (&Satan);
540 defsubr (&Scos);
541 defsubr (&Ssin);
542 defsubr (&Stan);
543 defsubr (&Sisnan);
544 #ifdef HAVE_COPYSIGN
545 defsubr (&Scopysign);
546 #endif
547 defsubr (&Sfrexp);
548 defsubr (&Sldexp);
549 defsubr (&Sfceiling);
550 defsubr (&Sffloor);
551 defsubr (&Sfround);
552 defsubr (&Sftruncate);
553 defsubr (&Sexp);
554 defsubr (&Sexpt);
555 defsubr (&Slog);
556 defsubr (&Ssqrt);
558 defsubr (&Sabs);
559 defsubr (&Sfloat);
560 defsubr (&Slogb);
561 defsubr (&Sceiling);
562 defsubr (&Sfloor);
563 defsubr (&Sround);
564 defsubr (&Struncate);