Hacker News new | ask | show | jobs
by pkhuong 4010 days ago
For values far from 1, log_2 will be dominated by the integer portion, which IEEE FP representation gives you for free, so it's not surprising that outputs eventually look right. I think the problem is that the approximation works with (x - 1), but only the "significand >= 1.5" branch subtracts that offset. Really, if you need speed, I'd probably start with an approximation over [1, 2)...

Re-implementation in SBCL:

  CL-USER> (defun approx (x)
             (let ((o (1- x)))
               (values (/ (+ (* o o 0.338953) (* o 2.198599)) (+ o 1.523692))
                       (log x 2.0))))

  CL-USER> (loop for x from 0.70 upto 1.51 by 0.1
                 do (format t "~,3F: ~{~F ~F~}~%" x (multiple-value-list (approx x))))
  0.700: -0.51407874 -0.5145732
  0.800: -0.32194924 -0.32192805
  0.900: -0.15204856 -0.15200303
  1.000: 0.0 0.0
  1.100: 0.13749498 0.13750356
  1.200: 0.26296926 0.26303446
  1.300: 0.3784003 0.37851173
  1.400: 0.48535433 0.48542693
  1.500: 0.585088 0.5849626
That's in the right ballpark for [0,7, 1.5]

  CL-USER> (defun almost-log2 (x)
             (let* ((bits (sb-kernel:single-float-bits x))
                    (exp (ash bits -23)))
               (values 
                (if (logbitp 22 bits)
                    (+ (- exp 126) (approx (sb-kernel:make-single-float
                                            (dpb 126 (byte 9 23) bits))))
                    (+ (- exp 127) (approx (sb-kernel:make-single-float
                                            (dpb 127 (byte 9 23) bits)))))
                (log x 2.0))))

  CL-USER> (loop for x from 0.1 upto 10.1 by 0.1
                 do (format t "~,3F: ~{~,8F ~,8F~}~%" x (multiple-value-list (almost-log2 x))))
  0.100: -3.32194920 -3.32192800
  0.200: -2.32194920 -2.32192800
  0.300: -1.73703070 -1.73696570
  0.400: -1.32194920 -1.32192800
  0.500: -1.00000000 -1.00000000
  0.600: -0.73703074 -0.73696554
  0.700: -0.51464570 -0.51457310
  0.800: -0.32194915 -0.32192796
  0.900: -0.15204845 -0.15200295
  1.000: 0.00000017 0.00000017
  1.100: 0.13749513 0.13750371
  1.200: 0.26296943 0.26303460
  1.300: 0.37840047 0.37851185
  1.400: 0.48535448 0.48542705
  1.500: 0.58509207 0.58496270
  1.600: 0.67805094 0.67807215
  1.700: 0.76547647 0.76553500
  1.800: 0.84795165 0.84799720
  1.900: 0.92598030 0.92599964
  2.000: 1.00000010 1.00000010
  2.100: 1.07039330 1.07038940
  2.200: 1.13749500 1.13750360
  [...]
  2.800: 1.48535400 1.48542650
  2.900: 1.53605470 1.53605260
  3.000: 1.58508750 1.58496210
  3.100: 1.63230250 1.63226800
  3.200: 1.67805030 1.67807150
  [...]
  3.800: 1.92597950 1.92599880
  3.900: 1.96346550 1.96347340
  4.000: 1.99999940 1.99999930
  4.100: 2.03562740 2.03562330
  4.200: 2.07039260 2.07038880
  [...]
  5.000: 2.32183340 2.32192750
  [...]
  6.000: 2.58508730 2.58496170
  [...]
  9.900: 3.30734100 3.30742860
  10.000: 3.32183430 3.32192850
1 comments

in [part I](http://www.ebaytechblog.com/2015/05/01/fast-approximate-loga...) the author argues against using [1,2) because it causes large relative error near x=1 (since log x = 0 there). I would disagree (why would you ever care about relative error in a logarithm -- IMHO absolute error makes much more sense... unless for some crazy reason you're dividing by log x near x=1) but that's his argument.
It's an interesting trade off; you get a nice bump in precision for a small increase in code complexity.