Raymond Toy pushed to branch master at cmucl / cmucl Commits: 209af323 by Raymond Toy at 2026-06-30T18:08:00-07:00 Address #435: Integrate single-float core-math into lisp - - - - - 7be88418 by Raymond Toy at 2026-06-30T18:08:00-07:00 Merge branch 'issue-435-add-core-math-lisp-support' into 'master' Address #435: Integrate single-float core-math into lisp See merge request cmucl/cmucl!357 - - - - - 7 changed files: - src/code/exports.lisp - src/code/irrat.lisp - src/compiler/float-tran.lisp - src/lisp/fdlibm.h - src/lisp/irrat.c - src/lisp/log2.c - src/lisp/openlibm/e_expf.c Changes: ===================================== src/code/exports.lisp ===================================== @@ -2191,10 +2191,33 @@ "%IEEE754-REM-PI/2" "%SINCOS" - "STANDARD-CHAR-TYPE" "MAKE-STANDARD-CHAR-TYPE" - "STANDARD-CHAR-TYPE-P") + "STANDARD-CHAR-TYPE-P" + + ;; Single-float functions + "%ACOSF" + "%ACOSHF" + "%ASINF" + "%ASINHF" + "%ATAN2F" + "%ATANF" + "%ATANHF" + "%COSF" + "%COSHF" + "%EXPF" + "%EXPM1F" + "%HYPOTF" + "%LOG10F" + "%LOG1PF" + "%LOG2F" + "%LOGF" + "%POWF" + "%SINF" + "%SINHF" + "%TANF" + "%TANHF" + ) #+heap-overflow-check (:export "DYNAMIC-SPACE-OVERFLOW-WARNING-HIT" "DYNAMIC-SPACE-OVERFLOW-ERROR-HIT" ===================================== src/code/irrat.lisp ===================================== @@ -40,28 +40,37 @@ ;;; call, and saves number consing to boot. ;;; (defmacro def-math-rtn (name num-args) - (multiple-value-bind (c-name lisp-name) - (if (listp name) - (values (first name) (second name)) - (values name - (intern (concatenate 'simple-string - "%" - (string-upcase name))))) - `(progn - (declaim (inline ,lisp-name)) - (export ',lisp-name) - (alien:def-alien-routine (,c-name ,lisp-name) double-float - ,@(let ((results nil)) - (dotimes (i num-args (nreverse results)) - (push (list (intern (format nil "ARG-~D" i)) - 'double-float) - results))))))) + (flet ((def-rtn (c-name lisp-name type) + `(progn + (declaim (inline ,lisp-name)) + (export ',lisp-name) + (alien:def-alien-routine (,c-name ,lisp-name) ,type + ,@(let ((results nil)) + (dotimes (i num-args (nreverse results)) + (push (list (intern (format nil "ARG-~D" i)) + type) + results))))))) + (multiple-value-bind (c-name lisp-name) + (if (listp name) + (values (first name) (second name)) + (values name + (intern (concatenate 'simple-string + "%" + (string-upcase name))))) + `(progn + ,(def-rtn c-name lisp-name 'double-float) + ;; The C99 convention is for the float version to have the same + ;; name as the double version but with an "f" appended. We do + ;; the same here and append an "F" to the lisp name too. + ,(def-rtn (concatenate 'string c-name "f") + (symbolicate lisp-name "F") + 'single-float))))) (eval-when (compile load eval) (defun handle-reals (function var) `((((foreach fixnum single-float bignum ratio)) - (coerce (,function (coerce ,var 'double-float)) 'single-float)) + (,(symbolicate function "F") (coerce ,var 'single-float))) ((double-float) (,function ,var)) #+double-double @@ -213,21 +222,34 @@ (y0 double-float :out) (y1 double-float :out)) -(declaim (inline %%sincos)) +(declaim (inline %%sincos %%sincosf)) (alien:def-alien-routine ("lisp_sincos" %%sincos) c-call:void (x double-float) (s double-float :out) (c double-float :out)) +(alien:def-alien-routine ("lisp_sincosf" %%sincosf) + c-call:void + (x single-float) + (s single-float :out) + (c single-float :out)) + (declaim (inline %sincos)) (export '%sincos) (defun %sincos (x) - (declare (double-float x)) - (multiple-value-bind (ign s c) - (%%sincos x) - (declare (ignore ign)) - (values s c))) + (declare (float x)) + (etypecase x + (single-float + (multiple-value-bind (ign s c) + (%%sincosf x) + (declare (ignore ign)) + (values s c))) + (double-float + (multiple-value-bind (ign s c) + (%%sincos x) + (declare (ignore ign)) + (values s c)))))  ;;;; Power functions. @@ -993,17 +1015,17 @@ (if (complexp theta) (error (intl:gettext "Argument to CIS is complex: ~S") theta) (number-dispatch ((theta real)) - ((rational) - (let ((arg (coerce theta 'double-float))) - (multiple-value-bind (s c) - (%sincos arg) - (complex (coerce c 'single-float) - (coerce s 'single-float))))) - (((foreach single-float double-float)) - (multiple-value-bind (s c) - (%sincos (coerce theta 'double-float)) - (complex (coerce c '(dispatch-type theta)) - (coerce s '(dispatch-type theta))))) + (((foreach rational single-float)) + (let ((arg (coerce theta 'single-float))) + (multiple-value-bind (ign s c) + (%%sincosf arg) + (declare (ignore ign)) + (complex c s)))) + ((double-float) + (multiple-value-bind (ign s c) + (%%sincos theta) + (declare (ignore ign)) + (complex c s))) #+double-double ((double-double-float) (multiple-value-bind (s c) @@ -1016,14 +1038,19 @@ ((rational) (if (or (> number 1) (< number -1)) (complex-asin number) - (coerce (%asin (coerce number 'double-float)) 'single-float))) - (((foreach single-float double-float)) + (%asinf (coerce number 'single-float)))) + ((single-float) (if (or (float-nan-p number) (and (<= number (coerce 1 '(dispatch-type number))) (>= number (coerce -1 '(dispatch-type number))))) - (coerce (%asin (coerce number 'double-float)) - '(dispatch-type number)) + (%asinf number) (complex-asin number))) + ((double-float) + (if (or (float-nan-p number) + (and (<= number (coerce 1 '(dispatch-type number))) + (>= number (coerce -1 '(dispatch-type number))))) + (%asin number)) + (complex-asin number)) #+double-double ((double-double-float) (if (or (float-nan-p number) @@ -1040,13 +1067,18 @@ ((rational) (if (or (> number 1) (< number -1)) (complex-acos number) - (coerce (%acos (coerce number 'double-float)) 'single-float))) - (((foreach single-float double-float)) + (%acosf (coerce number 'single-float)))) + ((single-float) (if (or (float-nan-p number) (and (<= number (coerce 1 '(dispatch-type number))) (>= number (coerce -1 '(dispatch-type number))))) - (coerce (%acos (coerce number 'double-float)) - '(dispatch-type number)) + (%acosf number) + (complex-acos number))) + ((double-float) + (if (or (float-nan-p number) + (and (<= number (coerce 1 '(dispatch-type number))) + (>= number (coerce -1 '(dispatch-type number))))) + (%acos number) (complex-acos number))) #+double-double ((double-double-float) @@ -1138,12 +1170,15 @@ ;; acosh is complex if number < 1 (if (< number 1) (complex-acosh number) - (coerce (%acosh (coerce number 'double-float)) 'single-float))) - (((foreach single-float double-float)) + (%acoshf (coerce number 'single-float)))) + ((single-float) (if (< number (coerce 1 '(dispatch-type number))) (complex-acosh number) - (coerce (%acosh (coerce number 'double-float)) - '(dispatch-type number)))) + (%acoshf number))) + ((double-float) + (if (< number (coerce 1 '(dispatch-type number))) + (complex-acosh number) + (%acosh number))) #+double-double ((double-double-float) (if (< number 1w0) @@ -1159,13 +1194,17 @@ ;; atanh is complex if |number| > 1 (if (or (> number 1) (< number -1)) (complex-atanh number) - (coerce (%atanh (coerce number 'double-float)) 'single-float))) - (((foreach single-float double-float)) + (%atanhf (coerce number 'single-float)))) + ((single-float) (if (or (> number (coerce 1 '(dispatch-type number))) (< number (coerce -1 '(dispatch-type number)))) (complex-atanh number) - (coerce (%atanh (coerce number 'double-float)) - '(dispatch-type number)))) + (%atanhf number))) + ((double-float) + (if (or (> number (coerce 1 '(dispatch-type number))) + (< number (coerce -1 '(dispatch-type number)))) + (complex-atanh number) + (%atanh number))) #+double-double ((double-double-float) (if (or (> number 1w0) ===================================== src/compiler/float-tran.lisp ===================================== @@ -683,25 +683,41 @@ (asinh %asinh *) (acosh %acosh float) (atanh %atanh float))) + (destructuring-bind (name prim rtype) + stuff + (let ((primf (intern (string (symbolicate prim "F")) + "KERNEL"))) + (deftransform name ((x) '(single-float) rtype :eval-name t) + `(,primf x)) + (deftransform name ((x) '(double-float) rtype :eval-name t :when :both) + `(,prim x))))) + +#-core-math +(dolist (stuff '((sqrt %sqrt float))) (destructuring-bind (name prim rtype) stuff (deftransform name ((x) '(single-float) rtype :eval-name t) `(coerce (,prim (coerce x 'double-float)) 'single-float)) (deftransform name ((x) '(double-float) rtype :eval-name t :when :both) - `(,prim x)))) + `(,prim x)))) -(defknown (%sincos) - (double-float) (values double-float double-float) +(defknown (kernel::%%sincos) + (double-float) (values null double-float double-float) (movable foldable flushable)) -(deftransform cis ((x) (single-float) * :when :both) - `(multiple-value-bind (s c) - (%sincos (coerce x 'double-float)) - (complex (coerce c 'single-float) - (coerce s 'single-float)))) +(defknown (kernel::%%sincosf) + (single-float) (values null single-float single-float) + (movable foldable flushable)) (deftransform cis ((x) (double-float) * :when :both) - `(multiple-value-bind (s c) - (%sincos x) + `(multiple-value-bind (ign s c) + (kernel::%%sincos x) + (declare (ignore ign)) + (complex c s))) + +(deftransform cis ((x) (single-float) * :when :both) + `(multiple-value-bind (ign s c) + (kernel::%%sincosf x) + (declare (ignore ign)) (complex c s))) ;;; The argument range is limited on the x86 FP trig. functions. A ===================================== src/lisp/fdlibm.h ===================================== @@ -69,6 +69,7 @@ extern double __ieee754_log10(double x); extern double __ieee754_pow(double x, double y); extern double __ieee754_hypot(double x, double y); extern double fdlibm_scalbn(double x, int n); +extern double fdlibm_log2(double x); enum FDLIBM_EXCEPTION { FDLIBM_DIVIDE_BY_ZERO, ===================================== src/lisp/irrat.c ===================================== @@ -32,6 +32,7 @@ extern double cr_hypot(double, double); extern double cr_log1p(double); extern double cr_expm1(double); extern void cr_sincos(double, double *, double *); +extern double cr_log2(double); extern float cr_sinf(float); extern float cr_cosf(float); @@ -54,6 +55,7 @@ extern float cr_hypotf(float, float); extern float cr_log1pf(float); extern float cr_expm1f(float); extern void cr_sincosf(float, float *, float *); +extern float cr_log2f(float); #else #include "openlibm_math.h" /* @@ -72,6 +74,7 @@ extern void openlibm_sincosf(float, float*, float*); extern float __ieee754_acosf(float); extern float __ieee754_acoshf(float); extern float __ieee754_logf(float); +extern float __ieee754_log2f(float); extern float __ieee754_atanhf(float); extern float __ieee754_asinf(float); extern float __ieee754_atan2f(float, float); @@ -102,7 +105,7 @@ extern float __ieee754_hypotf(float, float); return fdlibm_setexception(x, FDLIBM_OVERFLOW); \ } \ } while (0) - + double lisp_sin(double x) @@ -363,6 +366,16 @@ lisp_sincos(double x, double *s, double *c) #endif } +double +lisp_log2(double x) +{ +#ifdef FEATURE_CORE_MATH + return cr_log2(x); +#else + return fdlibm_log2(x); +#endif +} + /* * The single-float versions of the special functions. When core-math * is set, we use the single-float core-math functions. Otherwise, we @@ -592,3 +605,13 @@ lisp_sincosf(float x, float *s, float *c) openlibm_sincosf(x, s, c); #endif } + +float +lisp_log2f(float x) +{ +#ifdef FEATURE_CORE_MATH + return cr_log2f(x); +#else + return __ieee754_log2f(x); +#endif +} ===================================== src/lisp/log2.c ===================================== @@ -39,7 +39,7 @@ cp = 9.61796693925975554329e-01, /* 0x3FEEC709, 0xDC3A03FD =2/(3ln2) */ cp_h = 9.61796700954437255859e-01, /* 0x3FEEC709, 0xE0000000 =(float)cp */ cp_l = -7.02846165095275826516e-09; /* 0xBE3E2FE0, 0x145B01F5 =tail of cp_h*/ -double lisp_log2(double x) +double fdlibm_log2(double x) { double ax; int k, hx, lx, ix; ===================================== src/lisp/openlibm/e_expf.c ===================================== @@ -58,7 +58,15 @@ __ieee754_expf(float x) return x+x; /* NaN */ if(hx==0x7f800000) return (xsb==0)? x:0.0; /* exp(+-inf)={inf,0} */ - if(x > o_threshold) return huge*huge; /* overflow */ + /* + * Signal overflow. This keeps the compiler from + * constant-folding the result that then no longer signals + * an overflow. + */ + if(x > o_threshold) { + volatile float fhuge = huge; + return fhuge*fhuge; /* overflow */ + } if(x < u_threshold) return twom100*twom100; /* underflow */ } View it on GitLab: https://gitlab.common-lisp.net/cmucl/cmucl/-/compare/b767e8f7fb0e35c82a9b445... -- View it on GitLab: https://gitlab.common-lisp.net/cmucl/cmucl/-/compare/b767e8f7fb0e35c82a9b445... You're receiving this email because of your account on gitlab.common-lisp.net. Manage all notifications: https://gitlab.common-lisp.net/-/profile/notifications | Help: https://gitlab.common-lisp.net/help