logb now works correctly on large integers

* admin/merge-gnulib (GNULIB_MODULES): Add count-leading-zeros.
* etc/NEWS: Document the change.
* lib/count-leading-zeros.c, lib/count-leading-zeros.h:
* m4/count-leading-zeros.m4: New files, copied from Gnulib.
* lib/gnulib.mk, m4/gnulib-comp.m4: Regenerate.
* src/floatfns.c: Include count-leading-zeros.h.
(Flogb): Do not convert fixnum to float before taking the log,
as the rounding error can cause the answer to be off by 1.
* src/lisp.h (EMACS_UINT_WIDTH): New constant.
* test/src/floatfns-tests.el (logb-extreme-fixnum): New test.

10 files changed, 185 insertions, 19 deletions

diff --git a/admin/merge-gnulib b/admin/merge-gnulib
index 20a3240ea89..7e7971fe881 100755
--- a/admin/merge-gnulib
+++ b/admin/merge-gnulib
@@ -27,7 +27,8 @@ GNULIB_URL=git://git.savannah.gnu.org/gnulib.git
 
 GNULIB_MODULES='
   alloca-opt binary-io byteswap c-ctype c-strcase
-  careadlinkat close-stream count-one-bits count-trailing-zeros
+  careadlinkat close-stream
+  count-leading-zeros count-one-bits count-trailing-zeros
   crypto/md5 crypto/sha1 crypto/sha256 crypto/sha512
   dtoastr dtotimespec dup2 environ execinfo faccessat
   fcntl fcntl-h fdatasync fdopendir
diff --git a/etc/NEWS b/etc/NEWS
index 17353936e7f..a8db54c51ef 100644
--- a/etc/NEWS
+++ b/etc/NEWS
@@ -903,9 +903,9 @@ compares their numerical values.  According to this predicate,
 "foo2.png" is smaller than "foo12.png".
 
 ---
-** Numeric comparisons no longer return incorrect answers due to
-internal rounding errors.  For example, (< most-positive-fixnum (+ 1.0
-most-positive-fixnum)) now correctly returns t on 64-bit hosts.
+** Numeric comparisons and 'logb' no longer return incorrect answers
+due to internal rounding errors.  For example, (< most-positive-fixnum
+(+ 1.0 most-positive-fixnum)) now correctly returns t on 64-bit hosts.
 
 +++
 ** The new function 'char-from-name' converts a Unicode name string
