libstdc++
|
00001 // -*- C++ -*- 00002 00003 // Copyright (C) 2007, 2008, 2009 Free Software Foundation, Inc. 00004 // 00005 // This file is part of the GNU ISO C++ Library. This library is free 00006 // software; you can redistribute it and/or modify it under the terms 00007 // of the GNU General Public License as published by the Free Software 00008 // Foundation; either version 3, or (at your option) any later 00009 // version. 00010 00011 // This library is distributed in the hope that it will be useful, but 00012 // WITHOUT ANY WARRANTY; without even the implied warranty of 00013 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU 00014 // General Public License for more details. 00015 00016 // Under Section 7 of GPL version 3, you are granted additional 00017 // permissions described in the GCC Runtime Library Exception, version 00018 // 3.1, as published by the Free Software Foundation. 00019 00020 // You should have received a copy of the GNU General Public License and 00021 // a copy of the GCC Runtime Library Exception along with this program; 00022 // see the files COPYING3 and COPYING.RUNTIME respectively. If not, see 00023 // <http://www.gnu.org/licenses/>. 00024 00025 /** @file parallel/random_number.h 00026 * @brief Random number generator based on the Mersenne twister. 00027 * This file is a GNU parallel extension to the Standard C++ Library. 00028 */ 00029 00030 // Written by Johannes Singler. 00031 00032 #ifndef _GLIBCXX_PARALLEL_RANDOM_NUMBER_H 00033 #define _GLIBCXX_PARALLEL_RANDOM_NUMBER_H 1 00034 00035 #include <parallel/types.h> 00036 #include <tr1/random> 00037 00038 namespace __gnu_parallel 00039 { 00040 /** @brief Random number generator, based on the Mersenne twister. */ 00041 class random_number 00042 { 00043 private: 00044 std::tr1::mt19937 mt; 00045 uint64 supremum; 00046 uint64 RAND_SUP; 00047 double supremum_reciprocal; 00048 double RAND_SUP_REC; 00049 00050 // Assumed to be twice as long as the usual random number. 00051 uint64 cache; 00052 00053 // Bit results. 00054 int bits_left; 00055 00056 static uint32 00057 scale_down(uint64 x, 00058 #if _GLIBCXX_SCALE_DOWN_FPU 00059 uint64 /*supremum*/, double supremum_reciprocal) 00060 #else 00061 uint64 supremum, double /*supremum_reciprocal*/) 00062 #endif 00063 { 00064 #if _GLIBCXX_SCALE_DOWN_FPU 00065 return uint32(x * supremum_reciprocal); 00066 #else 00067 return static_cast<uint32>(x % supremum); 00068 #endif 00069 } 00070 00071 public: 00072 /** @brief Default constructor. Seed with 0. */ 00073 random_number() 00074 : mt(0), supremum(0x100000000ULL), 00075 RAND_SUP(1ULL << (sizeof(uint32) * 8)), 00076 supremum_reciprocal(double(supremum) / double(RAND_SUP)), 00077 RAND_SUP_REC(1.0 / double(RAND_SUP)), 00078 cache(0), bits_left(0) { } 00079 00080 /** @brief Constructor. 00081 * @param seed Random seed. 00082 * @param supremum Generate integer random numbers in the 00083 * interval @c [0,supremum). */ 00084 random_number(uint32 seed, uint64 supremum = 0x100000000ULL) 00085 : mt(seed), supremum(supremum), 00086 RAND_SUP(1ULL << (sizeof(uint32) * 8)), 00087 supremum_reciprocal(double(supremum) / double(RAND_SUP)), 00088 RAND_SUP_REC(1.0 / double(RAND_SUP)), 00089 cache(0), bits_left(0) { } 00090 00091 /** @brief Generate unsigned random 32-bit integer. */ 00092 uint32 00093 operator()() 00094 { return scale_down(mt(), supremum, supremum_reciprocal); } 00095 00096 /** @brief Generate unsigned random 32-bit integer in the 00097 interval @c [0,local_supremum). */ 00098 uint32 00099 operator()(uint64 local_supremum) 00100 { 00101 return scale_down(mt(), local_supremum, 00102 double(local_supremum * RAND_SUP_REC)); 00103 } 00104 00105 /** @brief Generate a number of random bits, run-time parameter. 00106 * @param bits Number of bits to generate. */ 00107 unsigned long 00108 genrand_bits(int bits) 00109 { 00110 unsigned long res = cache & ((1 << bits) - 1); 00111 cache = cache >> bits; 00112 bits_left -= bits; 00113 if (bits_left < 32) 00114 { 00115 cache |= ((uint64(mt())) << bits_left); 00116 bits_left += 32; 00117 } 00118 return res; 00119 } 00120 }; 00121 00122 } // namespace __gnu_parallel 00123 00124 #endif /* _GLIBCXX_PARALLEL_RANDOM_NUMBER_H */