large_int.hpp 7.0 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255
  1. // Copyright 2018-2019 Hans Dembinski
  2. //
  3. // Distributed under the Boost Software License, Version 1.0.
  4. // (See accompanying file LICENSE_1_0.txt
  5. // or copy at http://www.boost.org/LICENSE_1_0.txt)
  6. #ifndef BOOST_HISTOGRAM_DETAIL_LARGE_INT_HPP
  7. #define BOOST_HISTOGRAM_DETAIL_LARGE_INT_HPP
  8. #include <boost/histogram/detail/operators.hpp>
  9. #include <boost/histogram/detail/safe_comparison.hpp>
  10. #include <boost/mp11/algorithm.hpp>
  11. #include <boost/mp11/function.hpp>
  12. #include <boost/mp11/list.hpp>
  13. #include <boost/mp11/utility.hpp>
  14. #include <cassert>
  15. #include <cmath>
  16. #include <cstdint>
  17. #include <limits>
  18. #include <type_traits>
  19. #include <utility>
  20. #include <vector>
  21. namespace boost {
  22. namespace histogram {
  23. namespace detail {
  24. template <class T>
  25. using is_unsigned_integral = mp11::mp_and<std::is_integral<T>, std::is_unsigned<T>>;
  26. template <class T>
  27. bool safe_increment(T& t) {
  28. if (t < (std::numeric_limits<T>::max)()) {
  29. ++t;
  30. return true;
  31. }
  32. return false;
  33. }
  34. template <class T, class U>
  35. bool safe_radd(T& t, const U& u) {
  36. static_assert(is_unsigned_integral<T>::value, "T must be unsigned integral type");
  37. static_assert(is_unsigned_integral<U>::value, "T must be unsigned integral type");
  38. if (static_cast<T>((std::numeric_limits<T>::max)() - t) >= u) {
  39. t += static_cast<T>(u); // static_cast to suppress conversion warning
  40. return true;
  41. }
  42. return false;
  43. }
  44. // An integer type which can grow arbitrarily large (until memory is exhausted).
  45. // Use boost.multiprecision.cpp_int in your own code, it is much more sophisticated.
  46. // We use it only to reduce coupling between boost libs.
  47. template <class Allocator>
  48. struct large_int : totally_ordered<large_int<Allocator>, large_int<Allocator>>,
  49. partially_ordered<large_int<Allocator>, void> {
  50. explicit large_int(const Allocator& a = {}) : data(1, 0, a) {}
  51. explicit large_int(std::uint64_t v, const Allocator& a = {}) : data(1, v, a) {}
  52. large_int(const large_int&) = default;
  53. large_int& operator=(const large_int&) = default;
  54. large_int(large_int&&) = default;
  55. large_int& operator=(large_int&&) = default;
  56. large_int& operator=(std::uint64_t o) {
  57. data = decltype(data)(1, o);
  58. return *this;
  59. }
  60. large_int& operator++() {
  61. assert(data.size() > 0u);
  62. std::size_t i = 0;
  63. while (!safe_increment(data[i])) {
  64. data[i] = 0;
  65. ++i;
  66. if (i == data.size()) {
  67. data.push_back(1);
  68. break;
  69. }
  70. }
  71. return *this;
  72. }
  73. large_int& operator+=(const large_int& o) {
  74. if (this == &o) {
  75. auto tmp{o};
  76. return operator+=(tmp);
  77. }
  78. bool carry = false;
  79. std::size_t i = 0;
  80. for (std::uint64_t oi : o.data) {
  81. auto& di = maybe_extend(i);
  82. if (carry) {
  83. if (safe_increment(oi))
  84. carry = false;
  85. else {
  86. ++i;
  87. continue;
  88. }
  89. }
  90. if (!safe_radd(di, oi)) {
  91. add_remainder(di, oi);
  92. carry = true;
  93. }
  94. ++i;
  95. }
  96. while (carry) {
  97. auto& di = maybe_extend(i);
  98. if (safe_increment(di)) break;
  99. di = 0;
  100. ++i;
  101. }
  102. return *this;
  103. }
  104. large_int& operator+=(std::uint64_t o) {
  105. assert(data.size() > 0u);
  106. if (safe_radd(data[0], o)) return *this;
  107. add_remainder(data[0], o);
  108. // carry the one, data may grow several times
  109. std::size_t i = 1;
  110. while (true) {
  111. auto& di = maybe_extend(i);
  112. if (safe_increment(di)) break;
  113. di = 0;
  114. ++i;
  115. }
  116. return *this;
  117. }
  118. explicit operator double() const noexcept {
  119. assert(data.size() > 0u);
  120. double result = static_cast<double>(data[0]);
  121. std::size_t i = 0;
  122. while (++i < data.size())
  123. result += static_cast<double>(data[i]) * std::pow(2.0, i * 64);
  124. return result;
  125. }
  126. bool operator<(const large_int& o) const noexcept {
  127. assert(data.size() > 0u);
  128. assert(o.data.size() > 0u);
  129. // no leading zeros allowed
  130. assert(data.size() == 1 || data.back() > 0u);
  131. assert(o.data.size() == 1 || o.data.back() > 0u);
  132. if (data.size() < o.data.size()) return true;
  133. if (data.size() > o.data.size()) return false;
  134. auto s = data.size();
  135. while (s > 0u) {
  136. --s;
  137. if (data[s] < o.data[s]) return true;
  138. if (data[s] > o.data[s]) return false;
  139. }
  140. return false; // args are equal
  141. }
  142. bool operator==(const large_int& o) const noexcept {
  143. assert(data.size() > 0u);
  144. assert(o.data.size() > 0u);
  145. // no leading zeros allowed
  146. assert(data.size() == 1 || data.back() > 0u);
  147. assert(o.data.size() == 1 || o.data.back() > 0u);
  148. if (data.size() != o.data.size()) return false;
  149. return std::equal(data.begin(), data.end(), o.data.begin());
  150. }
  151. template <class U>
  152. std::enable_if_t<std::is_integral<U>::value, bool> operator<(
  153. const U& o) const noexcept {
  154. assert(data.size() > 0u);
  155. return data.size() == 1 && safe_less()(data[0], o);
  156. }
  157. template <class U>
  158. std::enable_if_t<std::is_integral<U>::value, bool> operator>(
  159. const U& o) const noexcept {
  160. assert(data.size() > 0u);
  161. assert(data.size() == 1 || data.back() > 0u); // no leading zeros allowed
  162. return data.size() > 1 || safe_less()(o, data[0]);
  163. }
  164. template <class U>
  165. std::enable_if_t<std::is_integral<U>::value, bool> operator==(
  166. const U& o) const noexcept {
  167. assert(data.size() > 0u);
  168. return data.size() == 1 && safe_equal()(data[0], o);
  169. }
  170. template <class U>
  171. std::enable_if_t<std::is_floating_point<U>::value, bool> operator<(
  172. const U& o) const noexcept {
  173. return operator double() < o;
  174. }
  175. template <class U>
  176. std::enable_if_t<std::is_floating_point<U>::value, bool> operator>(
  177. const U& o) const noexcept {
  178. return operator double() > o;
  179. }
  180. template <class U>
  181. std::enable_if_t<std::is_floating_point<U>::value, bool> operator==(
  182. const U& o) const noexcept {
  183. return operator double() == o;
  184. }
  185. template <class U>
  186. std::enable_if_t<
  187. (!std::is_arithmetic<U>::value && std::is_convertible<U, double>::value), bool>
  188. operator<(const U& o) const noexcept {
  189. return operator double() < o;
  190. }
  191. template <class U>
  192. std::enable_if_t<
  193. (!std::is_arithmetic<U>::value && std::is_convertible<U, double>::value), bool>
  194. operator>(const U& o) const noexcept {
  195. return operator double() > o;
  196. }
  197. template <class U>
  198. std::enable_if_t<
  199. (!std::is_arithmetic<U>::value && std::is_convertible<U, double>::value), bool>
  200. operator==(const U& o) const noexcept {
  201. return operator double() == o;
  202. }
  203. std::uint64_t& maybe_extend(std::size_t i) {
  204. while (i >= data.size()) data.push_back(0);
  205. return data[i];
  206. }
  207. static void add_remainder(std::uint64_t& d, const std::uint64_t o) noexcept {
  208. assert(d > 0u);
  209. // in decimal system it would look like this:
  210. // 8 + 8 = 6 = 8 - (9 - 8) - 1
  211. // 9 + 1 = 0 = 9 - (9 - 1) - 1
  212. auto tmp = (std::numeric_limits<std::uint64_t>::max)();
  213. tmp -= o;
  214. --d -= tmp;
  215. }
  216. template <class Archive>
  217. void serialize(Archive& ar, unsigned /* version */) {
  218. ar& make_nvp("data", data);
  219. }
  220. std::vector<std::uint64_t, Allocator> data;
  221. };
  222. } // namespace detail
  223. } // namespace histogram
  224. } // namespace boost
  225. #endif