(gitlab is uusable; i have to report here.) ecl 21.2.1 gives:
(log 1/6319748715279270675921934218987893281199411530039296)
Debugger received error of type: DIVISION-BY-ZERO #<a DIVISION-BY-ZERO 0x7ff5afe6c540> Error flushed. whereas other cl's (i testyed sbcl and ccl) give results like: ? (log 1/6319748715279270675921934218987893281199411530039296) -119.27552 I tested on amd64 (gentoo) and arm64 (debian and netbsd) with identical results. i did not have a musl box or other bsd to test on. run with --no-trap-fpe, the result is #<single-float negative infinity>. another example is: (truncate (log 1/6319748715279270675921934218987893281199418867)) Debugger received error of type: ARITHMETIC-ERROR #<a ARITHMETIC-ERROR 0xffff8f688a40> whereas this works: (truncate (log 1/631974871527927067592193421898789328119941867)) -103 -0.27893066 which of course suggests that the issue is c's long double's precision. it looks like ecl could use an mpq log function; https://github.com/linas/anant might work. -JimC -- James Cloos <cloos@jhcloos.com> OpenPGP: 0x997A9F17ED7DAEA6
This is permitted by http://www.lispworks.com/documentation/HyperSpec/Body/12_acc.htm On Sun, Jul 3, 2022 at 12:10 AM James Cloos <cloos@jhcloos.com> wrote:
(gitlab is uusable; i have to report here.)
ecl 21.2.1 gives:
(log 1/6319748715279270675921934218987893281199411530039296)
Debugger received error of type: DIVISION-BY-ZERO #<a DIVISION-BY-ZERO 0x7ff5afe6c540> Error flushed.
whereas other cl's (i testyed sbcl and ccl) give results like:
? (log 1/6319748715279270675921934218987893281199411530039296) -119.27552
I tested on amd64 (gentoo) and arm64 (debian and netbsd) with identical results. i did not have a musl box or other bsd to test on.
run with --no-trap-fpe, the result is #<single-float negative infinity>.
another example is:
(truncate (log 1/6319748715279270675921934218987893281199418867))
Debugger received error of type: ARITHMETIC-ERROR #<a ARITHMETIC-ERROR 0xffff8f688a40>
whereas this works:
(truncate (log 1/631974871527927067592193421898789328119941867))
-103 -0.27893066
which of course suggests that the issue is c's long double's precision.
it looks like ecl could use an mpq log function; https://github.com/linas/anant might work.
-JimC -- James Cloos <cloos@jhcloos.com> OpenPGP: 0x997A9F17ED7DAEA6
Out of curiosity -- the spec mentions "floating point approximations", and in this case I understand it is reasonable (a very small number is treated like zero). But is floating-point infinity also considered an approximation to very large numbers? (setq x (expt 10 5000)) (log x) => 11512.926 ECL, SBCL, Clisp, CCL, CMUCL return 11512.926. ABCL doesn't (it actually invokes a TYPE-ERROR condition, because the number is too large to be converted to DOUBLE-FLOAT) -- but it does not return float infinity. GCL's behavior is a bit more complicated: (setq x (expt 10 5000)) (progn (setq y (log x)) 'whatever) ;; don't let GCL print it yet! => WHATEVER (type-of y) => LONG-FLOAT But If I ask for the value of y on the REPL: Error: Fast links are on: do (si::use-fast-links nil) for debugging Signalled by SYSTEM::GCL-TOP-LEVEL. Condition in SYSTEM::GCL-TOP-LEVEL [or a callee]: INTERNAL-SIMPLE-ERROR: Can't print a non-number. (So I produced a non-printable LONG-FLOAT object by just calling `log` on an integer.) I did not test Clasp. Anyway - I was just wondering how to interpret the standard in this respect. J. Stas Boukarev <stassats@gmail.com> writes:
This is permitted by http://www.lispworks.com/documentation/HyperSpec/Body/12_acc.htm
On Sun, Jul 3, 2022 at 12:10 AM James Cloos <cloos@jhcloos.com> wrote:
(gitlab is uusable; i have to report here.)
ecl 21.2.1 gives:
(log 1/6319748715279270675921934218987893281199411530039296)
Debugger received error of type: DIVISION-BY-ZERO #<a DIVISION-BY-ZERO 0x7ff5afe6c540> Error flushed.
whereas other cl's (i testyed sbcl and ccl) give results like:
? (log 1/6319748715279270675921934218987893281199411530039296) -119.27552
I tested on amd64 (gentoo) and arm64 (debian and netbsd) with identical results. i did not have a musl box or other bsd to test on.
run with --no-trap-fpe, the result is #<single-float negative infinity>.
another example is:
(truncate (log 1/6319748715279270675921934218987893281199418867))
Debugger received error of type: ARITHMETIC-ERROR #<a ARITHMETIC-ERROR 0xffff8f688a40>
whereas this works:
(truncate (log 1/631974871527927067592193421898789328119941867))
-103 -0.27893066
which of course suggests that the issue is c's long double's precision.
it looks like ecl could use an mpq log function; https://github.com/linas/anant might work.
-JimC -- James Cloos <cloos@jhcloos.com> OpenPGP: 0x997A9F17ED7DAEA6
it looks like that bug only harms ratio. bignum integer seems to be ok. so aworkarount is, in the caseof a ratio, to take the difference of the logs of the numerator and the denominator. -JimC -- James Cloos <cloos@jhcloos.com> OpenPGP: 0x997A9F17ED7DAEA6
On Sun, Jul 3, 2022 at 8:10 AM James Cloos <cloos@jhcloos.com> wrote:
it looks like that bug only harms ratio. bignum integer seems to be ok.
so aworkarount is, in the caseof a ratio, to take the difference of the logs of the numerator and the denominator.
This is basically what cmucl does. But you have to be careful if the numerator and denominator are very close to each other. Then you'll take the difference of two logs that are essentially equal and lose lots or precision. IIRC, cmucl doesn't do this if the ratio is close to 1 because then converting to double-float isn't a problem. There are other ways to handle this though, like implementing a log2 function that returns the integer part and the fraction separately so you don't lose too many fraction bits when the numbers are huge.
-JimC -- James Cloos <cloos@jhcloos.com> OpenPGP: 0x997A9F17ED7DAEA6
-- Ray
i should have realized this even before noticing ecl's limitation for logs of ratios, but at least i finally did.... since i'monly using log to find the exponent for an arbitrary precision replacement for format's ~e, it works to truncate first. so i ended up with this: (defun ilog (n &optional (base nil base-p)) "return floor of base b log of n" (let ((s (if (< n 1) -1 1)) (trunc-n (truncate (if (< n 1) (/ n) n)))) (floor (* s (if base-p (log trunc-n base) (log trunc-n)))))) One could use let* and then test s when choosing between n and (/ n), but i choose instead to let the compiler optimize the two calls to (< n 1). that makes my ~e replacement work exactly as desired on all lisps i've tried. (I have a ~// function which chooses output baed on the arg's type, with ratios done via that print-e-significant-digits function and everything else using write.) -JimC -- James Cloos <cloos@jhcloos.com> OpenPGP: 0x997A9F17ED7DAEA6
Hey, sorry for taking my time with the response. While indeed conforming, I've added a separate code path for ratios. I'm conditioning on the integer length (like SBCL does). static cl_object ecl_log1_ratio(cl_object x) { cl_object num = x->ratio.num; cl_object den = x->ratio.den; if (ecl_integer_length(num) != ecl_integer_length(den)) { return ecl_minus(ecl_log1(num), ecl_log1(den)); } else { return ecl_log1_fixnum(x); } } Thanks all for the input. The change is now submitted as a merge request against develop (and will be part of the next release). Best regards, Daniel -- Daniel Kochmański ;; aka jackdaniel | Przemyśl, Poland TurtleWare - Daniel Kochmański | www.turtleware.eu "Be the change that you wish to see in the world." - Mahatma Gandhi ------- Original Message ------- On Sunday, July 10th, 2022 at 5:39 PM, James Cloos <cloos@jhcloos.com> wrote:
i should have realized this even before noticing ecl's limitation for logs of ratios, but at least i finally did....
since i'monly using log to find the exponent for an arbitrary precision replacement for format's ~e, it works to truncate first.
so i ended up with this:
(defun ilog (n &optional (base nil base-p)) "return floor of base b log of n" (let ((s (if (< n 1) -1 1)) (trunc-n (truncate (if (< n 1) (/ n) n)))) (floor (* s (if base-p (log trunc-n base) (log trunc-n))))))
One could use let* and then test s when choosing between n and (/ n), but i choose instead to let the compiler optimize the two calls to (< n 1).
that makes my ~e replacement work exactly as desired on all lisps i've tried.
(I have a ~// function which chooses output baed on the arg's type, with ratios done via that print-e-significant-digits function and everything else using write.)
-JimC -- James Cloos cloos@jhcloos.com OpenPGP: 0x997A9F17ED7DAEA6
participants (5)
-
Daniel Kochmański -
James Cloos -
Jeronimo Pellegrini -
Raymond Toy -
Stas Boukarev