1 /* Primitive operations on floating point for GNU Emacs Lisp interpreter.
3 Copyright (C) 1988, 1993-1994, 1999, 2001-2013 Free Software Foundation,
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, *modf, pow, sin, *sinh, sqrt, tan, *tanh.
30 C99 and C11 require the following math.h functions in addition to
31 the C89 functions. Of these, Emacs currently exports only the
32 starred ones to Lisp, since we haven't found a use for the others:
33 acosh, atanh, cbrt, *copysign, erf, erfc, exp2, expm1, fdim, fma,
34 fmax, fmin, fpclassify, hypot, ilogb, isfinite, isgreater,
35 isgreaterequal, isinf, isless, islessequal, islessgreater, *isnan,
36 isnormal, isunordered, lgamma, log1p, log2, *logb (approximately),
37 lrint/llrint, lround/llround, nan, nearbyint, nextafter,
38 nexttoward, remainder, remquo, *rint, round, scalbln, scalbn,
39 signbit, tgamma, trunc.
49 # define isfinite(x) ((x) - (x) == 0)
52 # define isnan(x) ((x) != (x))
55 /* Check that X is a floating point number. */
58 CHECK_FLOAT (Lisp_Object x
)
60 CHECK_TYPE (FLOATP (x
), Qfloatp
, x
);
63 /* Extract a Lisp number as a `double', or signal an error. */
66 extract_float (Lisp_Object num
)
68 CHECK_NUMBER_OR_FLOAT (num
);
71 return XFLOAT_DATA (num
);
72 return (double) XINT (num
);
77 DEFUN ("acos", Facos
, Sacos
, 1, 1, 0,
78 doc
: /* Return the inverse cosine of ARG. */)
81 double d
= extract_float (arg
);
83 return make_float (d
);
86 DEFUN ("asin", Fasin
, Sasin
, 1, 1, 0,
87 doc
: /* Return the inverse sine of ARG. */)
90 double d
= extract_float (arg
);
92 return make_float (d
);
95 DEFUN ("atan", Fatan
, Satan
, 1, 2, 0,
96 doc
: /* Return the inverse tangent of the arguments.
97 If only one argument Y is given, return the inverse tangent of Y.
98 If two arguments Y and X are given, return the inverse tangent of Y
99 divided by X, i.e. the angle in radians between the vector (X, Y)
101 (Lisp_Object y
, Lisp_Object x
)
103 double d
= extract_float (y
);
109 double d2
= extract_float (x
);
112 return make_float (d
);
115 DEFUN ("cos", Fcos
, Scos
, 1, 1, 0,
116 doc
: /* Return the cosine of ARG. */)
119 double d
= extract_float (arg
);
121 return make_float (d
);
124 DEFUN ("sin", Fsin
, Ssin
, 1, 1, 0,
125 doc
: /* Return the sine of ARG. */)
128 double d
= extract_float (arg
);
130 return make_float (d
);
133 DEFUN ("tan", Ftan
, Stan
, 1, 1, 0,
134 doc
: /* Return the tangent of ARG. */)
137 double d
= extract_float (arg
);
139 return make_float (d
);
142 DEFUN ("isnan", Fisnan
, Sisnan
, 1, 1, 0,
143 doc
: /* Return non nil iff argument X is a NaN. */)
147 return isnan (XFLOAT_DATA (x
)) ? Qt
: Qnil
;
151 DEFUN ("copysign", Fcopysign
, Scopysign
, 2, 2, 0,
152 doc
: /* Copy sign of X2 to value of X1, and return the result.
153 Cause an error if X1 or X2 is not a float. */)
154 (Lisp_Object x1
, Lisp_Object x2
)
161 f1
= XFLOAT_DATA (x1
);
162 f2
= XFLOAT_DATA (x2
);
164 return make_float (copysign (f1
, f2
));
168 DEFUN ("frexp", Ffrexp
, Sfrexp
, 1, 1, 0,
169 doc
: /* Get significand and exponent of a floating point number.
170 Breaks the floating point number X into its binary significand SGNFCAND
171 \(a floating point value between 0.5 (included) and 1.0 (excluded))
172 and an integral exponent EXP for 2, such that:
176 The function returns the cons cell (SGNFCAND . EXP).
177 If X is zero, both parts (SGNFCAND and EXP) are zero. */)
180 double f
= XFLOATINT (x
);
182 double sgnfcand
= frexp (f
, &exponent
);
183 return Fcons (make_float (sgnfcand
), make_number (exponent
));
186 DEFUN ("ldexp", Fldexp
, Sldexp
, 1, 2, 0,
187 doc
: /* Construct number X from significand SGNFCAND and exponent EXP.
188 Returns the floating point value resulting from multiplying SGNFCAND
189 (the significand) by 2 raised to the power of EXP (the exponent). */)
190 (Lisp_Object sgnfcand
, Lisp_Object exponent
)
192 CHECK_NUMBER (exponent
);
193 return make_float (ldexp (XFLOATINT (sgnfcand
), XINT (exponent
)));
196 DEFUN ("exp", Fexp
, Sexp
, 1, 1, 0,
197 doc
: /* Return the exponential base e of ARG. */)
200 double d
= extract_float (arg
);
202 return make_float (d
);
205 DEFUN ("expt", Fexpt
, Sexpt
, 2, 2, 0,
206 doc
: /* Return the exponential ARG1 ** ARG2. */)
207 (Lisp_Object arg1
, Lisp_Object arg2
)
211 CHECK_NUMBER_OR_FLOAT (arg1
);
212 CHECK_NUMBER_OR_FLOAT (arg2
);
213 if (INTEGERP (arg1
) /* common lisp spec */
214 && INTEGERP (arg2
) /* don't promote, if both are ints, and */
215 && XINT (arg2
) >= 0) /* we are sure the result is not fractional */
216 { /* this can be improved by pre-calculating */
217 EMACS_INT y
; /* some binary powers of x then accumulating */
218 EMACS_UINT acc
, x
; /* Unsigned so that overflow is well defined. */
223 acc
= (y
& 1 ? x
: 1);
225 while ((y
>>= 1) != 0)
234 f1
= FLOATP (arg1
) ? XFLOAT_DATA (arg1
) : XINT (arg1
);
235 f2
= FLOATP (arg2
) ? XFLOAT_DATA (arg2
) : XINT (arg2
);
237 return make_float (f3
);
240 DEFUN ("log", Flog
, Slog
, 1, 2, 0,
241 doc
: /* Return the natural logarithm of ARG.
242 If the optional argument BASE is given, return log ARG using that base. */)
243 (Lisp_Object arg
, Lisp_Object base
)
245 double d
= extract_float (arg
);
251 double b
= extract_float (base
);
256 d
= log (d
) / log (b
);
258 return make_float (d
);
261 DEFUN ("log10", Flog10
, Slog10
, 1, 1, 0,
262 doc
: /* Return the logarithm base 10 of ARG. */)
265 double d
= extract_float (arg
);
267 return make_float (d
);
270 DEFUN ("sqrt", Fsqrt
, Ssqrt
, 1, 1, 0,
271 doc
: /* Return the square root of ARG. */)
274 double d
= extract_float (arg
);
276 return make_float (d
);
279 DEFUN ("abs", Fabs
, Sabs
, 1, 1, 0,
280 doc
: /* Return the absolute value of ARG. */)
281 (register Lisp_Object arg
)
283 CHECK_NUMBER_OR_FLOAT (arg
);
286 arg
= make_float (fabs (XFLOAT_DATA (arg
)));
287 else if (XINT (arg
) < 0)
288 XSETINT (arg
, - XINT (arg
));
293 DEFUN ("float", Ffloat
, Sfloat
, 1, 1, 0,
294 doc
: /* Return the floating point number equal to ARG. */)
295 (register Lisp_Object arg
)
297 CHECK_NUMBER_OR_FLOAT (arg
);
300 return make_float ((double) XINT (arg
));
301 else /* give 'em the same float back */
305 DEFUN ("logb", Flogb
, Slogb
, 1, 1, 0,
306 doc
: /* Returns largest integer <= the base 2 log of the magnitude of ARG.
307 This is the same as the exponent of a float. */)
312 double f
= extract_float (arg
);
315 value
= MOST_NEGATIVE_FIXNUM
;
316 else if (isfinite (f
))
323 value
= MOST_POSITIVE_FIXNUM
;
325 XSETINT (val
, value
);
330 /* the rounding functions */
333 rounding_driver (Lisp_Object arg
, Lisp_Object divisor
,
334 double (*double_round
) (double),
335 EMACS_INT (*int_round2
) (EMACS_INT
, EMACS_INT
),
338 CHECK_NUMBER_OR_FLOAT (arg
);
340 if (! NILP (divisor
))
344 CHECK_NUMBER_OR_FLOAT (divisor
);
346 if (FLOATP (arg
) || FLOATP (divisor
))
350 f1
= FLOATP (arg
) ? XFLOAT_DATA (arg
) : XINT (arg
);
351 f2
= (FLOATP (divisor
) ? XFLOAT_DATA (divisor
) : XINT (divisor
));
352 if (! IEEE_FLOATING_POINT
&& f2
== 0)
353 xsignal0 (Qarith_error
);
355 f1
= (*double_round
) (f1
/ f2
);
356 if (FIXNUM_OVERFLOW_P (f1
))
357 xsignal3 (Qrange_error
, build_string (name
), arg
, divisor
);
358 arg
= make_number (f1
);
366 xsignal0 (Qarith_error
);
368 XSETINT (arg
, (*int_round2
) (i1
, i2
));
374 double d
= (*double_round
) (XFLOAT_DATA (arg
));
375 if (FIXNUM_OVERFLOW_P (d
))
376 xsignal2 (Qrange_error
, build_string (name
), arg
);
377 arg
= make_number (d
);
383 /* With C's /, the result is implementation-defined if either operand
384 is negative, so take care with negative operands in the following
385 integer functions. */
388 ceiling2 (EMACS_INT i1
, EMACS_INT i2
)
391 ? (i1
< 0 ? ((-1 - i1
) / -i2
) + 1 : - (i1
/ -i2
))
392 : (i1
<= 0 ? - (-i1
/ i2
) : ((i1
- 1) / i2
) + 1));
396 floor2 (EMACS_INT i1
, EMACS_INT i2
)
399 ? (i1
<= 0 ? -i1
/ -i2
: -1 - ((i1
- 1) / -i2
))
400 : (i1
< 0 ? -1 - ((-1 - i1
) / i2
) : i1
/ i2
));
404 truncate2 (EMACS_INT i1
, EMACS_INT i2
)
407 ? (i1
< 0 ? -i1
/ -i2
: - (i1
/ -i2
))
408 : (i1
< 0 ? - (-i1
/ i2
) : i1
/ i2
));
412 round2 (EMACS_INT i1
, EMACS_INT i2
)
414 /* The C language's division operator gives us one remainder R, but
415 we want the remainder R1 on the other side of 0 if R1 is closer
416 to 0 than R is; because we want to round to even, we also want R1
417 if R and R1 are the same distance from 0 and if C's quotient is
419 EMACS_INT q
= i1
/ i2
;
420 EMACS_INT r
= i1
% i2
;
421 EMACS_INT abs_r
= eabs (r
);
422 EMACS_INT abs_r1
= eabs (i2
) - abs_r
;
423 return q
+ (abs_r
+ (q
& 1) <= abs_r1
? 0 : (i2
^ r
) < 0 ? -1 : 1);
426 /* The code uses emacs_rint, so that it works to undefine HAVE_RINT
427 if `rint' exists but does not work right. */
429 #define emacs_rint rint
432 emacs_rint (double d
)
434 return floor (d
+ 0.5);
439 double_identity (double d
)
444 DEFUN ("ceiling", Fceiling
, Sceiling
, 1, 2, 0,
445 doc
: /* Return the smallest integer no less than ARG.
446 This rounds the value towards +inf.
447 With optional DIVISOR, return the smallest integer no less than ARG/DIVISOR. */)
448 (Lisp_Object arg
, Lisp_Object divisor
)
450 return rounding_driver (arg
, divisor
, ceil
, ceiling2
, "ceiling");
453 DEFUN ("floor", Ffloor
, Sfloor
, 1, 2, 0,
454 doc
: /* Return the largest integer no greater than ARG.
455 This rounds the value towards -inf.
456 With optional DIVISOR, return the largest integer no greater than ARG/DIVISOR. */)
457 (Lisp_Object arg
, Lisp_Object divisor
)
459 return rounding_driver (arg
, divisor
, floor
, floor2
, "floor");
462 DEFUN ("round", Fround
, Sround
, 1, 2, 0,
463 doc
: /* Return the nearest integer to ARG.
464 With optional DIVISOR, return the nearest integer to ARG/DIVISOR.
466 Rounding a value equidistant between two integers may choose the
467 integer closer to zero, or it may prefer an even integer, depending on
468 your machine. For example, \(round 2.5\) can return 3 on some
469 systems, but 2 on others. */)
470 (Lisp_Object arg
, Lisp_Object divisor
)
472 return rounding_driver (arg
, divisor
, emacs_rint
, round2
, "round");
475 DEFUN ("truncate", Ftruncate
, Struncate
, 1, 2, 0,
476 doc
: /* Truncate a floating point number to an int.
477 Rounds ARG toward zero.
478 With optional DIVISOR, truncate ARG/DIVISOR. */)
479 (Lisp_Object arg
, Lisp_Object divisor
)
481 return rounding_driver (arg
, divisor
, double_identity
, truncate2
,
487 fmod_float (Lisp_Object x
, Lisp_Object y
)
491 f1
= FLOATP (x
) ? XFLOAT_DATA (x
) : XINT (x
);
492 f2
= FLOATP (y
) ? XFLOAT_DATA (y
) : XINT (y
);
496 /* If the "remainder" comes out with the wrong sign, fix it. */
497 if (f2
< 0 ? f1
> 0 : f1
< 0)
500 return make_float (f1
);
503 DEFUN ("fceiling", Ffceiling
, Sfceiling
, 1, 1, 0,
504 doc
: /* Return the smallest integer no less than ARG, as a float.
505 \(Round toward +inf.\) */)
508 double d
= extract_float (arg
);
510 return make_float (d
);
513 DEFUN ("ffloor", Fffloor
, Sffloor
, 1, 1, 0,
514 doc
: /* Return the largest integer no greater than ARG, as a float.
515 \(Round towards -inf.\) */)
518 double d
= extract_float (arg
);
520 return make_float (d
);
523 DEFUN ("fround", Ffround
, Sfround
, 1, 1, 0,
524 doc
: /* Return the nearest integer to ARG, as a float. */)
527 double d
= extract_float (arg
);
529 return make_float (d
);
532 DEFUN ("ftruncate", Fftruncate
, Sftruncate
, 1, 1, 0,
533 doc
: /* Truncate a floating point number to an integral float value.
534 Rounds the value toward zero. */)
537 double d
= extract_float (arg
);
542 return make_float (d
);
546 syms_of_floatfns (void)
556 defsubr (&Scopysign
);
560 defsubr (&Sfceiling
);
563 defsubr (&Sftruncate
);
576 defsubr (&Struncate
);