Free cookie consent management tool by TermsFeed Policy Generator

source: trunk/sources/HeuristicLab.Random/NormalDistributedRandom.cs @ 143

Last change on this file since 143 was 138, checked in by gkronber, 17 years ago

fixed serialization number-format of double values in DoubleBoundedConstraint and NormalDistributedRandom (ticket #102)

File size: 13.6 KB
Line 
1#region License Information
2/* HeuristicLab
3 * Copyright (C) 2002-2008 Heuristic and Evolutionary Algorithms Laboratory (HEAL)
4 *
5 * This file is part of HeuristicLab.
6 *
7 * HeuristicLab is free software: you can redistribute it and/or modify
8 * it under the terms of the GNU General Public License as published by
9 * the Free Software Foundation, either version 3 of the License, or
10 * (at your option) any later version.
11 *
12 * HeuristicLab is distributed in the hope that it will be useful,
13 * but WITHOUT ANY WARRANTY; without even the implied warranty of
14 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
15 * GNU General Public License for more details.
16 *
17 * You should have received a copy of the GNU General Public License
18 * along with HeuristicLab. If not, see <http://www.gnu.org/licenses/>.
19 */
20#endregion
21
22using System;
23using System.Collections.Generic;
24using System.Text;
25using HeuristicLab.Core;
26using System.Xml;
27using System.Globalization;
28
29namespace HeuristicLab.Random {
30
31  /// <summary>
32  /// Normally distributed random variable.
33  /// Uses the Ziggurat transform to efficiently generate normal distributed random numbers.
34  /// See "The Ziggurat Method for Generating Random Variables" (G. Marsaglia and W.W. Tsang 2000).
35  /// </summary>
36  public class NormalDistributedRandom : ItemBase, IRandom {
37
38    private double mu;
39
40    public double Mu {
41      get { return mu; }
42      set { mu = value; }
43    }
44    private double sigma;
45
46    public double Sigma {
47      get { return sigma; }
48      set { sigma = value; }
49    }
50    private IRandom uniform;
51    private double[] w = new double[] {
52      1.7290404664e-09,
53      1.2680928529e-10,
54      1.6897518107e-10,
55      1.9862687883e-10,
56      2.2232431174e-10,
57      2.4244936614e-10,
58      2.6016130916e-10,
59      2.7611987696e-10,
60      2.9073962682e-10,
61      3.0429969655e-10,
62      3.1699795566e-10,
63      3.2898020419e-10,
64      3.4035738117e-10,
65      3.5121602848e-10,
66      3.6162509098e-10,
67      3.7164057942e-10,
68      3.8130856805e-10,
69      3.9066758162e-10,
70      3.9975012189e-10,
71      4.0858399997e-10,
72      4.1719308563e-10,
73      4.2559822333e-10,
74      4.3381759296e-10,
75      4.4186720949e-10,
76      4.4976131153e-10,
77      4.5751258337e-10,
78      4.6513240481e-10,
79      4.7263104541e-10,
80      4.8001774777e-10,
81      4.8730097735e-10,
82      4.9448850570e-10,
83      5.0158732723e-10,
84      5.0860404777e-10,
85      5.1554460700e-10,
86      5.2241466708e-10,
87      5.2921933502e-10,
88      5.3596349581e-10,
89      5.4265170135e-10,
90      5.4928817050e-10,
91      5.5587695558e-10,
92      5.6242188684e-10,
93      5.6892646150e-10,
94      5.7539412124e-10,
95      5.8182819673e-10,
96      5.8823168558e-10,
97      5.9460769641e-10,
98      6.0095900478e-10,
99      6.0728838625e-10,
100      6.1359850534e-10,
101      6.1989202660e-10,
102      6.2617133700e-10,
103      6.3243904558e-10,
104      6.3869737277e-10,
105      6.4494881657e-10,
106      6.5119559745e-10,
107      6.5744004685e-10,
108      6.6368432972e-10,
109      6.6993072201e-10,
110      6.7618144417e-10,
111      6.8243871665e-10,
112      6.8870464887e-10,
113      6.9498151678e-10,
114      7.0127148533e-10,
115      7.0757677495e-10,
116      7.1389966161e-10,
117      7.2024242126e-10,
118      7.2660727435e-10,
119      7.3299660786e-10,
120      7.3941280876e-10,
121      7.4585826404e-10,
122      7.5233547170e-10,
123      7.5884698525e-10,
124      7.6539541372e-10,
125      7.7198347714e-10,
126      7.7861395109e-10,
127      7.8528972214e-10,
128      7.9201378789e-10,
129      7.9878920145e-10,
130      8.0561923799e-10,
131      8.1250728368e-10,
132      8.1945689123e-10,
133      8.2647166888e-10,
134      8.3355555791e-10,
135      8.4071272166e-10,
136      8.4794732347e-10,
137      8.5526402627e-10,
138      8.6266754851e-10,
139      8.7016316375e-10,
140      8.7775620106e-10,
141      8.8545243360e-10,
142      8.9325818964e-10,
143      9.0117996399e-10,
144      9.0922497309e-10,
145      9.1740082198e-10,
146      9.2571583732e-10,
147      9.3417884539e-10,
148      9.4279972718e-10,
149      9.5158891877e-10,
150      9.6055785548e-10,
151      9.6971930486e-10,
152      9.7908692265e-10,
153      9.8867602993e-10,
154      9.9850361313e-10,
155      1.0085882129e-09,
156      1.0189509236e-09,
157      1.0296150599e-09,
158      1.0406069340e-09,
159      1.0519566329e-09,
160      1.0636980186e-09,
161      1.0758701707e-09,
162      1.0885182755e-09,
163      1.1016947354e-09,
164      1.1154610569e-09,
165      1.1298901814e-09,
166      1.1450695947e-09,
167      1.1611052120e-09,
168      1.1781275955e-09,
169      1.1962995039e-09,
170      1.2158286600e-09,
171      1.2369856250e-09,
172      1.2601323318e-09,
173      1.2857697129e-09,
174      1.3146201905e-09,
175      1.3477839955e-09,
176      1.3870635751e-09,
177      1.4357403044e-09,
178      1.5008658760e-09,
179      1.6030947680e-09
180    };
181    private long[] k = new long[] {
182      1991057938,
183      0,
184      1611602771,
185      1826899878,
186      1918584482,
187      1969227037,
188      2001281515,
189      2023368125,
190      2039498179,
191      2051788381,
192      2061460127,
193      2069267110,
194      2075699398,
195      2081089314,
196      2085670119,
197      2089610331,
198      2093034710,
199      2096037586,
200      2098691595,
201      2101053571,
202      2103168620,
203      2105072996,
204      2106796166,
205      2108362327,
206      2109791536,
207      2111100552,
208      2112303493,
209      2113412330,
210      2114437283,
211      2115387130,
212      2116269447,
213      2117090813,
214      2117856962,
215      2118572919,
216      2119243101,
217      2119871411,
218      2120461303,
219      2121015852,
220      2121537798,
221      2122029592,
222      2122493434,
223      2122931299,
224      2123344971,
225      2123736059,
226      2124106020,
227      2124456175,
228      2124787725,
229      2125101763,
230      2125399283,
231      2125681194,
232      2125948325,
233      2126201433,
234      2126441213,
235      2126668298,
236      2126883268,
237      2127086657,
238      2127278949,
239      2127460589,
240      2127631985,
241      2127793506,
242      2127945490,
243      2128088244,
244      2128222044,
245      2128347141,
246      2128463758,
247      2128572095,
248      2128672327,
249      2128764606,
250      2128849065,
251      2128925811,
252      2128994934,
253      2129056501,
254      2129110560,
255      2129157136,
256      2129196237,
257      2129227847,
258      2129251929,
259      2129268426,
260      2129277255,
261      2129278312,
262      2129271467,
263      2129256561,
264      2129233410,
265      2129201800,
266      2129161480,
267      2129112170,
268      2129053545,
269      2128985244,
270      2128906855,
271      2128817916,
272      2128717911,
273      2128606255,
274      2128482298,
275      2128345305,
276      2128194452,
277      2128028813,
278      2127847342,
279      2127648860,
280      2127432031,
281      2127195339,
282      2126937058,
283      2126655214,
284      2126347546,
285      2126011445,
286      2125643893,
287      2125241376,
288      2124799783,
289      2124314271,
290      2123779094,
291      2123187386,
292      2122530867,
293      2121799464,
294      2120980787,
295      2120059418,
296      2119015917,
297      2117825402,
298      2116455471,
299      2114863093,
300      2112989789,
301      2110753906,
302      2108037662,
303      2104664315,
304      2100355223,
305      2094642347,
306      2086670106,
307      2074676188,
308      2054300022,
309      2010539237
310   };
311    private double[] f = new double[] {
312      1.0000000000e+00,
313      9.6359968185e-01,
314      9.3628269434e-01,
315      9.1304361820e-01,
316      8.9228165150e-01,
317      8.7324303389e-01,
318      8.5550057888e-01,
319      8.3878362179e-01,
320      8.2290720940e-01,
321      8.0773830414e-01,
322      7.9317700863e-01,
323      7.7914607525e-01,
324      7.6558417082e-01,
325      7.5244158506e-01,
326      7.3967725039e-01,
327      7.2725689411e-01,
328      7.1515148878e-01,
329      7.0333611965e-01,
330      6.9178915024e-01,
331      6.8049186468e-01,
332      6.6942769289e-01,
333      6.5858197212e-01,
334      6.4794182777e-01,
335      6.3749545813e-01,
336      6.2723249197e-01,
337      6.1714339256e-01,
338      6.0721951723e-01,
339      5.9745317698e-01,
340      5.8783704042e-01,
341      5.7836467028e-01,
342      5.6902998686e-01,
343      5.5982738733e-01,
344      5.5075180531e-01,
345      5.4179835320e-01,
346      5.3296267986e-01,
347      5.2424055338e-01,
348      5.1562821865e-01,
349      5.0712203979e-01,
350      4.9871864915e-01,
351      4.9041482806e-01,
352      4.8220765591e-01,
353      4.7409430146e-01,
354      4.6607214212e-01,
355      4.5813870430e-01,
356      4.5029163361e-01,
357      4.4252872467e-01,
358      4.3484783173e-01,
359      4.2724698782e-01,
360      4.1972434521e-01,
361      4.1227802634e-01,
362      4.0490642190e-01,
363      3.9760786295e-01,
364      3.9038079977e-01,
365      3.8322380185e-01,
366      3.7613546848e-01,
367      3.6911445856e-01,
368      3.6215949059e-01,
369      3.5526937246e-01,
370      3.4844297171e-01,
371      3.4167915583e-01,
372      3.3497685194e-01,
373      3.2833510637e-01,
374      3.2175290585e-01,
375      3.1522938609e-01,
376      3.0876362324e-01,
377      3.0235484242e-01,
378      2.9600214958e-01,
379      2.8970485926e-01,
380      2.8346219659e-01,
381      2.7727350593e-01,
382      2.7113807201e-01,
383      2.6505529881e-01,
384      2.5902456045e-01,
385      2.5304529071e-01,
386      2.4711695313e-01,
387      2.4123899639e-01,
388      2.3541094363e-01,
389      2.2963231802e-01,
390      2.2390270233e-01,
391      2.1822164953e-01,
392      2.1258877218e-01,
393      2.0700371265e-01,
394      2.0146611333e-01,
395      1.9597564638e-01,
396      1.9053204358e-01,
397      1.8513499200e-01,
398      1.7978426814e-01,
399      1.7447963357e-01,
400      1.6922089458e-01,
401      1.6400785744e-01,
402      1.5884037316e-01,
403      1.5371830761e-01,
404      1.4864157140e-01,
405      1.4361007512e-01,
406      1.3862377405e-01,
407      1.3368265331e-01,
408      1.2878671288e-01,
409      1.2393598258e-01,
410      1.1913054436e-01,
411      1.1437050998e-01,
412      1.0965602100e-01,
413      1.0498725623e-01,
414      1.0036443919e-01,
415      9.5787845552e-02,
416      9.1257803142e-02,
417      8.6774669588e-02,
418      8.2338899374e-02,
419      7.7950984240e-02,
420      7.3611505330e-02,
421      6.9321118295e-02,
422      6.5080583096e-02,
423      6.0890771449e-02,
424      5.6752663106e-02,
425      5.2667401731e-02,
426      4.8636294901e-02,
427      4.4660862535e-02,
428      4.0742866695e-02,
429      3.6884389818e-02,
430      3.3087886870e-02,
431      2.9356317595e-02,
432      2.5693291798e-02,
433      2.2103304043e-02,
434      1.8592102453e-02,
435      1.5167297795e-02,
436      1.1839478277e-02,
437      8.6244847625e-03,
438      5.5489949882e-03,
439      2.6696291752e-03
440    };
441
442    public NormalDistributedRandom() {
443      this.mu = 0.0;
444      this.sigma = 1.0;
445      this.uniform = new MersenneTwister();
446    }
447
448    public NormalDistributedRandom(IRandom uniformRandom, double mu, double sigma) {
449      this.mu = mu;
450      this.sigma = sigma;
451      this.uniform = uniformRandom;
452    }
453
454    #region IRandom Members
455
456    public void Reset() {
457      uniform.Reset();
458    }
459
460    public void Reset(int seed) {
461      uniform.Reset(seed);
462    }
463
464    public int Next() {
465      throw new Exception("The method or operation is not implemented.");
466    }
467
468    public int Next(int maxVal) {
469      throw new Exception("The method or operation is not implemented.");
470    }
471
472    public int Next(int minVal, int maxVal) {
473      throw new Exception("The method or operation is not implemented.");
474    }
475
476    public double NextDouble() {
477      double signFactor = uniform.Next()%2==0?1.0:-1.0;
478      return sigma* signFactor * NextPositiveDouble() + mu;
479    }
480
481    private double NextPositiveDouble() {
482      int j = uniform.Next();
483      int i = (j & 127);
484      if(Math.Abs(j) < k[i]) {
485        return j * w[i];
486      } else {
487        double r = 3.442620;
488        double x, y;
489        x = j * w[i];
490        if(i == 0) {
491          do {
492            x = -Math.Log(ScaledUniform()) * 0.2904764;
493            y = -Math.Log(ScaledUniform());
494          } while(y + y < x * x);
495          return (j > 0) ? r + x : -r - x;
496        }
497        if(f[i] + ScaledUniform() * (f[i - 1] - f[i]) < Math.Exp(-0.5 * x * x)) {
498          return x;
499        } else {
500          // recurse
501          return (NextPositiveDouble());
502        }
503      }
504    }
505
506    private double ScaledUniform() {
507      return 0.5 + uniform.Next() * 0.2328306e-9;
508    }
509
510
511    #endregion
512
513    #region persistence
514    public override XmlNode GetXmlNode(string name, XmlDocument document, IDictionary<Guid, IStorable> persistedObjects) {
515      XmlNode node = base.GetXmlNode(name, document, persistedObjects);
516
517      XmlNode muNode = document.CreateNode(XmlNodeType.Element, "Mu", null);
518      muNode.InnerText = mu.ToString(CultureInfo.InvariantCulture);
519      node.AppendChild(muNode);
520
521      XmlNode sigmaNode = document.CreateNode(XmlNodeType.Element, "Sigma", null);
522      sigmaNode.InnerText = sigma.ToString(CultureInfo.InvariantCulture);
523      node.AppendChild(sigmaNode);
524
525      node.AppendChild(PersistenceManager.Persist("UniformRandom", uniform, document, persistedObjects));
526
527      return node;
528    }
529
530    public override void Populate(XmlNode node, IDictionary<Guid, IStorable> restoredObjects) {
531      base.Populate(node, restoredObjects);
532
533      mu = double.Parse(node.SelectSingleNode("Mu").InnerText, CultureInfo.InvariantCulture);
534      sigma = double.Parse(node.SelectSingleNode("Sigma").InnerText, CultureInfo.InvariantCulture);
535      uniform = (IRandom)PersistenceManager.Restore(node.SelectSingleNode("UniformRandom"), restoredObjects);
536    }
537
538    public override object Clone(IDictionary<Guid, object> clonedObjects) {
539      NormalDistributedRandom clone = new NormalDistributedRandom((IRandom)Auxiliary.Clone(uniform, clonedObjects), mu, sigma);
540      clonedObjects.Add(Guid, clone);
541      return clone;
542    }
543
544    #endregion
545  }
546}
Note: See TracBrowser for help on using the repository browser.