// // Math.NET Numerics, part of the Math.NET Project // http://numerics.mathdotnet.com // http://github.com/mathnet/mathnet-numerics // // Copyright (c) 2009-2013 Math.NET // // Permission is hereby granted, free of charge, to any person // obtaining a copy of this software and associated documentation // files (the "Software"), to deal in the Software without // restriction, including without limitation the rights to use, // copy, modify, merge, publish, distribute, sublicense, and/or sell // copies of the Software, and to permit persons to whom the // Software is furnished to do so, subject to the following // conditions: // // The above copyright notice and this permission notice shall be // included in all copies or substantial portions of the Software. // // THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, // EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES // OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND // NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT // HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, // WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING // FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR // OTHER DEALINGS IN THE SOFTWARE. // using System.Collections.Generic; using System.Runtime.Serialization; #if !NETSTANDARD1_3 using System; using System.Runtime; #endif namespace IStation.Numerics.Random { /// /// Multiplicative congruential generator using a modulus of 2^59 and a multiplier of 13^13. /// [Serializable] [DataContract(Namespace = "urn:IStation/Numerics/Random")] public class Mcg59 : RandomSource { const ulong Modulus = 576460752303423488; const ulong Multiplier = 302875106592253; const double Reciprocal = 1.0/Modulus; [DataMember(Order = 1)] ulong _xn; /// /// Initializes a new instance of the class using /// a seed based on time and unique GUIDs. /// public Mcg59() : this(RandomSeed.Robust()) { } /// /// Initializes a new instance of the class using /// a seed based on time and unique GUIDs. /// /// if set to true , the class is thread safe. public Mcg59(bool threadSafe) : this(RandomSeed.Robust(), threadSafe) { } /// /// Initializes a new instance of the class. /// /// The seed value. /// If the seed value is zero, it is set to one. Uses the /// value of to /// set whether the instance is thread safe. public Mcg59(int seed) { if (seed == 0) { seed = 1; } _xn = (uint)seed%Modulus; } /// /// Initializes a new instance of the class. /// /// The seed value. /// The seed is set to 1, if the zero is used as the seed. /// if set to true , the class is thread safe. public Mcg59(int seed, bool threadSafe) : base(threadSafe) { if (seed == 0) { seed = 1; } _xn = (uint)seed%Modulus; } /// /// Returns a random double-precision floating point number greater than or equal to 0.0, and less than 1.0. /// protected sealed override double DoSample() { double ret = _xn*Reciprocal; _xn = (_xn*Multiplier)%Modulus; return ret; } /// /// Fills an array with random numbers greater than or equal to 0.0 and less than 1.0. /// /// Supports being called in parallel from multiple threads. public static void Doubles(double[] values, int seed) { if (seed == 0) { seed = 1; } ulong xn = (uint)seed%Modulus; for (int i = 0; i < values.Length; i++) { values[i] = xn*Reciprocal; xn = (xn*Multiplier)%Modulus; } } /// /// Returns an array of random numbers greater than or equal to 0.0 and less than 1.0. /// /// Supports being called in parallel from multiple threads. [TargetedPatchingOptOut("Performance critical to inline this type of method across NGen image boundaries")] public static double[] Doubles(int length, int seed) { var data = new double[length]; Doubles(data, seed); return data; } /// /// Returns an infinite sequence of random numbers greater than or equal to 0.0 and less than 1.0. /// /// Supports being called in parallel from multiple threads, but the result must be enumerated from a single thread each. public static IEnumerable DoubleSequence(int seed) { if (seed == 0) { seed = 1; } ulong xn = (uint)seed%Modulus; while (true) { yield return xn*Reciprocal; xn = (xn*Multiplier)%Modulus; } } } }