lib/integer_length.c

   1 /* integer_length - find most significant bit in an 'unsigned int'.
   2    Copyright (C) 2011-2020 Free Software Foundation, Inc.
   3
   4    This program is free software: you can redistribute it and/or modify
   5    it under the terms of the GNU General Public License as published by
   6    the Free Software Foundation; either version 3 of the License, or
   7    (at your option) any later version.
   8
   9    This program is distributed in the hope that it will be useful,
  10    but WITHOUT ANY WARRANTY; without even the implied warranty of
  11    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  12    GNU General Public License for more details.
  13
  14    You should have received a copy of the GNU General Public License
  15    along with this program.  If not, see <https://www.gnu.org/licenses/>.  */
  16
  17 /* Written by Bruno Haible <bruno@clisp.org>, 2011.  */
  18
  19 #include <config.h>
  20
  21 /* Specification.  */
  22 #include "integer_length.h"
  23
  24 #include <limits.h>
  25
  26 #include "float+.h"
  27
  28 #if defined _MSC_VER && !(__clang_major__ >= 4)
  29 # include <intrin.h>
  30 #endif
  31
  32 #define NBITS (sizeof (unsigned int) * CHAR_BIT)
  33
  34 int
  35 integer_length (unsigned int x)
  36 {
  37 #if __GNUC__ > 3 || (__GNUC__ == 3 && __GNUC_MINOR__ >= 4) || (__clang_major__ >= 4)
  38   if (x == 0)
  39     return 0;
  40   else
  41     return NBITS - __builtin_clz (x);
  42 #elif defined _MSC_VER
  43   /* _BitScanReverse
  44      <https://docs.microsoft.com/en-us/cpp/intrinsics/bitscanreverse-bitscanreverse64> */
  45   unsigned long bit;
  46   if (_BitScanReverse (&bit, x))
  47     return bit + 1;
  48   else
  49     return 0;
  50 #else
  51 # if defined DBL_EXPBIT0_WORD && defined DBL_EXPBIT0_BIT
  52   if (NBITS <= DBL_MANT_BIT)
  53     {
  54       /* Use 'double' operations.
  55          Assumes an IEEE 754 'double' implementation.  */
  56 #  define DBL_EXP_MASK ((DBL_MAX_EXP - DBL_MIN_EXP) | 7)
  57 #  define DBL_EXP_BIAS (DBL_EXP_MASK / 2 - 1)
  58 #  define NWORDS \
  59     ((sizeof (double) + sizeof (unsigned int) - 1) / sizeof (unsigned int))
  60       typedef union { double value; unsigned int word[NWORDS]; }
  61               memory_double;
  62
  63       if (x == 0)
  64         return 0;
  65       else
  66         {
  67           memory_double m;
  68           unsigned int exponent;
  69
  70           if (1)
  71             {
  72               /* Use a single integer to floating-point conversion.  */
  73               m.value = x;
  74             }
  75           else
  76             {
  77               /* Use a single floating-point subtraction.  */
  78               /* 2^(DBL_MANT_DIG-1).  */
  79               static const double TWO_DBL_MANT_DIG =
  80                 /* Assume DBL_MANT_DIG <= 5 * 31.
  81                    Use the identity
  82                    n = floor(n/5) + floor((n+1)/5) + ... + floor((n+4)/5).  */
  83                 (double) (1U << ((DBL_MANT_DIG - 1) / 5))
  84                 * (double) (1U << ((DBL_MANT_DIG - 1 + 1) / 5))
  85                 * (double) (1U << ((DBL_MANT_DIG - 1 + 2) / 5))
  86                 * (double) (1U << ((DBL_MANT_DIG - 1 + 3) / 5))
  87                 * (double) (1U << ((DBL_MANT_DIG - 1 + 4) / 5));
  88
  89               /* Construct 2^(DBL_MANT_DIG-1) + x by hand.  */
  90               m.word[DBL_EXPBIT0_WORD] =
  91                 (DBL_MANT_DIG + DBL_EXP_BIAS) << DBL_EXPBIT0_BIT;
  92               m.word[1 - DBL_EXPBIT0_WORD] = x;
  93
  94               /* Subtract 2^(DBL_MANT_DIG-1).  */
  95               m.value = m.value - TWO_DBL_MANT_DIG;
  96             }
  97
  98           exponent =
  99             (m.word[DBL_EXPBIT0_WORD] >> DBL_EXPBIT0_BIT) & DBL_EXP_MASK;
 100           return exponent - DBL_EXP_BIAS;
 101         }
 102     }
 103   else
 104 # endif
 105     if (NBITS == 32)
 106       {
 107         /* 6 comparisons.  */
 108         if (x != 0)
 109           {
 110             int result = 1;
 111             if (x >= 0x10000)
 112               {
 113                 x = x >> 16;
 114                 result += 16;
 115               }
 116             if (x >= 0x100)
 117               {
 118                 x = x >> 8;
 119                 result += 8;
 120               }
 121             if (x >= 0x10)
 122               {
 123                 x = x >> 4;
 124                 result += 4;
 125               }
 126             if (x >= 0x4)
 127               {
 128                 x = x >> 2;
 129                 result += 2;
 130               }
 131             if (x >= 0x2)
 132               {
 133                 x = x >> 1;
 134                 result += 1;
 135               }
 136             return result;
 137           }
 138         else
 139           return 0;
 140       }
 141     else
 142       {
 143         /* Naive loop.
 144            Works for any value of NBITS.  */
 145         int j;
 146
 147         for (j = NBITS - 1; j >= 0; j--)
 148           if (x & (1U << j))
 149             return j + 1;
 150         return 0;
 151       }
 152 #endif
 153 }