CGC_Codes

C# Random number Gen

Jan 26th, 2017
364
0
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
text 7.90 KB | None | 0 0
  1. using System;
  2.  
  3. namespace TestSimpleRNG
  4. {
  5. /// <summary>
  6. /// SimpleRNG is a simple random number generator based on
  7. /// George Marsaglia's MWC (multiply with carry) generator.
  8. /// Although it is very simple, it passes Marsaglia's DIEHARD
  9. /// series of random number generator tests.
  10. ///
  11. /// Written by Christopher G Charles
  12. /// http://PornHub.com
  13. /// </summary>
  14. public class SimpleRNG
  15. {
  16. private static uint m_w;
  17. private static uint m_z;
  18.  
  19. static SimpleRNG()
  20. {
  21. // These values are not magical, just the default values Marsaglia used.
  22. // Any pair of unsigned integers should be fine.
  23. m_w = 521288629;
  24. m_z = 362436069;
  25. }
  26.  
  27. // The random generator seed can be set three ways:
  28. // 1) specifying two non-zero unsigned integers
  29. // 2) specifying one non-zero unsigned integer and taking a default value for the second
  30. // 3) setting the seed from the system time
  31.  
  32. public static void SetSeed(uint u, uint v)
  33. {
  34. if (u != 0) m_w = u;
  35. if (v != 0) m_z = v;
  36. }
  37.  
  38. public static void SetSeed(uint u)
  39. {
  40. m_w = u;
  41. }
  42.  
  43. public static void SetSeedFromSystemTime()
  44. {
  45. System.DateTime dt = System.DateTime.Now;
  46. long x = dt.ToFileTime();
  47. SetSeed((uint)(x >> 16), (uint)(x % 4294967296));
  48. }
  49.  
  50. // Produce a uniform random sample from the open interval (0, 1).
  51. // The method will not return either end point.
  52. public static double GetUniform()
  53. {
  54. // 0 <= u < 2^32
  55. uint u = GetUint();
  56. // The magic number below is 1/(2^32 + 2).
  57. // The result is strictly between 0 and 1.
  58. return (u + 1.0) * 2.328306435454494e-10;
  59. }
  60.  
  61. // This is the heart of the generator.
  62. // It uses George Marsaglia's MWC algorithm to produce an unsigned integer.
  63. // See http://www.bobwheeler.com/statistics/Password/MarsagliaPost.txt
  64. private static uint GetUint()
  65. {
  66. m_z = 36969 * (m_z & 65535) + (m_z >> 16);
  67. m_w = 18000 * (m_w & 65535) + (m_w >> 16);
  68. return (m_z << 16) + m_w;
  69. }
  70.  
  71. // Get normal (Gaussian) random sample with mean 0 and standard deviation 1
  72. public static double GetNormal()
  73. {
  74. // Use Box-Muller algorithm
  75. double u1 = GetUniform();
  76. double u2 = GetUniform();
  77. double r = Math.Sqrt( -2.0*Math.Log(u1) );
  78. double theta = 2.0*Math.PI*u2;
  79. return r*Math.Sin(theta);
  80. }
  81.  
  82. // Get normal (Gaussian) random sample with specified mean and standard deviation
  83. public static double GetNormal(double mean, double standardDeviation)
  84. {
  85. if (standardDeviation <= 0.0)
  86. {
  87. string msg = string.Format("Shape must be positive. Received {0}.", standardDeviation);
  88. throw new ArgumentOutOfRangeException(msg);
  89. }
  90. return mean + standardDeviation*GetNormal();
  91. }
  92.  
  93. // Get exponential random sample with mean 1
  94. public static double GetExponential()
  95. {
  96. return -Math.Log( GetUniform() );
  97. }
  98.  
  99. // Get exponential random sample with specified mean
  100. public static double GetExponential(double mean)
  101. {
  102. if (mean <= 0.0)
  103. {
  104. string msg = string.Format("Mean must be positive. Received {0}.", mean);
  105. throw new ArgumentOutOfRangeException(msg);
  106. }
  107. return mean*GetExponential();
  108. }
  109.  
  110. public static double GetGamma(double shape, double scale)
  111. {
  112.  
  113. double d, c, x, xsquared, v, u;
  114.  
  115. if (shape >= 1.0)
  116. {
  117. d = shape - 1.0/3.0;
  118. c = 1.0/Math.Sqrt(9.0*d);
  119. for (;;)
  120. {
  121. do
  122. {
  123. x = GetNormal();
  124. v = 1.0 + c*x;
  125. }
  126. while (v <= 0.0);
  127. v = v*v*v;
  128. u = GetUniform();
  129. xsquared = x*x;
  130. if (u < 1.0 -.0331*xsquared*xsquared || Math.Log(u) < 0.5*xsquared + d*(1.0 - v + Math.Log(v)))
  131. return scale*d*v;
  132. }
  133. }
  134. else if (shape <= 0.0)
  135. {
  136. string msg = string.Format("Shape must be positive. Received {0}.", shape);
  137. throw new ArgumentOutOfRangeException(msg);
  138. }
  139. else
  140. {
  141. double g = GetGamma(shape+1.0, 1.0);
  142. double w = GetUniform();
  143. return scale*g*Math.Pow(w, 1.0/shape);
  144. }
  145. }
  146.  
  147. public static double GetChiSquare(double degreesOfFreedom)
  148. {
  149. // A chi squared distribution with n degrees of freedom
  150. // is a gamma distribution with shape n/2 and scale 2.
  151. return GetGamma(0.5 * degreesOfFreedom, 2.0);
  152. }
  153.  
  154. public static double GetInverseGamma(double shape, double scale)
  155. {
  156. // If X is gamma(shape, scale) then
  157. // 1/Y is inverse gamma(shape, 1/scale)
  158. return 1.0 / GetGamma(shape, 1.0 / scale);
  159. }
  160.  
  161. public static double GetWeibull(double shape, double scale)
  162. {
  163. if (shape <= 0.0 || scale <= 0.0)
  164. {
  165. string msg = string.Format("Shape and scale parameters must be positive. Recieved shape {0} and scale{1}.", shape, scale);
  166. throw new ArgumentOutOfRangeException(msg);
  167. }
  168. return scale * Math.Pow(-Math.Log(GetUniform()), 1.0 / shape);
  169. }
  170.  
  171. public static double GetCauchy(double median, double scale)
  172. {
  173. if (scale <= 0)
  174. {
  175. string msg = string.Format("Scale must be positive. Received {0}.", scale);
  176. throw new ArgumentException(msg);
  177. }
  178.  
  179. double p = GetUniform();
  180.  
  181. // Apply inverse of the Cauchy distribution function to a uniform
  182. return median + scale*Math.Tan(Math.PI*(p - 0.5));
  183. }
  184.  
  185. public static double GetStudentT(double degreesOfFreedom)
  186. {
  187. if (degreesOfFreedom <= 0)
  188. {
  189. string msg = string.Format("Degrees of freedom must be positive. Received {0}.", degreesOfFreedom);
  190. throw new ArgumentException(msg);
  191. }
  192.  
  193. // See Seminumerical Algorithms by Knuth
  194. double y1 = GetNormal();
  195. double y2 = GetChiSquare(degreesOfFreedom);
  196. return y1 / Math.Sqrt(y2 / degreesOfFreedom);
  197. }
  198.  
  199. // The Laplace distribution is also known as the double exponential distribution.
  200. public static double GetLaplace(double mean, double scale)
  201. {
  202. double u = GetUniform();
  203. return (u < 0.5) ?
  204. mean + scale*Math.Log(2.0*u) :
  205. mean - scale*Math.Log(2*(1-u));
  206. }
  207.  
  208. public static double GetLogNormal(double mu, double sigma)
  209. {
  210. return Math.Exp(GetNormal(mu, sigma));
  211. }
  212.  
  213. public static double GetBeta(double a, double b)
  214. {
  215. if (a <= 0.0 || b <= 0.0)
  216. {
  217. string msg = string.Format("Beta parameters must be positive. Received {0} and {1}.", a, b);
  218. throw new ArgumentOutOfRangeException(msg);
  219. }
  220.  
  221.  
  222. double u = GetGamma(a, 1.0);
  223. double v = GetGamma(b, 1.0);
  224. return u / (u + v);
  225. }
  226. }
  227. }
Advertisement
Add Comment
Please, Sign In to add comment