diff --git a/lib/count-leading-zeros.c b/lib/count-leading-zeros.c
new file mode 100644
index 00000000000..d0c0704f582
--- /dev/null
+++ b/lib/count-leading-zeros.c
@@ -0,0 +1,3 @@
+#include <config.h>
+#define COUNT_LEADING_ZEROS_INLINE _GL_EXTERN_INLINE
+#include "count-leading-zeros.h"
diff --git a/lib/count-leading-zeros.h b/lib/count-leading-zeros.h
new file mode 100644
index 00000000000..e197137e66e
--- /dev/null
+++ b/lib/count-leading-zeros.h
@@ -0,0 +1,114 @@
+/* count-leading-zeros.h -- counts the number of leading 0 bits in a word.
+   Copyright (C) 2012-2017 Free Software Foundation, Inc.
+
+   This program is free software: you can redistribute it and/or modify
+   it under the terms of the GNU General Public License as published by
+   the Free Software Foundation; either version 3 of the License, or
+   (at your option) any later version.
+
+   This program is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+   GNU General Public License for more details.
+
+   You should have received a copy of the GNU General Public License
+   along with this program.  If not, see <http://www.gnu.org/licenses/>.  */
+
+/* Written by Eric Blake.  */
+
+#ifndef COUNT_LEADING_ZEROS_H
+#define COUNT_LEADING_ZEROS_H 1
+
+#include <limits.h>
+#include <stdlib.h>
+
+#ifndef _GL_INLINE_HEADER_BEGIN
+ #error "Please include config.h first."
+#endif
+_GL_INLINE_HEADER_BEGIN
+#ifndef COUNT_LEADING_ZEROS_INLINE
+# define COUNT_LEADING_ZEROS_INLINE _GL_INLINE
+#endif
+
+/* Assuming the GCC builtin is BUILTIN and the MSC builtin is MSC_BUILTIN,
+   expand to code that computes the number of leading zeros of the local
+   variable 'x' of type TYPE (an unsigned integer type) and return it
+   from the current function.  */
+#if __GNUC__ > 3 || (__GNUC__ == 3 && __GNUC_MINOR__ >= 4)
+# define COUNT_LEADING_ZEROS(BUILTIN, MSC_BUILTIN, TYPE)                \
+  return x ? BUILTIN (x) : CHAR_BIT * sizeof x;
+#elif _MSC_VER
+# pragma intrinsic _BitScanReverse
+# pragma intrinsic _BitScanReverse64
+# define COUNT_LEADING_ZEROS(BUILTIN, MSC_BUILTIN, TYPE)                \
+    do                                                                  \
+      {                                                                 \
+        unsigned long result;                                           \
+        return MSC_BUILTIN (&result, x) ? result : CHAR_BIT * sizeof x; \
+      }                                                                 \
+    while (0)
+#else
+# define COUNT_LEADING_ZEROS(BUILTIN, MSC_BUILTIN, TYPE)                \
+    do                                                                  \
+      {                                                                 \
+        int count;                                                      \
+        unsigned int leading_32;                                        \
+        if (! x)                                                        \
+          return CHAR_BIT * sizeof x;                                   \
+        for (count = 0;                                                 \
+             (leading_32 = ((x >> (sizeof (TYPE) * CHAR_BIT - 32))      \
+                            & 0xffffffffU),                             \
+              count < CHAR_BIT * sizeof x - 32 && !leading_32);         \
+             count += 32)                                               \
+          x = x << 31 << 1;                                             \
+        return count + count_leading_zeros_32 (leading_32);             \
+      }                                                                 \
+    while (0)
+
+/* Compute and return the number of leading zeros in X,
+   where 0 < X < 2**32.  */
+COUNT_LEADING_ZEROS_INLINE int
+count_leading_zeros_32 (unsigned int x)
+{
+  /* http://graphics.stanford.edu/~seander/bithacks.html */
+  static const char de_Bruijn_lookup[32] = {
+    31, 22, 30, 21, 18, 10, 29, 2, 20, 17, 15, 13, 9, 6, 28, 1,
+    23, 19, 11, 3, 16, 14, 7, 24, 12, 4, 8, 25, 5, 26, 27, 0
+  };
+
+  x |= x >> 1;
+  x |= x >> 2;
+  x |= x >> 4;
+  x |= x >> 8;
+  x |= x >> 16;
+  return de_Bruijn_lookup[((x * 0x07c4acddU) & 0xffffffffU) >> 27];
+}
+#endif
+
+/* Compute and return the number of leading zeros in X. */
+COUNT_LEADING_ZEROS_INLINE int
+count_leading_zeros (unsigned int x)
+{
+  COUNT_LEADING_ZEROS (__builtin_clz, _BitScanReverse, unsigned int);
+}
+
+/* Compute and return the number of leading zeros in X. */
+COUNT_LEADING_ZEROS_INLINE int
+count_leading_zeros_l (unsigned long int x)
+{
+  COUNT_LEADING_ZEROS (__builtin_clzl, _BitScanReverse, unsigned long int);
+}
+
+#if HAVE_UNSIGNED_LONG_LONG_INT
+/* Compute and return the number of leading zeros in X. */
+COUNT_LEADING_ZEROS_INLINE int
+count_leading_zeros_ll (unsigned long long int x)
+{
+  COUNT_LEADING_ZEROS (__builtin_clzll, _BitScanReverse64,
+                       unsigned long long int);
+}
+#endif
+
+_GL_INLINE_HEADER_END
+
+#endif /* COUNT_LEADING_ZEROS_H */
diff --git a/lib/gnulib.mk b/lib/gnulib.mk
index 4398fe37173..e4aa90ecac9 100644
--- a/lib/gnulib.mk
+++ b/lib/gnulib.mk
@@ -21,7 +21,7 @@
 # the same distribution terms as the rest of that program.
 #
 # Generated by gnulib-tool.
