372 lines · c
1//===----------------------------------------------------------------------===//2//3// Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions.4// See https://llvm.org/LICENSE.txt for license information.5// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception6//7//===----------------------------------------------------------------------===//8 9#ifndef _LIBCPP___RANDOM_LINEAR_CONGRUENTIAL_ENGINE_H10#define _LIBCPP___RANDOM_LINEAR_CONGRUENTIAL_ENGINE_H11 12#include <__config>13#include <__random/is_seed_sequence.h>14#include <__type_traits/enable_if.h>15#include <__type_traits/integral_constant.h>16#include <__type_traits/is_unsigned.h>17#include <cstdint>18#include <iosfwd>19 20#if !defined(_LIBCPP_HAS_NO_PRAGMA_SYSTEM_HEADER)21# pragma GCC system_header22#endif23 24_LIBCPP_PUSH_MACROS25#include <__undef_macros>26 27_LIBCPP_BEGIN_NAMESPACE_STD28 29enum __lce_alg_type {30 _LCE_Full,31 _LCE_Part,32 _LCE_Schrage,33 _LCE_Promote,34};35 36template <unsigned long long __a,37 unsigned long long __c,38 unsigned long long __m,39 unsigned long long _Mp,40 bool _HasOverflow = (__a != 0ull && (__m & (__m - 1ull)) != 0ull), // a != 0, m != 0, m != 2^n41 bool _Full = (!_HasOverflow || __m - 1ull <= (_Mp - __c) / __a), // (a * x + c) % m works42 bool _Part = (!_HasOverflow || __m - 1ull <= _Mp / __a), // (a * x) % m works43 bool _Schrage = (_HasOverflow && __m % __a <= __m / __a)> // r <= q44struct __lce_alg_picker {45 static _LIBCPP_CONSTEXPR const __lce_alg_type __mode =46 _Full ? _LCE_Full47 : _Part ? _LCE_Part48 : _Schrage ? _LCE_Schrage49 : _LCE_Promote;50 51#if !_LIBCPP_HAS_INT12852 static_assert(_Mp != (unsigned long long)(-1) || _Full || _Part || _Schrage,53 "The current values for a, c, and m are not currently supported on platforms without __int128");54#endif55};56 57template <unsigned long long __a,58 unsigned long long __c,59 unsigned long long __m,60 unsigned long long _Mp,61 __lce_alg_type _Mode = __lce_alg_picker<__a, __c, __m, _Mp>::__mode>62struct __lce_ta;63 64// 6465 66#if _LIBCPP_HAS_INT12867template <unsigned long long _Ap, unsigned long long _Cp, unsigned long long _Mp>68struct __lce_ta<_Ap, _Cp, _Mp, (unsigned long long)(-1), _LCE_Promote> {69 typedef unsigned long long result_type;70 _LIBCPP_HIDE_FROM_ABI static result_type next(result_type __xp) {71 __extension__ using __calc_type = unsigned __int128;72 const __calc_type __a = static_cast<__calc_type>(_Ap);73 const __calc_type __c = static_cast<__calc_type>(_Cp);74 const __calc_type __m = static_cast<__calc_type>(_Mp);75 const __calc_type __x = static_cast<__calc_type>(__xp);76 return static_cast<result_type>((__a * __x + __c) % __m);77 }78};79#endif80 81template <unsigned long long __a, unsigned long long __c, unsigned long long __m>82struct __lce_ta<__a, __c, __m, (unsigned long long)(-1), _LCE_Schrage> {83 typedef unsigned long long result_type;84 _LIBCPP_HIDE_FROM_ABI static result_type next(result_type __x) {85 // Schrage's algorithm86 const result_type __q = __m / __a;87 const result_type __r = __m % __a;88 const result_type __t0 = __a * (__x % __q);89 const result_type __t1 = __r * (__x / __q);90 __x = __t0 + (__t0 < __t1) * __m - __t1;91 __x += __c - (__x >= __m - __c) * __m;92 return __x;93 }94};95 96template <unsigned long long __a, unsigned long long __m>97struct __lce_ta<__a, 0ull, __m, (unsigned long long)(-1), _LCE_Schrage> {98 typedef unsigned long long result_type;99 _LIBCPP_HIDE_FROM_ABI static result_type next(result_type __x) {100 // Schrage's algorithm101 const result_type __q = __m / __a;102 const result_type __r = __m % __a;103 const result_type __t0 = __a * (__x % __q);104 const result_type __t1 = __r * (__x / __q);105 __x = __t0 + (__t0 < __t1) * __m - __t1;106 return __x;107 }108};109 110template <unsigned long long __a, unsigned long long __c, unsigned long long __m>111struct __lce_ta<__a, __c, __m, (unsigned long long)(-1), _LCE_Part> {112 typedef unsigned long long result_type;113 _LIBCPP_HIDE_FROM_ABI static result_type next(result_type __x) {114 // Use (((a*x) % m) + c) % m115 __x = (__a * __x) % __m;116 __x += __c - (__x >= __m - __c) * __m;117 return __x;118 }119};120 121template <unsigned long long __a, unsigned long long __c, unsigned long long __m>122struct __lce_ta<__a, __c, __m, (unsigned long long)(-1), _LCE_Full> {123 typedef unsigned long long result_type;124 _LIBCPP_HIDE_FROM_ABI static result_type next(result_type __x) { return (__a * __x + __c) % __m; }125};126 127template <unsigned long long __a, unsigned long long __c>128struct __lce_ta<__a, __c, 0ull, (unsigned long long)(-1), _LCE_Full> {129 typedef unsigned long long result_type;130 _LIBCPP_HIDE_FROM_ABI static result_type next(result_type __x) { return __a * __x + __c; }131};132 133// 32134 135template <unsigned long long __a, unsigned long long __c, unsigned long long __m>136struct __lce_ta<__a, __c, __m, unsigned(-1), _LCE_Promote> {137 typedef unsigned result_type;138 _LIBCPP_HIDE_FROM_ABI static result_type next(result_type __x) {139 return static_cast<result_type>(__lce_ta<__a, __c, __m, (unsigned long long)(-1)>::next(__x));140 }141};142 143template <unsigned long long _Ap, unsigned long long _Cp, unsigned long long _Mp>144struct __lce_ta<_Ap, _Cp, _Mp, unsigned(-1), _LCE_Schrage> {145 typedef unsigned result_type;146 _LIBCPP_HIDE_FROM_ABI static result_type next(result_type __x) {147 const result_type __a = static_cast<result_type>(_Ap);148 const result_type __c = static_cast<result_type>(_Cp);149 const result_type __m = static_cast<result_type>(_Mp);150 // Schrage's algorithm151 const result_type __q = __m / __a;152 const result_type __r = __m % __a;153 const result_type __t0 = __a * (__x % __q);154 const result_type __t1 = __r * (__x / __q);155 __x = __t0 + (__t0 < __t1) * __m - __t1;156 __x += __c - (__x >= __m - __c) * __m;157 return __x;158 }159};160 161template <unsigned long long _Ap, unsigned long long _Mp>162struct __lce_ta<_Ap, 0ull, _Mp, unsigned(-1), _LCE_Schrage> {163 typedef unsigned result_type;164 _LIBCPP_HIDE_FROM_ABI static result_type next(result_type __x) {165 const result_type __a = static_cast<result_type>(_Ap);166 const result_type __m = static_cast<result_type>(_Mp);167 // Schrage's algorithm168 const result_type __q = __m / __a;169 const result_type __r = __m % __a;170 const result_type __t0 = __a * (__x % __q);171 const result_type __t1 = __r * (__x / __q);172 __x = __t0 + (__t0 < __t1) * __m - __t1;173 return __x;174 }175};176 177template <unsigned long long _Ap, unsigned long long _Cp, unsigned long long _Mp>178struct __lce_ta<_Ap, _Cp, _Mp, unsigned(-1), _LCE_Part> {179 typedef unsigned result_type;180 _LIBCPP_HIDE_FROM_ABI static result_type next(result_type __x) {181 const result_type __a = static_cast<result_type>(_Ap);182 const result_type __c = static_cast<result_type>(_Cp);183 const result_type __m = static_cast<result_type>(_Mp);184 // Use (((a*x) % m) + c) % m185 __x = (__a * __x) % __m;186 __x += __c - (__x >= __m - __c) * __m;187 return __x;188 }189};190 191template <unsigned long long _Ap, unsigned long long _Cp, unsigned long long _Mp>192struct __lce_ta<_Ap, _Cp, _Mp, unsigned(-1), _LCE_Full> {193 typedef unsigned result_type;194 _LIBCPP_HIDE_FROM_ABI static result_type next(result_type __x) {195 const result_type __a = static_cast<result_type>(_Ap);196 const result_type __c = static_cast<result_type>(_Cp);197 const result_type __m = static_cast<result_type>(_Mp);198 return (__a * __x + __c) % __m;199 }200};201 202template <unsigned long long _Ap, unsigned long long _Cp>203struct __lce_ta<_Ap, _Cp, 0ull, unsigned(-1), _LCE_Full> {204 typedef unsigned result_type;205 _LIBCPP_HIDE_FROM_ABI static result_type next(result_type __x) {206 const result_type __a = static_cast<result_type>(_Ap);207 const result_type __c = static_cast<result_type>(_Cp);208 return __a * __x + __c;209 }210};211 212// 16213 214template <unsigned long long __a, unsigned long long __c, unsigned long long __m, __lce_alg_type __mode>215struct __lce_ta<__a, __c, __m, (unsigned short)(-1), __mode> {216 typedef unsigned short result_type;217 _LIBCPP_HIDE_FROM_ABI static result_type next(result_type __x) {218 return static_cast<result_type>(__lce_ta<__a, __c, __m, unsigned(-1)>::next(__x));219 }220};221 222template <class _UIntType, _UIntType __a, _UIntType __c, _UIntType __m>223class linear_congruential_engine;224 225template <class _CharT, class _Traits, class _Up, _Up _Ap, _Up _Cp, _Up _Np>226_LIBCPP_HIDE_FROM_ABI basic_ostream<_CharT, _Traits>&227operator<<(basic_ostream<_CharT, _Traits>& __os, const linear_congruential_engine<_Up, _Ap, _Cp, _Np>&);228 229template <class _CharT, class _Traits, class _Up, _Up _Ap, _Up _Cp, _Up _Np>230_LIBCPP_HIDE_FROM_ABI basic_istream<_CharT, _Traits>&231operator>>(basic_istream<_CharT, _Traits>& __is, linear_congruential_engine<_Up, _Ap, _Cp, _Np>& __x);232 233template <class _UIntType, _UIntType __a, _UIntType __c, _UIntType __m>234class linear_congruential_engine {235public:236 // types237 typedef _UIntType result_type;238 239private:240 result_type __x_;241 242 static _LIBCPP_CONSTEXPR const result_type _Mp = result_type(-1);243 244 static_assert(__m == 0 || __a < __m, "linear_congruential_engine invalid parameters");245 static_assert(__m == 0 || __c < __m, "linear_congruential_engine invalid parameters");246 static_assert(is_unsigned<_UIntType>::value, "_UIntType must be unsigned type");247 248public:249 static _LIBCPP_CONSTEXPR const result_type _Min = __c == 0u ? 1u : 0u;250 static _LIBCPP_CONSTEXPR const result_type _Max = __m - _UIntType(1u);251 static_assert(_Min < _Max, "linear_congruential_engine invalid parameters");252 253 // engine characteristics254 static inline _LIBCPP_CONSTEXPR const result_type multiplier = __a;255 static inline _LIBCPP_CONSTEXPR const result_type increment = __c;256 static inline _LIBCPP_CONSTEXPR const result_type modulus = __m;257 _LIBCPP_HIDE_FROM_ABI static _LIBCPP_CONSTEXPR result_type min() { return _Min; }258 _LIBCPP_HIDE_FROM_ABI static _LIBCPP_CONSTEXPR result_type max() { return _Max; }259 static inline _LIBCPP_CONSTEXPR const result_type default_seed = 1u;260 261 // constructors and seeding functions262#ifndef _LIBCPP_CXX03_LANG263 _LIBCPP_HIDE_FROM_ABI linear_congruential_engine() : linear_congruential_engine(default_seed) {}264 _LIBCPP_HIDE_FROM_ABI explicit linear_congruential_engine(result_type __s) { seed(__s); }265#else266 _LIBCPP_HIDE_FROM_ABI explicit linear_congruential_engine(result_type __s = default_seed) { seed(__s); }267#endif268 template <class _Sseq, __enable_if_t<__is_seed_sequence<_Sseq, linear_congruential_engine>::value, int> = 0>269 _LIBCPP_HIDE_FROM_ABI explicit linear_congruential_engine(_Sseq& __q) {270 seed(__q);271 }272 _LIBCPP_HIDE_FROM_ABI void seed(result_type __s = default_seed) {273 seed(integral_constant<bool, __m == 0>(), integral_constant<bool, __c == 0>(), __s);274 }275 template <class _Sseq, __enable_if_t<__is_seed_sequence<_Sseq, linear_congruential_engine>::value, int> = 0>276 _LIBCPP_HIDE_FROM_ABI void seed(_Sseq& __q) {277 __seed(278 __q,279 integral_constant<unsigned,280 1 + (__m == 0 ? (sizeof(result_type) * __CHAR_BIT__ - 1) / 32 : (__m > 0x100000000ull))>());281 }282 283 // generating functions284 _LIBCPP_HIDE_FROM_ABI result_type operator()() {285 return __x_ = static_cast<result_type>(__lce_ta<__a, __c, __m, _Mp>::next(__x_));286 }287 _LIBCPP_HIDE_FROM_ABI void discard(unsigned long long __z) {288 for (; __z; --__z)289 operator()();290 }291 292 friend _LIBCPP_HIDE_FROM_ABI bool293 operator==(const linear_congruential_engine& __x, const linear_congruential_engine& __y) {294 return __x.__x_ == __y.__x_;295 }296 friend _LIBCPP_HIDE_FROM_ABI bool297 operator!=(const linear_congruential_engine& __x, const linear_congruential_engine& __y) {298 return !(__x == __y);299 }300 301private:302 _LIBCPP_HIDE_FROM_ABI void seed(true_type, true_type, result_type __s) { __x_ = __s == 0 ? 1 : __s; }303 _LIBCPP_HIDE_FROM_ABI void seed(true_type, false_type, result_type __s) { __x_ = __s; }304 _LIBCPP_HIDE_FROM_ABI void seed(false_type, true_type, result_type __s) { __x_ = __s % __m == 0 ? 1 : __s % __m; }305 _LIBCPP_HIDE_FROM_ABI void seed(false_type, false_type, result_type __s) { __x_ = __s % __m; }306 307 template <class _Sseq>308 _LIBCPP_HIDE_FROM_ABI void __seed(_Sseq& __q, integral_constant<unsigned, 1>);309 template <class _Sseq>310 _LIBCPP_HIDE_FROM_ABI void __seed(_Sseq& __q, integral_constant<unsigned, 2>);311 312 template <class _CharT, class _Traits, class _Up, _Up _Ap, _Up _Cp, _Up _Np>313 friend basic_ostream<_CharT, _Traits>&314 operator<<(basic_ostream<_CharT, _Traits>& __os, const linear_congruential_engine<_Up, _Ap, _Cp, _Np>&);315 316 template <class _CharT, class _Traits, class _Up, _Up _Ap, _Up _Cp, _Up _Np>317 friend basic_istream<_CharT, _Traits>&318 operator>>(basic_istream<_CharT, _Traits>& __is, linear_congruential_engine<_Up, _Ap, _Cp, _Np>& __x);319};320 321template <class _UIntType, _UIntType __a, _UIntType __c, _UIntType __m>322template <class _Sseq>323void linear_congruential_engine<_UIntType, __a, __c, __m>::__seed(_Sseq& __q, integral_constant<unsigned, 1>) {324 const unsigned __k = 1;325 uint32_t __ar[__k + 3];326 __q.generate(__ar, __ar + __k + 3);327 result_type __s = static_cast<result_type>(__ar[3] % __m);328 __x_ = __c == 0 && __s == 0 ? result_type(1) : __s;329}330 331template <class _UIntType, _UIntType __a, _UIntType __c, _UIntType __m>332template <class _Sseq>333void linear_congruential_engine<_UIntType, __a, __c, __m>::__seed(_Sseq& __q, integral_constant<unsigned, 2>) {334 const unsigned __k = 2;335 uint32_t __ar[__k + 3];336 __q.generate(__ar, __ar + __k + 3);337 result_type __s = static_cast<result_type>((__ar[3] + ((uint64_t)__ar[4] << 32)) % __m);338 __x_ = __c == 0 && __s == 0 ? result_type(1) : __s;339}340 341template <class _CharT, class _Traits, class _UIntType, _UIntType __a, _UIntType __c, _UIntType __m>342inline _LIBCPP_HIDE_FROM_ABI basic_ostream<_CharT, _Traits>&343operator<<(basic_ostream<_CharT, _Traits>& __os, const linear_congruential_engine<_UIntType, __a, __c, __m>& __x) {344 __save_flags<_CharT, _Traits> __lx(__os);345 typedef basic_ostream<_CharT, _Traits> _Ostream;346 __os.flags(_Ostream::dec | _Ostream::left);347 __os.fill(__os.widen(' '));348 return __os << __x.__x_;349}350 351template <class _CharT, class _Traits, class _UIntType, _UIntType __a, _UIntType __c, _UIntType __m>352_LIBCPP_HIDE_FROM_ABI basic_istream<_CharT, _Traits>&353operator>>(basic_istream<_CharT, _Traits>& __is, linear_congruential_engine<_UIntType, __a, __c, __m>& __x) {354 __save_flags<_CharT, _Traits> __lx(__is);355 typedef basic_istream<_CharT, _Traits> _Istream;356 __is.flags(_Istream::dec | _Istream::skipws);357 _UIntType __t;358 __is >> __t;359 if (!__is.fail())360 __x.__x_ = __t;361 return __is;362}363 364typedef linear_congruential_engine<uint_fast32_t, 16807, 0, 2147483647> minstd_rand0;365typedef linear_congruential_engine<uint_fast32_t, 48271, 0, 2147483647> minstd_rand;366 367_LIBCPP_END_NAMESPACE_STD368 369_LIBCPP_POP_MACROS370 371#endif // _LIBCPP___RANDOM_LINEAR_CONGRUENTIAL_ENGINE_H372