Scopira 20080306

linconrandom.h

00001 
00002 /*
00003  *  Copyright (c) 2002-2006    National Research Council
00004  *
00005  *  All rights reserved.
00006  *
00007  *  This material is confidential and proprietary information of
00008  *  National Research Council Canada ("Confidential Information").
00009  *  This Confidential Information may only be used and reproduced
00010  *  in accordance with the terms of the license agreement.
00011  *
00012  */
00013 
00014 #ifndef __INCLUDED__SCOPIRA_TOOL_LINCONRANDOM_H__
00015 #define __INCLUDED__SCOPIRA_TOOL_LINCONRANDOM_H__
00016 
00017 #include <assert.h>
00018 
00019 #include <scopira/tool/constmod.h>
00020 #include <scopira/tool/platform.h>
00021 
00022 namespace scopira
00023 {
00024   namespace tool
00025   {
00026     template<class IntType, IntType a, IntType c, IntType m, IntType val>
00027       class lincon_gen;
00028 
00031     typedef lincon_gen<int, 16807, 0, 2147483647, 1043618065> minstd_rand0;
00034     typedef lincon_gen<int, 48271, 0, 2147483647, 399268537> minstd_rand;
00035 
00037     typedef minstd_rand0 lincon_core;
00038 
00039     class rand48;
00040   }
00041 }
00042 
00050 template<class IntType, IntType a, IntType c, IntType m, IntType val>
00051 class scopira::tool::lincon_gen
00052 {
00053   public:
00054     typedef IntType result_type;
00055 
00056     static const bool has_fixed_range = true;
00057     static const result_type min_value = ( c == 0 ? 1 : 0 );
00058     static const result_type max_value = m-1;
00059 
00061     explicit lincon_gen(IntType x0 = 1)
00062       : _x(x0) { 
00063       assert(c || x0); // if c == 0 and x(0) == 0 then x(n) = 0 for all n
00064       // overflow check
00065       // disabled because it gives spurious "divide by zero" gcc warnings
00066       // assert(m == 0 || (a*(m-1)+c) % m == (c < a ? c-a+m : c-a)); 
00067     }
00068 
00070     void seed(IntType x0) { assert(c || x0); _x = x0; }
00071     result_type min(void) const { return c == 0 ? 1 : 0; }
00072     result_type max(void) const { return m-1; }
00073 
00075     IntType operator()(void) {
00076       _x = const_mod<IntType, m>::mult_add(a, _x, c);
00077       return _x;
00078     }
00079 
00080     bool validation(IntType x) const { return val == x; }
00081   private:
00082     IntType _x;
00083 };
00084 
00093 class scopira::tool::rand48 
00094 {
00095   public:
00096     typedef int32_t result_type;
00097     static const bool has_fixed_range = true;
00098     static const int32_t min_value = 0;
00099     //static const int32_t max_value = integer_traits<int32_t>::const_max;
00100 
00101     int32_t min(void) const { return 0; }
00102     int32_t max(void) const { return std::numeric_limits<int32_t>::max(); }
00103     
00104     explicit rand48(int32_t x0 = 1) : lcf(cnv(x0)) { }
00105     explicit rand48(uint64_t x0) : lcf(x0) { }
00106     //template<class It> rand48(It& first, It last) : lcf(first, last) { }
00107 
00108     // compiler-generated copy ctor and assignment operator are fine
00109 
00110     void seed(int32_t x0 = 1) { lcf.seed(cnv(x0)); }
00111     void seed(uint64_t x0) { lcf.seed(x0); }
00112     //template<class It> void seed(It& first, It last) { lcf.seed(first,last); }
00113 
00114     int32_t operator()(void) { return static_cast<int32_t>(lcf() >> 17); }
00115 
00116   private:
00117     scopira::tool::lincon_gen<uint64_t, uint64_t(0xDEECE66DUL) | (uint64_t(0x5) << 32),
00118       0xB, uint64_t(1)<<48, 0> lcf;
00119 
00120   private:
00121     static uint64_t cnv(int32_t x) { return (static_cast<uint64_t>(x) << 16) | 0x330e;  }
00122 };
00123 
00124 #endif
00125