random/uniform.h

Go to the documentation of this file.
00001 #ifndef BZ_RANDOM_UNIFORM_H
00002 #define BZ_RANDOM_UNIFORM_H
00003 
00004 #include <random/default.h>
00005 
00006 #ifndef FLT_MANT_DIG
00007  #include <float.h>
00008 #endif
00009 
00010 BZ_NAMESPACE(ranlib)
00011 
00012 /*****************************************************************************
00013  * UniformClosedOpen generator: uniform random numbers in [0,1).
00014  *****************************************************************************/
00015 
00016 template<typename T = defaultType, typename IRNG = defaultIRNG, 
00017     typename stateTag = defaultState>
00018 class UniformClosedOpen { };
00019 
00020 // These constants are 1/2^32, 1/2^64, 1/2^96, 1/2^128
00021 const long double norm32open = .2328306436538696289062500000000000000000E-9L;
00022 const long double norm64open = .5421010862427522170037264004349708557129E-19L;
00023 const long double norm96open = .1262177448353618888658765704452457967477E-28L;
00024 const long double norm128open = .2938735877055718769921841343055614194547E-38L;
00025 
00026 
00027 template<typename IRNG, typename stateTag>
00028 class UniformClosedOpen<float,IRNG,stateTag> 
00029   : public IRNGWrapper<IRNG,stateTag> 
00030 {
00031 
00032 public:
00033     typedef float T_numtype;
00034 
00035     float random()
00036     {
00037 #if FLT_MANT_DIG > 96
00038  #if (FLT_MANT_DIG > 128) && !defined(BZ_IGNORE_RNG_ERRORS)
00039   #error RNG code assumes float mantissa is <= 128 bits (not true for your platform).  Use -DBZ_IGNORE_RNG_ERRORS to ignore this warning.
00040  #endif
00041         IRNG_int i1 = this->irng_.random();
00042         IRNG_int i2 = this->irng_.random();
00043         IRNG_int i3 = this->irng_.random();
00044         IRNG_int i4 = this->irng_.random();
00045         return i1 * norm128open + i2 * norm96open + i3 * norm64open
00046             + i4 * norm32open;
00047 #elif FLT_MANT_DIG > 64
00048         IRNG_int i1 = this->irng_.random();
00049         IRNG_int i2 = this->irng_.random();
00050         IRNG_int i3 = this->irng_.random();
00051         return i1 * norm96open + i2 * norm64open + i3 * norm32open;
00052 #elif FLT_MANT_DIG > 32
00053         IRNG_int i1 = this->irng_.random();
00054         IRNG_int i2 = this->irng_.random();
00055         return i1 * norm64open + i2 * norm32open;
00056 #else
00057         IRNG_int i1 = this->irng_.random();
00058         return i1 * norm32open;
00059 #endif
00060     }    
00061 
00062     float getUniform()
00063     { return random(); }
00064 };
00065 
00066 template<typename IRNG, typename stateTag>
00067 class UniformClosedOpen<double,IRNG,stateTag> 
00068   : public IRNGWrapper<IRNG,stateTag> 
00069 {
00070 public:
00071     typedef double T_numtype;
00072 
00073     double random()
00074     {    
00075 #if DBL_MANT_DIG > 96
00076  #if (DBL_MANT_DIG > 128) && !defined(BZ_IGNORE_RNG_ERRORS)
00077   #error RNG code assumes double mantissa is <= 128 bits (not true for your platform).  Use -DBZ_IGNORE_RNG_ERRORS to ignore this warning.
00078  #endif
00079 
00080         IRNG_int i1 = this->irng_.random();
00081         IRNG_int i2 = this->irng_.random();
00082         IRNG_int i3 = this->irng_.random();
00083         IRNG_int i4 = this->irng_.random();
00084         return i1 * norm128open + i2 * norm96open + i3 * norm64open
00085             + i4 * norm32open;
00086 #elif DBL_MANT_DIG > 64
00087         IRNG_int i1 = this->irng_.random();
00088         IRNG_int i2 = this->irng_.random();
00089         IRNG_int i3 = this->irng_.random();
00090         return i1 * norm96open + i2 * norm64open + i3 * norm32open;
00091 #elif DBL_MANT_DIG > 32
00092         IRNG_int i1 = this->irng_.random();
00093         IRNG_int i2 = this->irng_.random();
00094         return i1 * norm64open + i2 * norm32open;
00095 #else
00096         IRNG_int i1 = this->irng_.random();
00097         return i1 * norm32open;
00098 #endif
00099 
00100     }
00101 
00102     double getUniform() { return random(); }
00103 };
00104 
00105 template<typename IRNG, typename stateTag>
00106 class UniformClosedOpen<long double,IRNG,stateTag>
00107   : public IRNGWrapper<IRNG,stateTag> 
00108 {
00109 public:
00110     typedef long double T_numtype;
00111 
00112     long double random()
00113     {
00114 #if LDBL_MANT_DIG > 96
00115  #if (LDBL_MANT_DIG > 128) && !defined(BZ_IGNORE_RNG_ERRORS)
00116   #error RNG code assumes long double mantissa is <= 128 bits (not true for your platform).  Use -DBZ_IGNORE_RNG_ERRORS to ignore this warning.  
00117 #endif
00118 
00119         IRNG_int i1 = this->irng_.random();
00120         IRNG_int i2 = this->irng_.random();
00121         IRNG_int i3 = this->irng_.random();
00122         IRNG_int i4 = this->irng_.random();
00123         return i1 * norm128open + i2 * norm96open + i3 * norm64open
00124             + i4 * norm32open;
00125 #elif LDBL_MANT_DIG > 64
00126         IRNG_int i1 = this->irng_.random();
00127         IRNG_int i2 = this->irng_.random();
00128         IRNG_int i3 = this->irng_.random();
00129         return i1 * norm96open + i2 * norm64open + i3 * norm32open;
00130 #elif LDBL_MANT_DIG > 32
00131         IRNG_int i1 = this->irng_.random();
00132         IRNG_int i2 = this->irng_.random();
00133         return i1 * norm64open + i2 * norm32open;
00134 #else
00135         IRNG_int i1 = this->irng_.random();
00136         return i1 * norm32open;
00137 #endif
00138     }
00139 
00140     long double getUniform() { return random(); }
00141 };
00142 
00143 // For people who don't care about open or closed: just give them
00144 // ClosedOpen (this is the fastest).
00145 
00146 template<class T, typename IRNG = defaultIRNG, 
00147     typename stateTag = defaultState>
00148 class Uniform : public UniformClosedOpen<T,IRNG,stateTag> 
00149 { };
00150 
00151 /*****************************************************************************
00152  * UniformClosed generator: uniform random numbers in [0,1].
00153  *****************************************************************************/
00154 
00155 // This constant is 1/(2^32-1)
00156 const long double norm32closed = .2328306437080797375431469961868475648078E-9L;
00157 
00158 // These constants are 2^32/(2^64-1) and 1/(2^64-1), respectively.
00159 
00160 const long double norm64closed1 =
00161     .23283064365386962891887177448353618888727188481031E-9L;
00162 const long double norm64closed2 =
00163     .54210108624275221703311375920552804341370213034169E-19L;
00164 
00165 // These constants are 2^64/(2^96-1), 2^32/(2^96-1), and 1/(2^96-1)
00166 const long double norm96closed1 = .2328306436538696289062500000029387358771E-9L;
00167 const long double norm96closed2 =
00168     .5421010862427522170037264004418131333707E-19L;
00169 const long double norm96closed3 =
00170     .1262177448353618888658765704468388886588E-28L;
00171 
00172 // These constants are 2^96/(2^128-1), 2^64/(2^128-1), 2^32/(2^128-1) and
00173 // 1/(2^128-1).
00174 const long double norm128clos1 = .2328306436538696289062500000000000000007E-9L;
00175 const long double norm128clos2 = .5421010862427522170037264004349708557145E-19L;
00176 const long double norm128clos3 = .1262177448353618888658765704452457967481E-28L;
00177 const long double norm128clos4 = .2938735877055718769921841343055614194555E-38L;
00178 
00179 
00180 template<typename T = double, typename IRNG = defaultIRNG, 
00181     typename stateTag = defaultState>
00182 class UniformClosed { };
00183 
00184 template<typename IRNG, typename stateTag>
00185 class UniformClosed<float,IRNG,stateTag> : public IRNGWrapper<IRNG,stateTag> {
00186 
00187 public:
00188     typedef float T_numtype;
00189 
00190     float random()
00191     {
00192 #if FLTMANT_DIG > 96
00193  #if (FLT_MANT_DIG > 128) && !defined(BZ_IGNORE_RNG_ERRORS)
00194   #error RNG code assumes float mantissa is <= 128 bits (not true for your platform).  Use -DBZ_IGNORE_RNG_ERRORS to ignore this warning.
00195  #endif
00196         IRNG_int i1 = this->irng_.random();
00197         IRNG_int i2 = this->irng_.random();
00198         IRNG_int i3 = this->irng_.random();
00199         IRNG_int i4 = this->irng_.random();
00200 
00201         return i1 * norm128clos1 + i2 * norm128clos2
00202           + i3 * norm128clos3 + i4 * norm128clos4;
00203 #elif FLT_MANT_DIG > 64
00204         IRNG_int i1 = this->irng_.random();
00205         IRNG_int i2 = this->irng_.random();
00206         IRNG_int i3 = this->irng_.random();
00207 
00208         return i1 * norm96closed1 + i2 * norm96closed2
00209           + i3 * norm96closed3;
00210 #elif FLT_MANT_DIG > 32
00211         IRNG_int i1 = this->irng_.random();
00212         IRNG_int i2 = this->irng_.random();
00213 
00214         return i1 * norm64closed1 + i2 * norm64closed2;
00215 #else
00216         IRNG_int i = this->irng_.random();
00217         return i * norm32closed;
00218 #endif
00219 
00220     }
00221 
00222     float getUniform()
00223     { return random(); }
00224 };
00225 
00226 template<typename IRNG, typename stateTag>
00227 class UniformClosed<double,IRNG,stateTag> : public IRNGWrapper<IRNG,stateTag> {
00228 
00229 public:
00230     typedef double T_numtype;
00231 
00232     double random()
00233     {
00234 #if DBL_MANT_DIG > 96
00235  #if (DBL_MANT_DIG > 128) && !defined(BZ_IGNORE_RNG_ERRORS)
00236   #error RNG code assumes double mantissa is <= 128 bits (not true for your platform).  Use -DBZ_IGNORE_RNG_ERRORS to ignore this warning.
00237  #endif
00238         IRNG_int i1 = this->irng_.random();
00239         IRNG_int i2 = this->irng_.random();
00240         IRNG_int i3 = this->irng_.random();
00241         IRNG_int i4 = this->irng_.random();
00242 
00243         return i1 * norm128clos1 + i2 * norm128clos2
00244           + i3 * norm128clos3 + i4 * norm128clos4;
00245 #elif DBL_MANT_DIG > 64
00246         IRNG_int i1 = this->irng_.random();
00247         IRNG_int i2 = this->irng_.random();
00248         IRNG_int i3 = this->irng_.random();
00249 
00250         return i1 * norm96closed1 + i2 * norm96closed2
00251           + i3 * norm96closed3;
00252 #elif DBL_MANT_DIG > 32
00253         IRNG_int i1 = this->irng_.random();
00254         IRNG_int i2 = this->irng_.random();
00255 
00256         return i1 * norm64closed1 + i2 * norm64closed2;
00257 #else
00258         IRNG_int i = this->irng_.random();
00259         return i * norm32closed;
00260 #endif
00261 
00262     }
00263 
00264     double getUniform()
00265     { return random(); }
00266 };
00267 
00268 template<typename IRNG, typename stateTag>
00269 class UniformClosed<long double,IRNG,stateTag>
00270   : public IRNGWrapper<IRNG,stateTag> {
00271 
00272 public:
00273     typedef long double T_numtype;
00274 
00275     long double random()
00276     {
00277 #if LDBL_MANT_DIG > 96
00278  #if (LDBL_MANT_DIG > 128) && !defined(BZ_IGNORE_RNG_ERRORS)
00279   #error RNG code assumes long double mantissa is <= 128 bits (not true for your platform).  Use -DBZ_IGNORE_RNG_ERRORS to ignore this warning.
00280  #endif
00281         IRNG_int i1 = this->irng_.random();
00282         IRNG_int i2 = this->irng_.random();
00283         IRNG_int i3 = this->irng_.random();
00284         IRNG_int i4 = this->irng_.random();
00285 
00286         return i1 * norm128clos1 + i2 * norm128clos2 
00287           + i3 * norm128clos3 + i4 * norm128clos4;
00288 #elif LDBL_MANT_DIG > 64
00289         IRNG_int i1 = this->irng_.random();
00290         IRNG_int i2 = this->irng_.random();
00291         IRNG_int i3 = this->irng_.random();
00292 
00293         return i1 * norm96closed1 + i2 * norm96closed2
00294           + i3 * norm96closed3;
00295 #elif LDBL_MANT_DIG > 32
00296         IRNG_int i1 = this->irng_.random();
00297         IRNG_int i2 = this->irng_.random();
00298 
00299         return i1 * norm64closed1 + i2 * norm64closed2;
00300 #else
00301         IRNG_int i = this->irng_.random();
00302         return i * norm32closed;
00303 #endif
00304     }
00305 
00306     long double getUniform()
00307     { return random(); }
00308   
00309 };
00310 
00311 /*****************************************************************************
00312  * UniformOpen generator: uniform random numbers in (0,1).
00313  *****************************************************************************/
00314 
00315 template<typename T = double, typename IRNG = defaultIRNG,
00316     typename stateTag = defaultState>
00317 class UniformOpen : public UniformClosedOpen<T,IRNG,stateTag> {
00318   public:
00319     typedef T T_numtype;
00320 
00321     T random()
00322     {
00323         // Turn a [0,1) into a (0,1) interval by weeding out
00324         // any zeros.
00325         T x;
00326         do {
00327           x = UniformClosedOpen<T,IRNG,stateTag>::random();
00328         } while (x == 0.0L);
00329 
00330         return x;
00331     }
00332 
00333     T getUniform()
00334     { return random(); }
00335 
00336 };
00337 
00338 /*****************************************************************************
00339  * UniformOpenClosed generator: uniform random numbers in (0,1]
00340  *****************************************************************************/
00341 
00342 template<typename T = double, typename IRNG = defaultIRNG,
00343     typename stateTag = defaultState>
00344 class UniformOpenClosed : public UniformClosedOpen<T,IRNG,stateTag> {
00345 
00346 public:
00347     typedef T T_numtype;
00348 
00349     T random()
00350     {
00351         // Antithetic value: taking 1-X where X is [0,1) results
00352         // in a (0,1] distribution.
00353         return 1.0 - UniformClosedOpen<T,IRNG,stateTag>::random();
00354     }
00355 
00356     T getUniform()
00357     { return random(); }
00358 };
00359 
00360 BZ_NAMESPACE_END
00361 
00362 #endif // BZ_RANDOM_UNIFORM_H

Generated on Wed Oct 17 17:52:59 2007 for blitz by  doxygen 1.5.2