-# Reproduce by: gnulib-tool --import --lib=libgnu --source-base=lib --m4-base=m4 --doc-base=doc --tests-base=tests --aux-dir=build-aux --avoid=close --avoid=dup --avoid=fchdir --avoid=fstat --avoid=malloc-posix --avoid=msvc-inval --avoid=msvc-nothrow --avoid=open --avoid=openat-die --avoid=opendir --avoid=raise --avoid=save-cwd --avoid=select --avoid=setenv --avoid=sigprocmask --avoid=stdarg --avoid=stdbool --avoid=threadlib --avoid=unsetenv --makefile-name=gnulib.mk --conditional-dependencies --no-libtool --macro-prefix=gl --no-vc-files alloca-opt binary-io byteswap c-ctype c-strcase careadlinkat close-stream count-one-bits count-trailing-zeros crypto/md5 crypto/sha1 crypto/sha256 crypto/sha512 dtoastr dtotimespec dup2 environ execinfo faccessat fcntl fcntl-h fdatasync fdopendir filemode filevercmp flexmember fstatat fsync getloadavg getopt-gnu gettime gettimeofday gitlog-to-changelog ignore-value intprops largefile lstat manywarnings memrchr mkostemp mktime pipe2 pselect pthread_sigmask putenv qcopy-acl readlink readlinkat sig2str socklen stat-time std-gnu11 stdalign stddef stdio stpcpy strftime strtoimax strtoumax symlink sys_stat sys_time time time_r time_rz timegm timer-time timespec-add timespec-sub unsetenv update-copyright utimens vla warnings
+# Reproduce by: gnulib-tool --import --lib=libgnu --source-base=lib --m4-base=m4 --doc-base=doc --tests-base=tests --aux-dir=build-aux --avoid=close --avoid=dup --avoid=fchdir --avoid=fstat --avoid=malloc-posix --avoid=msvc-inval --avoid=msvc-nothrow --avoid=open --avoid=openat-die --avoid=opendir --avoid=raise --avoid=save-cwd --avoid=select --avoid=setenv --avoid=sigprocmask --avoid=stdarg --avoid=stdbool --avoid=threadlib --avoid=unsetenv --makefile-name=gnulib.mk --conditional-dependencies --no-libtool --macro-prefix=gl --no-vc-files alloca-opt binary-io byteswap c-ctype c-strcase careadlinkat close-stream count-leading-zeros count-one-bits count-trailing-zeros crypto/md5 crypto/sha1 crypto/sha256 crypto/sha512 dtoastr dtotimespec dup2 environ execinfo faccessat fcntl fcntl-h fdatasync fdopendir filemode filevercmp flexmember fstatat fsync getloadavg getopt-gnu gettime gettimeofday gitlog-to-changelog ignore-value intprops largefile lstat manywarnings memrchr mkostemp mktime pipe2 pselect pthread_sigmask putenv qcopy-acl readlink readlinkat sig2str socklen stat-time std-gnu11 stdalign stddef stdio stpcpy strftime strtoimax strtoumax symlink sys_stat sys_time time time_r time_rz timegm timer-time timespec-add timespec-sub unsetenv update-copyright utimens vla warnings
 
 
 MOSTLYCLEANFILES += core *.stackdump
@@ -151,6 +151,14 @@ EXTRA_DIST += close-stream.h
 
 ## end   gnulib module close-stream
 
