math.hpp 7.9 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296
  1. /////////////////////////////////////////////////////////////////////////////
  2. //
  3. // (C) Copyright Ion Gaztanaga 2014-2014
  4. //
  5. // Distributed under the Boost Software License, Version 1.0.
  6. // (See accompanying file LICENSE_1_0.txt or copy at
  7. // http://www.boost.org/LICENSE_1_0.txt)
  8. //
  9. // See http://www.boost.org/libs/intrusive for documentation.
  10. //
  11. /////////////////////////////////////////////////////////////////////////////
  12. #ifndef BOOST_INTRUSIVE_DETAIL_MATH_HPP
  13. #define BOOST_INTRUSIVE_DETAIL_MATH_HPP
  14. #ifndef BOOST_CONFIG_HPP
  15. # include <boost/config.hpp>
  16. #endif
  17. #if defined(BOOST_HAS_PRAGMA_ONCE)
  18. # pragma once
  19. #endif
  20. #include <cstddef>
  21. #include <climits>
  22. #include <boost/intrusive/detail/mpl.hpp>
  23. namespace boost {
  24. namespace intrusive {
  25. namespace detail {
  26. ///////////////////////////
  27. // floor_log2 Dispatcher
  28. ////////////////////////////
  29. #if defined(_MSC_VER) && (_MSC_VER >= 1300)
  30. }}} //namespace boost::intrusive::detail
  31. //Use _BitScanReverseXX intrinsics
  32. #if defined(_M_X64) || defined(_M_AMD64) || defined(_M_IA64) //64 bit target
  33. #define BOOST_INTRUSIVE_BSR_INTRINSIC_64_BIT
  34. #endif
  35. #ifndef __INTRIN_H_ // Avoid including any windows system header
  36. #ifdef __cplusplus
  37. extern "C" {
  38. #endif // __cplusplus
  39. #if defined(BOOST_INTRUSIVE_BSR_INTRINSIC_64_BIT) //64 bit target
  40. unsigned char _BitScanReverse64(unsigned long *index, unsigned __int64 mask);
  41. #pragma intrinsic(_BitScanReverse64)
  42. #else //32 bit target
  43. unsigned char _BitScanReverse(unsigned long *index, unsigned long mask);
  44. #pragma intrinsic(_BitScanReverse)
  45. #endif
  46. #ifdef __cplusplus
  47. }
  48. #endif // __cplusplus
  49. #endif // __INTRIN_H_
  50. #ifdef BOOST_INTRUSIVE_BSR_INTRINSIC_64_BIT
  51. #define BOOST_INTRUSIVE_BSR_INTRINSIC _BitScanReverse64
  52. #undef BOOST_INTRUSIVE_BSR_INTRINSIC_64_BIT
  53. #else
  54. #define BOOST_INTRUSIVE_BSR_INTRINSIC _BitScanReverse
  55. #endif
  56. namespace boost {
  57. namespace intrusive {
  58. namespace detail {
  59. inline std::size_t floor_log2 (std::size_t x)
  60. {
  61. unsigned long log2;
  62. BOOST_INTRUSIVE_BSR_INTRINSIC( &log2, (unsigned long)x );
  63. return log2;
  64. }
  65. #undef BOOST_INTRUSIVE_BSR_INTRINSIC
  66. #elif defined(__GNUC__) && ((__GNUC__ >= 4) || (__GNUC__ == 3 && __GNUC_MINOR__ >= 4)) //GCC >=3.4
  67. //Compile-time error in case of missing specialization
  68. template<class Uint>
  69. struct builtin_clz_dispatch;
  70. #if defined(BOOST_HAS_LONG_LONG)
  71. template<>
  72. struct builtin_clz_dispatch< ::boost::ulong_long_type >
  73. {
  74. static ::boost::ulong_long_type call(::boost::ulong_long_type n)
  75. { return __builtin_clzll(n); }
  76. };
  77. #endif
  78. template<>
  79. struct builtin_clz_dispatch<unsigned long>
  80. {
  81. static unsigned long call(unsigned long n)
  82. { return __builtin_clzl(n); }
  83. };
  84. template<>
  85. struct builtin_clz_dispatch<unsigned int>
  86. {
  87. static unsigned int call(unsigned int n)
  88. { return __builtin_clz(n); }
  89. };
  90. inline std::size_t floor_log2(std::size_t n)
  91. {
  92. return sizeof(std::size_t)*CHAR_BIT - std::size_t(1) - builtin_clz_dispatch<std::size_t>::call(n);
  93. }
  94. #else //Portable methods
  95. ////////////////////////////
  96. // Generic method
  97. ////////////////////////////
  98. inline std::size_t floor_log2_get_shift(std::size_t n, true_ )//power of two size_t
  99. { return n >> 1; }
  100. inline std::size_t floor_log2_get_shift(std::size_t n, false_ )//non-power of two size_t
  101. { return (n >> 1) + ((n & 1u) & (n != 1)); }
  102. template<std::size_t N>
  103. inline std::size_t floor_log2 (std::size_t x, integral_constant<std::size_t, N>)
  104. {
  105. const std::size_t Bits = N;
  106. const bool Size_t_Bits_Power_2= !(Bits & (Bits-1));
  107. std::size_t n = x;
  108. std::size_t log2 = 0;
  109. std::size_t remaining_bits = Bits;
  110. std::size_t shift = floor_log2_get_shift(remaining_bits, bool_<Size_t_Bits_Power_2>());
  111. while(shift){
  112. std::size_t tmp = n >> shift;
  113. if (tmp){
  114. log2 += shift, n = tmp;
  115. }
  116. shift = floor_log2_get_shift(shift, bool_<Size_t_Bits_Power_2>());
  117. }
  118. return log2;
  119. }
  120. ////////////////////////////
  121. // DeBruijn method
  122. ////////////////////////////
  123. //Taken from:
  124. //http://stackoverflow.com/questions/11376288/fast-computing-of-log2-for-64-bit-integers
  125. //Thanks to Desmond Hume
  126. inline std::size_t floor_log2 (std::size_t v, integral_constant<std::size_t, 32>)
  127. {
  128. static const int MultiplyDeBruijnBitPosition[32] =
  129. {
  130. 0, 9, 1, 10, 13, 21, 2, 29, 11, 14, 16, 18, 22, 25, 3, 30,
  131. 8, 12, 20, 28, 15, 17, 24, 7, 19, 27, 23, 6, 26, 5, 4, 31
  132. };
  133. v |= v >> 1;
  134. v |= v >> 2;
  135. v |= v >> 4;
  136. v |= v >> 8;
  137. v |= v >> 16;
  138. return MultiplyDeBruijnBitPosition[(std::size_t)(v * 0x07C4ACDDU) >> 27];
  139. }
  140. inline std::size_t floor_log2 (std::size_t v, integral_constant<std::size_t, 64>)
  141. {
  142. static const std::size_t MultiplyDeBruijnBitPosition[64] = {
  143. 63, 0, 58, 1, 59, 47, 53, 2,
  144. 60, 39, 48, 27, 54, 33, 42, 3,
  145. 61, 51, 37, 40, 49, 18, 28, 20,
  146. 55, 30, 34, 11, 43, 14, 22, 4,
  147. 62, 57, 46, 52, 38, 26, 32, 41,
  148. 50, 36, 17, 19, 29, 10, 13, 21,
  149. 56, 45, 25, 31, 35, 16, 9, 12,
  150. 44, 24, 15, 8, 23, 7, 6, 5};
  151. v |= v >> 1;
  152. v |= v >> 2;
  153. v |= v >> 4;
  154. v |= v >> 8;
  155. v |= v >> 16;
  156. v |= v >> 32;
  157. return MultiplyDeBruijnBitPosition[((std::size_t)((v - (v >> 1))*0x07EDD5E59A4E28C2ULL)) >> 58];
  158. }
  159. inline std::size_t floor_log2 (std::size_t x)
  160. {
  161. const std::size_t Bits = sizeof(std::size_t)*CHAR_BIT;
  162. return floor_log2(x, integral_constant<std::size_t, Bits>());
  163. }
  164. #endif
  165. //Thanks to Laurent de Soras in
  166. //http://www.flipcode.com/archives/Fast_log_Function.shtml
  167. inline float fast_log2 (float val)
  168. {
  169. union caster_t
  170. {
  171. unsigned x;
  172. float val;
  173. } caster;
  174. caster.val = val;
  175. unsigned x = caster.x;
  176. const int log_2 = int((x >> 23) & 255) - 128;
  177. x &= ~(unsigned(255u) << 23u);
  178. x += unsigned(127) << 23u;
  179. caster.x = x;
  180. val = caster.val;
  181. //1+log2(m), m ranging from 1 to 2
  182. //3rd degree polynomial keeping first derivate continuity.
  183. //For less precision the line can be commented out
  184. val = ((-1.f/3.f) * val + 2.f) * val - (2.f/3.f);
  185. return val + static_cast<float>(log_2);
  186. }
  187. inline bool is_pow2(std::size_t x)
  188. { return (x & (x-1)) == 0; }
  189. template<std::size_t N>
  190. struct static_is_pow2
  191. {
  192. static const bool value = (N & (N-1)) == 0;
  193. };
  194. inline std::size_t ceil_log2 (std::size_t x)
  195. {
  196. return static_cast<std::size_t>(!(is_pow2)(x)) + floor_log2(x);
  197. }
  198. inline std::size_t ceil_pow2 (std::size_t x)
  199. {
  200. return std::size_t(1u) << (ceil_log2)(x);
  201. }
  202. inline std::size_t previous_or_equal_pow2(std::size_t x)
  203. {
  204. return std::size_t(1u) << floor_log2(x);
  205. }
  206. template<class SizeType, std::size_t N>
  207. struct numbits_eq
  208. {
  209. static const bool value = sizeof(SizeType)*CHAR_BIT == N;
  210. };
  211. template<class SizeType, class Enabler = void >
  212. struct sqrt2_pow_max;
  213. template <class SizeType>
  214. struct sqrt2_pow_max<SizeType, typename voider<typename enable_if< numbits_eq<SizeType, 32> >::type>::type>
  215. {
  216. static const SizeType value = 0xb504f334;
  217. static const std::size_t pow = 31;
  218. };
  219. #ifndef BOOST_NO_INT64_T
  220. template <class SizeType>
  221. struct sqrt2_pow_max<SizeType, typename voider<typename enable_if< numbits_eq<SizeType, 64> >::type>::type>
  222. {
  223. static const SizeType value = 0xb504f333f9de6484ull;
  224. static const std::size_t pow = 63;
  225. };
  226. #endif //BOOST_NO_INT64_T
  227. // Returns floor(pow(sqrt(2), x * 2 + 1)).
  228. // Defined for X from 0 up to the number of bits in size_t minus 1.
  229. inline std::size_t sqrt2_pow_2xplus1 (std::size_t x)
  230. {
  231. const std::size_t value = (std::size_t)sqrt2_pow_max<std::size_t>::value;
  232. const std::size_t pow = (std::size_t)sqrt2_pow_max<std::size_t>::pow;
  233. return (value >> (pow - x)) + 1;
  234. }
  235. } //namespace detail
  236. } //namespace intrusive
  237. } //namespace boost
  238. #endif //BOOST_INTRUSIVE_DETAIL_MATH_HPP