
今天要介紹的是Math.NET中擴充的其他随機數生成算法。真正意義上的随機數(或者随機事件)在某次産生過程中是按照實驗過程中表現的分布機率随機産生的,其結果是不可預測的,是不可見的。而計算機中的随機函數是按照一定算法模拟産生的,其結果是确定的,是可見的。我們可以這樣認為這個可預見的結果其出現的機率是100%。是以用計算機随機函數所産生的“随機數”并不随機,是僞随機數。僞随機數的作用在開發中的使用非常常見,是以.NET在System命名空間,提供了一個簡單的Random随機數生成類型。但這個類型并不能滿足所有的需求,本節開始就将陸續介紹Math.NET中有關随機數的擴充以及其他僞随機生成算法編寫的随機數生成器。
本部落格所有文章分類的總目錄:【總目錄】本部落格博文總目錄-實時更新
開源Math.NET基礎數學類庫使用總目錄:【目錄】開源Math.NET基礎數學類庫使用總目錄
前言
真正意義上的随機數(或者随機事件)在某次産生過程中是按照實驗過程中表現的分布機率随機産生的,其結果是不可預測的,是不可見的。而計算機中的随機函數是按照一定算法模拟産生的,其結果是确定的,是可見的。我們可以這樣認為這個可預見的結果其出現的機率是100%。是以用計算機随機函數所産生的“随機數”并不随機,是僞随機數。僞随機數的作用在開發中的使用非常常見,是以.NET在System命名空間,提供了一個簡單的Random随機數生成類型。但這個類型并不能滿足所有的需求,本節開始就将陸續介紹Math.NET中有關随機數的擴充以及其他僞随機生成算法編寫的随機數生成器。
今天要介紹的是Math.NET中擴充的其他随機數生成算法。
如果本文資源或者顯示有問題,請參考 本文原文位址:http://www.cnblogs.com/asxinyu/p/4301555.html
http://zh.wikipedia.org/wiki/随機數
随機數是專門的随機試驗的結果。
在統計學的不同技術中需要使用随機數,比如在從統計總體中抽取有代表性的樣本的時候,或者在将實驗動物配置設定到不同的試驗組的過程中,或者在進行蒙特卡羅模拟法計算的時候等等。産生随機數有多種不同的方法。這些方法被稱為随機數生成器。随機數最重要的特性是它在産生時後面的那個數與前面的那個數毫無關系。
真正的随機數是使用實體現象産生的:比如擲錢币、骰子、轉輪、使用電子元件的噪音、核裂變等等。這樣的随機數生成器叫做實體性随機數生成器,它們的缺點是技術要求比較高。在實際應用中往往使用僞随機數就足夠了。這些數列是“似乎”随機的數,實際上它們是通過一個固定的、可以重複的計算方法産生的。它們不真正地随機,因為它們實際上是可以計算出來的,但是它們具有類似于随機數的統計特征。這樣的生成器叫做僞随機數生成器。在真正關鍵性的應用中,比如在密碼學中,人們一般使用真正的随機數。
1.Math.NET的其他随機數生成算法
Math.NET在MathNet.Numerics.Random命名空間下包括了若幹個其他随機數生成算法,基本介紹如下:
1.Mcg31m1類 與 Mcg59 類,都屬于 矩陣同餘發生器,矩陣同餘發生器是乘同餘線性發生器的一個推廣;2者的參數有些差别;
3.MersenneTwister,Mersenne Twister算法譯為馬特賽特旋轉演算法,是僞随機數發生器之一,其主要作用是生成僞随機數。此算法是Makoto Matsumoto (松本)和Takuji Nishimura (西村)于1997年開發的,基于有限二進制字段上的矩陣線性再生。可以快速産生高品質的僞随機數,修正了古老随機數産生算法的很多缺陷。 Mersenne Twister這個名字來自周期長度通常取Mersenne質數這樣一個事實。常見的有兩個變種Mersenne Twister MT19937和Mersenne Twister MT19937-64。Math.NET實作的版本是前者(Mersenne Twister MT19937)。介紹
4.Mrg32k3a,又叫 素數模乘同餘法,素數模乘同餘發生器是提出的一個統計性質較好和周期較大的發生器,目前是使用最廣的一種均勻随機數發生器下面給出兩組經檢驗統計性質是良好的素數模乘同餘發生器;
5.Palf,是一個并行加法滞後的斐波那契僞随機數字生成器。
6.WH1982與WH2006類,是乘線性同餘法,也叫積式發生器,Math.NET實作的2個版本友善是1982和2006,代表作者釋出該算法論文的2個時間,可以參考作者的2篇論文:
1.Wichmann, B. A. & Hill, I. D. (1982), "Algorithm AS 183:An efficient and portable pseudo-random number generator". Applied Statistics 31 (1982) 188-190
2. Wichmann, B. A. & Hill, I. D. (2006), "Generating good pseudo-random numbers".Computational Statistics & Data Analysis 51:3 (2006) 1614-1622
8.Xorshift類,實作的是George Marsaglia,在論文“Xorshift RNGs”提出的一個算法,原理可以參考這篇論文:http://www.jstatsoft.org/v08/i14/paper;
2.Math.NET擴充随機數生成算法的實作
上面已經已經對Math.NET擴充的幾個随機數生成算法進行了介紹,對于一般人來說,直接用就可以了,但對于特殊的人來說,可能要用到其中一種,可以直接使用C#進行調用即可,當然為了節省大家的時間,這裡對Math.NET的實作做一個簡單的介紹,這樣大家可以更加快速的擴充自己的算法,同時也可以了解實作的原理,對Math.NET有一個更加深入的了解。
Math.NET對随機數的擴充也是借用了System.Random類,在它的基礎上實作了一個随機數發生器的基類:RandomSource,所有的擴充算法都實作該類,這樣使用RandomSource就非常友善。RandomSource的結構很簡單,對System.Random進行了簡單的封裝,增加了幾個直接生成其他類型随機數的方法,并都可以在繼承中使用。其源碼如下:
1 public abstract class RandomSource : System.Random
2 {
3 readonly bool _threadSafe;
4 readonly object _lock = new object();
5
6 /// <summary>
7 /// Initializes a new instance of the <see cref="RandomSource"/> class using
8 /// the value of <see cref="Control.ThreadSafeRandomNumberGenerators"/> to set whether
9 /// the instance is thread safe or not.
10 /// </summary>
11 protected RandomSource() : base(RandomSeed.Robust())
12 {
13 _threadSafe = Control.ThreadSafeRandomNumberGenerators;
14 }
15
16 /// <summary>
17 /// Initializes a new instance of the <see cref="RandomSource"/> class.
18 /// </summary>
19 /// <param name="threadSafe">if set to <c>true</c> , the class is thread safe.</param>
20 /// <remarks>Thread safe instances are two and half times slower than non-thread
21 /// safe classes.</remarks>
22 protected RandomSource(bool threadSafe) : base(RandomSeed.Robust())
23 {
24 _threadSafe = threadSafe;
25 }
26
27 /// <summary>
28 /// Fills an array with uniform random numbers greater than or equal to 0.0 and less than 1.0.
29 /// </summary>
30 /// <param name="values">The array to fill with random values.</param>
31 public void NextDoubles(double[] values)
32 {
33 if (_threadSafe)
34 {
35 lock (_lock)
36 {
37 for (var i = 0; i < values.Length; i++)
38 {
39 values[i] = DoSample();
40 }
41 }
42 }
43 else
44 {
45 for (var i = 0; i < values.Length; i++)
46 {
47 values[i] = DoSample();
48 }
49 }
50 }
51
52 /// <summary>
53 /// Returns an infinite sequence of uniform random numbers greater than or equal to 0.0 and less than 1.0.
54 /// </summary>
55 public IEnumerable<double> NextDoubleSequence()
56 {
57 for (int i = 0; i < 64; i++)
58 {
59 yield return NextDouble();
60 }
61
62 var buffer = new double[64];
63 while (true)
64 {
65 NextDoubles(buffer);
66 for (int i = 0; i < buffer.Length; i++)
67 {
68 yield return buffer[i];
69 }
70 }
71 }
72
73 /// <summary>
74 /// Returns a nonnegative random number.
75 /// </summary>
76 /// <returns>
77 /// A 32-bit signed integer greater than or equal to zero and less than <see cref="F:System.Int32.MaxValue"/>.
78 /// </returns>
79 public override sealed int Next()
80 {
81 if (_threadSafe)
82 {
83 lock (_lock)
84 {
85 return (int)(DoSample()*int.MaxValue);
86 }
87 }
88
89 return (int)(DoSample()*int.MaxValue);
90 }
91
92 /// <summary>
93 /// Returns a random number less then a specified maximum.
94 /// </summary>
95 /// <param name="maxValue">The exclusive upper bound of the random number returned.</param>
96 /// <returns>A 32-bit signed integer less than <paramref name="maxValue"/>.</returns>
97 /// <exception cref="T:System.ArgumentOutOfRangeException"><paramref name="maxValue"/> is negative. </exception>
98 public override sealed int Next(int maxValue)
99 {
100 if (maxValue <= 0)
101 {
102 throw new ArgumentException(Resources.ArgumentMustBePositive);
103 }
104
105 if (_threadSafe)
106 {
107 lock (_lock)
108 {
109 return (int)(DoSample()*maxValue);
110 }
111 }
112
113 return (int)(DoSample()*maxValue);
114 }
115
116 /// <summary>
117 /// Returns a random number within a specified range.
118 /// </summary>
119 /// <param name="minValue">The inclusive lower bound of the random number returned.</param>
120 /// <param name="maxValue">The exclusive upper bound of the random number returned. <paramref name="maxValue"/> must be greater than or equal to <paramref name="minValue"/>.</param>
121 /// <returns>
122 /// A 32-bit signed integer greater than or equal to <paramref name="minValue"/> and less than <paramref name="maxValue"/>; that is, the range of return values includes <paramref name="minValue"/> but not <paramref name="maxValue"/>. If <paramref name="minValue"/> equals <paramref name="maxValue"/>, <paramref name="minValue"/> is returned.
123 /// </returns>
124 /// <exception cref="T:System.ArgumentOutOfRangeException"><paramref name="minValue"/> is greater than <paramref name="maxValue"/>. </exception>
125 public override sealed int Next(int minValue, int maxValue)
126 {
127 if (minValue > maxValue)
128 {
129 throw new ArgumentException(Resources.ArgumentMinValueGreaterThanMaxValue);
130 }
131
132 if (_threadSafe)
133 {
134 lock (_lock)
135 {
136 return (int)(DoSample()*(maxValue - minValue)) + minValue;
137 }
138 }
139
140 return (int)(DoSample()*(maxValue - minValue)) + minValue;
141 }
142
143 /// <summary>
144 /// Fills an array with random numbers within a specified range.
145 /// </summary>
146 /// <param name="values">The array to fill with random values.</param>
147 /// <param name="minValue">The inclusive lower bound of the random number returned.</param>
148 /// <param name="maxValue">The exclusive upper bound of the random number returned. <paramref name="maxValue"/> must be greater than or equal to <paramref name="minValue"/>.</param>
149 public void NextInt32s(int[] values, int minValue, int maxValue)
150 {
151 if (_threadSafe)
152 {
153 lock (_lock)
154 {
155 for (var i = 0; i < values.Length; i++)
156 {
157 values[i] = (int)(DoSample()*(maxValue - minValue)) + minValue;
158 }
159 }
160 }
161 else
162 {
163 for (var i = 0; i < values.Length; i++)
164 {
165 values[i] = (int)(DoSample()*(maxValue - minValue)) + minValue;
166 }
167 }
168 }
169
170 /// <summary>
171 /// Returns an infinite sequence of random numbers within a specified range.
172 /// </summary>
173 /// <param name="minValue">The inclusive lower bound of the random number returned.</param>
174 /// <param name="maxValue">The exclusive upper bound of the random number returned. <paramref name="maxValue"/> must be greater than or equal to <paramref name="minValue"/>.</param>
175 public IEnumerable<int> NextInt32Sequence(int minValue, int maxValue)
176 {
177 for (int i = 0; i < 64; i++)
178 {
179 yield return Next(minValue, maxValue);
180 }
181
182 var buffer = new int[64];
183 while (true)
184 {
185 NextInt32s(buffer, minValue, maxValue);
186 for (int i = 0; i < buffer.Length; i++)
187 {
188 yield return buffer[i];
189 }
190 }
191 }
192
193 /// <summary>
194 /// Fills the elements of a specified array of bytes with random numbers.
195 /// </summary>
196 /// <param name="buffer">An array of bytes to contain random numbers.</param>
197 /// <exception cref="T:System.ArgumentNullException"><paramref name="buffer"/> is null. </exception>
198 public override void NextBytes(byte[] buffer)
199 {
200 if (buffer == null)
201 {
202 throw new ArgumentNullException("buffer");
203 }
204
205 if (_threadSafe)
206 {
207 lock (_lock)
208 {
209 for (var i = 0; i < buffer.Length; i++)
210 {
211 buffer[i] = (byte)(((int)(DoSample()*int.MaxValue))%256);
212 }
213 }
214
215 return;
216 }
217
218 for (var i = 0; i < buffer.Length; i++)
219 {
220 buffer[i] = (byte)(((int)(DoSample()*int.MaxValue))%256);
221 }
222 }
223
224 /// <summary>
225 /// Returns a random number between 0.0 and 1.0.
226 /// </summary>
227 /// <returns>A double-precision floating point number greater than or equal to 0.0, and less than 1.0.</returns>
228 protected override sealed double Sample()
229 {
230 if (_threadSafe)
231 {
232 lock (_lock)
233 {
234 return DoSample();
235 }
236 }
237
238 return DoSample();
239 }
240
241 /// <summary>
242 /// Returns a random number between 0.0 and 1.0.
243 /// </summary>
244 /// <returns>
245 /// A double-precision floating point number greater than or equal to 0.0, and less than 1.0.
246 /// </returns>
247 protected abstract double DoSample();
248 }
View Code
在擴充随機數生成算法的時候,直接繼承該類即可,看一下Mcg59的實作源碼:
1 public class Mcg59 : RandomSource
2 {
3 const ulong Modulus = 576460752303423488;
4 const ulong Multiplier = 302875106592253;
5 const double Reciprocal = 1.0/Modulus;
6 ulong _xn;
7
8 /// <summary>
9 /// Initializes a new instance of the <see cref="Mcg59"/> class using
10 /// a seed based on time and unique GUIDs.
11 /// </summary>
12 public Mcg59() : this(RandomSeed.Robust())
13 {
14 }
15
16 /// <summary>
17 /// Initializes a new instance of the <see cref="Mcg59"/> class using
18 /// a seed based on time and unique GUIDs.
19 /// </summary>
20 /// <param name="threadSafe">if set to <c>true</c> , the class is thread safe.</param>
21 public Mcg59(bool threadSafe) : this(RandomSeed.Robust(), threadSafe)
22 {
23 }
24
25 /// <summary>
26 /// Initializes a new instance of the <see cref="Mcg59"/> class.
27 /// </summary>
28 /// <param name="seed">The seed value.</param>
29 /// <remarks>If the seed value is zero, it is set to one. Uses the
30 /// value of <see cref="Control.ThreadSafeRandomNumberGenerators"/> to
31 /// set whether the instance is thread safe.</remarks>
32 public Mcg59(int seed)
33 {
34 if (seed == 0)
35 {
36 seed = 1;
37 }
38
39 _xn = (uint)seed%Modulus;
40 }
41
42 /// <summary>
43 /// Initializes a new instance of the <see cref="Mcg59"/> class.
44 /// </summary>
45 /// <param name="seed">The seed value.</param>
46 /// <remarks>The seed is set to 1, if the zero is used as the seed.</remarks>
47 /// <param name="threadSafe">if set to <c>true</c> , the class is thread safe.</param>
48 public Mcg59(int seed, bool threadSafe) : base(threadSafe)
49 {
50 if (seed == 0)
51 {
52 seed = 1;
53 }
54
55 _xn = (uint)seed%Modulus;
56 }
57
58 /// <summary>
59 /// Returns a random number between 0.0 and 1.0.
60 /// </summary>
61 /// <returns>
62 /// A double-precision floating point number greater than or equal to 0.0, and less than 1.0.
63 /// </returns>
64 protected override sealed double DoSample()
65 {
66 double ret = _xn*Reciprocal;
67 _xn = (_xn*Multiplier)%Modulus;
68 return ret;
69 }
70
71 /// <summary>
72 /// Fills an array with random numbers greater than or equal to 0.0 and less than 1.0.
73 /// </summary>
74 /// <remarks>Supports being called in parallel from multiple threads.</remarks>
75 public static void Doubles(double[] values, int seed)
76 {
77 if (seed == 0)
78 {
79 seed = 1;
80 }
81
82 ulong xn = (uint)seed%Modulus;
83
84 for (int i = 0; i < values.Length; i++)
85 {
86 values[i] = xn*Reciprocal;
87 xn = (xn*Multiplier)%Modulus;
88 }
89 }
90
91 /// <summary>
92 /// Returns an array of random numbers greater than or equal to 0.0 and less than 1.0.
93 /// </summary>
94 /// <remarks>Supports being called in parallel from multiple threads.</remarks>
95 [TargetedPatchingOptOut("Performance critical to inline this type of method across NGen image boundaries")]
96 public static double[] Doubles(int length, int seed)
97 {
98 var data = new double[length];
99 Doubles(data, seed);
100 return data;
101 }
102
103 /// <summary>
104 /// Returns an infinite sequence of random numbers greater than or equal to 0.0 and less than 1.0.
105 /// </summary>
106 /// <remarks>Supports being called in parallel from multiple threads, but the result must be enumerated from a single thread each.</remarks>
107 public static IEnumerable<double> DoubleSequence(int seed)
108 {
109 if (seed == 0)
110 {
111 seed = 1;
112 }
113
114 ulong xn = (uint)seed%Modulus;
115
116 while (true)
117 {
118 yield return xn*Reciprocal;
119 xn = (xn*Multiplier)%Modulus;
120 }
121 }
122 }
随機數的使用大家都很熟練,和Random類差不多,上述擴充的算法中,也包括了很多 靜态方法,可以直接使用。這裡不再舉例說明。
3.資源
源碼下載下傳:http://www.cnblogs.com/asxinyu/p/4264638.html
如果本文資源或者顯示有問題,請參考 本文原文位址:http://www.cnblogs.com/asxinyu/p/4301555.html
.NET資料挖掘與機器學習,作者部落格:
http://www.cnblogs.com/asxinyu
E-mail:[email protected]