+## begin gnulib module count-leading-zeros
+
+libgnu_a_SOURCES += count-leading-zeros.c
+
+EXTRA_DIST += count-leading-zeros.h
+
+## end   gnulib module count-leading-zeros
+
 ## begin gnulib module count-one-bits
 
 libgnu_a_SOURCES += count-one-bits.c
diff --git a/m4/count-leading-zeros.m4 b/m4/count-leading-zeros.m4
new file mode 100644
index 00000000000..b6e97cf13f8
--- /dev/null
+++ b/m4/count-leading-zeros.m4
@@ -0,0 +1,12 @@
+# count-leading-zeros.m4 serial 2
+dnl Copyright (C) 2012-2017 Free Software Foundation, Inc.
+dnl This file is free software; the Free Software Foundation
+dnl gives unlimited permission to copy and/or distribute it,
+dnl with or without modifications, as long as this notice is preserved.
+
+AC_DEFUN([gl_COUNT_LEADING_ZEROS],
+[
+  dnl We don't need (and can't compile) count_leading_zeros_ll
+  dnl unless the type 'unsigned long long int' exists.
+  AC_REQUIRE([AC_TYPE_UNSIGNED_LONG_LONG_INT])
+])
diff --git a/m4/gnulib-comp.m4 b/m4/gnulib-comp.m4
index 14af9fbd6bf..a3e30fd736e 100644
--- a/m4/gnulib-comp.m4
+++ b/m4/gnulib-comp.m4
@@ -54,6 +54,7 @@ AC_DEFUN([gl_EARLY],
   # Code from module careadlinkat:
   # Code from module clock-time:
   # Code from module close-stream:
+  # Code from module count-leading-zeros:
   # Code from module count-one-bits:
   # Code from module count-trailing-zeros:
   # Code from module crypto/md5:
@@ -190,6 +191,7 @@ AC_DEFUN([gl_INIT],
   gl_CLOCK_TIME
   gl_CLOSE_STREAM
   gl_MODULE_INDICATOR([close-stream])
+  gl_COUNT_LEADING_ZEROS
   gl_COUNT_ONE_BITS
   gl_COUNT_TRAILING_ZEROS
   gl_MD5
@@ -871,6 +873,8 @@ AC_DEFUN([gl_FILE_LIST], [
   lib/careadlinkat.h
   lib/close-stream.c
   lib/close-stream.h
+  lib/count-leading-zeros.c
+  lib/count-leading-zeros.h
   lib/count-one-bits.c
   lib/count-one-bits.h
   lib/count-trailing-zeros.c
@@ -1000,6 +1004,7 @@ AC_DEFUN([gl_FILE_LIST], [
   m4/c-strtod.m4
   m4/clock_time.m4
   m4/close-stream.m4
+  m4/count-leading-zeros.m4
   m4/count-one-bits.m4
   m4/count-trailing-zeros.m4
   m4/dirent_h.m4
diff --git a/src/floatfns.c b/src/floatfns.c
index dda03698093..94da22a3ba7 100644
--- a/src/floatfns.c
+++ b/src/floatfns.c
@@ -45,6 +45,8 @@ along with GNU Emacs.  If not, see <http://www.gnu.org/licenses/>.  */
 
 #include <math.h>
 
+#include <count-leading-zeros.h>
+
 /* 'isfinite' and 'isnan' cause build failures on Solaris 10 with the
    bundled GCC in c99 mode.  Work around the bugs with simple
    implementations that are good enough.  */
@@ -290,28 +292,46 @@ DEFUN ("float", Ffloat, Sfloat, 1, 1, 0,
     return arg;
 }
 
+static int
+ecount_leading_zeros (EMACS_UINT x)
+{
+  return (EMACS_UINT_WIDTH == UINT_WIDTH ? count_leading_zeros (x)
+	  : EMACS_UINT_WIDTH == ULONG_WIDTH ? count_leading_zeros_l (x)
+	  : count_leading_zeros_ll (x));
+}
+
 DEFUN ("logb", Flogb, Slogb, 1, 1, 0,
        doc: /* Returns largest integer <= the base 2 log of the magnitude of ARG.
 This is the same as the exponent of a float.  */)
   (Lisp_Object arg)
 {
-  Lisp_Object val;
   EMACS_INT value;
-  double f = extract_float (arg);
+  CHECK_NUMBER_OR_FLOAT (arg);
 
-  if (f == 0.0)
-    value = MOST_NEGATIVE_FIXNUM;
-  else if (isfinite (f))
+  if (FLOATP (arg))
     {
-      int ivalue;
-      frexp (f, &ivalue);
-      value = ivalue - 1;
+      double f = XFLOAT_DATA (arg);
+
+      if (f == 0)
+	value = MOST_NEGATIVE_FIXNUM;
+      else if (isfinite (f))
+	{
+	  int ivalue;
+	  frexp (f, &ivalue);
+	  value = ivalue - 1;
+	}
+      else
+	value = MOST_POSITIVE_FIXNUM;
     }
   else
-    value = MOST_POSITIVE_FIXNUM;
+    {
+      EMACS_INT i = eabs (XINT (arg));
+      value = (i == 0
+	       ? MOST_NEGATIVE_FIXNUM
+	       : EMACS_UINT_WIDTH - 1 - ecount_leading_zeros (i));
+    }
 
-  XSETINT (val, value);
-  return val;
+  return make_number (value);
 }
 
 
diff --git a/src/lisp.h b/src/lisp.h
index 220188cdb87..6d0b5283356 100644
--- a/src/lisp.h
+++ b/src/lisp.h
@@ -80,19 +80,19 @@ DEFINE_GDB_SYMBOL_END (GCTYPEBITS)
 # elif INTPTR_MAX <= INT_MAX && !defined WIDE_EMACS_INT
 typedef int EMACS_INT;
 typedef unsigned int EMACS_UINT;
-enum { EMACS_INT_WIDTH = INT_WIDTH };
+enum { EMACS_INT_WIDTH = INT_WIDTH, EMACS_UINT_WIDTH = UINT_WIDTH };
 #  define EMACS_INT_MAX INT_MAX
 #  define pI ""
 # elif INTPTR_MAX <= LONG_MAX && !defined WIDE_EMACS_INT
 typedef long int EMACS_INT;
 typedef unsigned long EMACS_UINT;
-enum { EMACS_INT_WIDTH = LONG_WIDTH };
+enum { EMACS_INT_WIDTH = LONG_WIDTH, EMACS_UINT_WIDTH = ULONG_WIDTH };
 #  define EMACS_INT_MAX LONG_MAX
 #  define pI "l"
 # elif INTPTR_MAX <= LLONG_MAX
 typedef long long int EMACS_INT;
 typedef unsigned long long int EMACS_UINT;
-enum { EMACS_INT_WIDTH = LLONG_WIDTH };
+enum { EMACS_INT_WIDTH = LLONG_WIDTH, EMACS_UINT_WIDTH = ULLONG_WIDTH };
 #  define EMACS_INT_MAX LLONG_MAX
 #  ifdef __MINGW32__
 #   define pI "I64"
diff --git a/test/src/floatfns-tests.el b/test/src/floatfns-tests.el
index a2116a59459..448d6167f25 100644
--- a/test/src/floatfns-tests.el
+++ b/test/src/floatfns-tests.el
@@ -25,4 +25,7 @@
   (should-error (round most-negative-fixnum -1.0))
   (should-error (truncate most-negative-fixnum -1.0)))
 
+(ert-deftest logb-extreme-fixnum ()
+  (should (= (logb most-negative-fixnum) (1+ (logb most-positive-fixnum)))))
+
 (provide 'floatfns-tests)