diff options
Diffstat (limited to 'rand/src')
| -rw-r--r-- | rand/src/distributions/exponential.rs | 124 | ||||
| -rw-r--r-- | rand/src/distributions/gamma.rs | 386 | ||||
| -rw-r--r-- | rand/src/distributions/mod.rs | 409 | ||||
| -rw-r--r-- | rand/src/distributions/normal.rs | 201 | ||||
| -rw-r--r-- | rand/src/distributions/range.rs | 241 | ||||
| -rw-r--r-- | rand/src/distributions/ziggurat_tables.rs | 280 | ||||
| -rw-r--r-- | rand/src/jitter.rs | 754 | ||||
| -rw-r--r-- | rand/src/lib.rs | 1214 | ||||
| -rw-r--r-- | rand/src/os.rs | 617 | ||||
| -rw-r--r-- | rand/src/prng/chacha.rs | 321 | ||||
| -rw-r--r-- | rand/src/prng/isaac.rs | 328 | ||||
| -rw-r--r-- | rand/src/prng/isaac64.rs | 340 | ||||
| -rw-r--r-- | rand/src/prng/mod.rs | 51 | ||||
| -rw-r--r-- | rand/src/prng/xorshift.rs | 101 | ||||
| -rw-r--r-- | rand/src/rand_impls.rs | 299 | ||||
| -rw-r--r-- | rand/src/read.rs | 123 | ||||
| -rw-r--r-- | rand/src/reseeding.rs | 229 | ||||
| -rw-r--r-- | rand/src/seq.rs | 337 |
18 files changed, 6355 insertions, 0 deletions
diff --git a/rand/src/distributions/exponential.rs b/rand/src/distributions/exponential.rs new file mode 100644 index 0000000..c3c924c --- /dev/null +++ b/rand/src/distributions/exponential.rs @@ -0,0 +1,124 @@ +// Copyright 2013 The Rust Project Developers. See the COPYRIGHT +// file at the top-level directory of this distribution and at +// http://rust-lang.org/COPYRIGHT. +// +// Licensed under the Apache License, Version 2.0 <LICENSE-APACHE or +// http://www.apache.org/licenses/LICENSE-2.0> or the MIT license +// <LICENSE-MIT or http://opensource.org/licenses/MIT>, at your +// option. This file may not be copied, modified, or distributed +// except according to those terms. + +//! The exponential distribution. + +use {Rng, Rand}; +use distributions::{ziggurat, ziggurat_tables, Sample, IndependentSample}; + +/// A wrapper around an `f64` to generate Exp(1) random numbers. +/// +/// See `Exp` for the general exponential distribution. +/// +/// Implemented via the ZIGNOR variant[1] of the Ziggurat method. The +/// exact description in the paper was adjusted to use tables for the +/// exponential distribution rather than normal. +/// +/// [1]: Jurgen A. Doornik (2005). [*An Improved Ziggurat Method to +/// Generate Normal Random +/// Samples*](http://www.doornik.com/research/ziggurat.pdf). Nuffield +/// College, Oxford +/// +/// # Example +/// +/// ```rust +/// use rand::distributions::exponential::Exp1; +/// +/// let Exp1(x) = rand::random(); +/// println!("{}", x); +/// ``` +#[derive(Clone, Copy, Debug)] +pub struct Exp1(pub f64); + +// This could be done via `-rng.gen::<f64>().ln()` but that is slower. +impl Rand for Exp1 { + #[inline] + fn rand<R:Rng>(rng: &mut R) -> Exp1 { + #[inline] + fn pdf(x: f64) -> f64 { + (-x).exp() + } + #[inline] + fn zero_case<R:Rng>(rng: &mut R, _u: f64) -> f64 { + ziggurat_tables::ZIG_EXP_R - rng.gen::<f64>().ln() + } + + Exp1(ziggurat(rng, false, + &ziggurat_tables::ZIG_EXP_X, + &ziggurat_tables::ZIG_EXP_F, + pdf, zero_case)) + } +} + +/// The exponential distribution `Exp(lambda)`. +/// +/// This distribution has density function: `f(x) = lambda * +/// exp(-lambda * x)` for `x > 0`. +/// +/// # Example +/// +/// ```rust +/// use rand::distributions::{Exp, IndependentSample}; +/// +/// let exp = Exp::new(2.0); +/// let v = exp.ind_sample(&mut rand::thread_rng()); +/// println!("{} is from a Exp(2) distribution", v); +/// ``` +#[derive(Clone, Copy, Debug)] +pub struct Exp { + /// `lambda` stored as `1/lambda`, since this is what we scale by. + lambda_inverse: f64 +} + +impl Exp { + /// Construct a new `Exp` with the given shape parameter + /// `lambda`. Panics if `lambda <= 0`. + #[inline] + pub fn new(lambda: f64) -> Exp { + assert!(lambda > 0.0, "Exp::new called with `lambda` <= 0"); + Exp { lambda_inverse: 1.0 / lambda } + } +} + +impl Sample<f64> for Exp { + fn sample<R: Rng>(&mut self, rng: &mut R) -> f64 { self.ind_sample(rng) } +} +impl IndependentSample<f64> for Exp { + fn ind_sample<R: Rng>(&self, rng: &mut R) -> f64 { + let Exp1(n) = rng.gen::<Exp1>(); + n * self.lambda_inverse + } +} + +#[cfg(test)] +mod test { + use distributions::{Sample, IndependentSample}; + use super::Exp; + + #[test] + fn test_exp() { + let mut exp = Exp::new(10.0); + let mut rng = ::test::rng(); + for _ in 0..1000 { + assert!(exp.sample(&mut rng) >= 0.0); + assert!(exp.ind_sample(&mut rng) >= 0.0); + } + } + #[test] + #[should_panic] + fn test_exp_invalid_lambda_zero() { + Exp::new(0.0); + } + #[test] + #[should_panic] + fn test_exp_invalid_lambda_neg() { + Exp::new(-10.0); + } +} diff --git a/rand/src/distributions/gamma.rs b/rand/src/distributions/gamma.rs new file mode 100644 index 0000000..2806495 --- /dev/null +++ b/rand/src/distributions/gamma.rs @@ -0,0 +1,386 @@ +// Copyright 2013 The Rust Project Developers. See the COPYRIGHT +// file at the top-level directory of this distribution and at +// http://rust-lang.org/COPYRIGHT. +// +// Licensed under the Apache License, Version 2.0 <LICENSE-APACHE or +// http://www.apache.org/licenses/LICENSE-2.0> or the MIT license +// <LICENSE-MIT or http://opensource.org/licenses/MIT>, at your +// option. This file may not be copied, modified, or distributed +// except according to those terms. +// +// ignore-lexer-test FIXME #15679 + +//! The Gamma and derived distributions. + +use self::GammaRepr::*; +use self::ChiSquaredRepr::*; + +use {Rng, Open01}; +use super::normal::StandardNormal; +use super::{IndependentSample, Sample, Exp}; + +/// The Gamma distribution `Gamma(shape, scale)` distribution. +/// +/// The density function of this distribution is +/// +/// ```text +/// f(x) = x^(k - 1) * exp(-x / θ) / (Γ(k) * θ^k) +/// ``` +/// +/// where `Γ` is the Gamma function, `k` is the shape and `θ` is the +/// scale and both `k` and `θ` are strictly positive. +/// +/// The algorithm used is that described by Marsaglia & Tsang 2000[1], +/// falling back to directly sampling from an Exponential for `shape +/// == 1`, and using the boosting technique described in [1] for +/// `shape < 1`. +/// +/// # Example +/// +/// ```rust +/// use rand::distributions::{IndependentSample, Gamma}; +/// +/// let gamma = Gamma::new(2.0, 5.0); +/// let v = gamma.ind_sample(&mut rand::thread_rng()); +/// println!("{} is from a Gamma(2, 5) distribution", v); +/// ``` +/// +/// [1]: George Marsaglia and Wai Wan Tsang. 2000. "A Simple Method +/// for Generating Gamma Variables" *ACM Trans. Math. Softw.* 26, 3 +/// (September 2000), +/// 363-372. DOI:[10.1145/358407.358414](http://doi.acm.org/10.1145/358407.358414) +#[derive(Clone, Copy, Debug)] +pub struct Gamma { + repr: GammaRepr, +} + +#[derive(Clone, Copy, Debug)] +enum GammaRepr { + Large(GammaLargeShape), + One(Exp), + Small(GammaSmallShape) +} + +// These two helpers could be made public, but saving the +// match-on-Gamma-enum branch from using them directly (e.g. if one +// knows that the shape is always > 1) doesn't appear to be much +// faster. + +/// Gamma distribution where the shape parameter is less than 1. +/// +/// Note, samples from this require a compulsory floating-point `pow` +/// call, which makes it significantly slower than sampling from a +/// gamma distribution where the shape parameter is greater than or +/// equal to 1. +/// +/// See `Gamma` for sampling from a Gamma distribution with general +/// shape parameters. +#[derive(Clone, Copy, Debug)] +struct GammaSmallShape { + inv_shape: f64, + large_shape: GammaLargeShape +} + +/// Gamma distribution where the shape parameter is larger than 1. +/// +/// See `Gamma` for sampling from a Gamma distribution with general +/// shape parameters. +#[derive(Clone, Copy, Debug)] +struct GammaLargeShape { + scale: f64, + c: f64, + d: f64 +} + +impl Gamma { + /// Construct an object representing the `Gamma(shape, scale)` + /// distribution. + /// + /// Panics if `shape <= 0` or `scale <= 0`. + #[inline] + pub fn new(shape: f64, scale: f64) -> Gamma { + assert!(shape > 0.0, "Gamma::new called with shape <= 0"); + assert!(scale > 0.0, "Gamma::new called with scale <= 0"); + + let repr = if shape == 1.0 { + One(Exp::new(1.0 / scale)) + } else if shape < 1.0 { + Small(GammaSmallShape::new_raw(shape, scale)) + } else { + Large(GammaLargeShape::new_raw(shape, scale)) + }; + Gamma { repr: repr } + } +} + +impl GammaSmallShape { + fn new_raw(shape: f64, scale: f64) -> GammaSmallShape { + GammaSmallShape { + inv_shape: 1. / shape, + large_shape: GammaLargeShape::new_raw(shape + 1.0, scale) + } + } +} + +impl GammaLargeShape { + fn new_raw(shape: f64, scale: f64) -> GammaLargeShape { + let d = shape - 1. / 3.; + GammaLargeShape { + scale: scale, + c: 1. / (9. * d).sqrt(), + d: d + } + } +} + +impl Sample<f64> for Gamma { + fn sample<R: Rng>(&mut self, rng: &mut R) -> f64 { self.ind_sample(rng) } +} +impl Sample<f64> for GammaSmallShape { + fn sample<R: Rng>(&mut self, rng: &mut R) -> f64 { self.ind_sample(rng) } +} +impl Sample<f64> for GammaLargeShape { + fn sample<R: Rng>(&mut self, rng: &mut R) -> f64 { self.ind_sample(rng) } +} + +impl IndependentSample<f64> for Gamma { + fn ind_sample<R: Rng>(&self, rng: &mut R) -> f64 { + match self.repr { + Small(ref g) => g.ind_sample(rng), + One(ref g) => g.ind_sample(rng), + Large(ref g) => g.ind_sample(rng), + } + } +} +impl IndependentSample<f64> for GammaSmallShape { + fn ind_sample<R: Rng>(&self, rng: &mut R) -> f64 { + let Open01(u) = rng.gen::<Open01<f64>>(); + + self.large_shape.ind_sample(rng) * u.powf(self.inv_shape) + } +} +impl IndependentSample<f64> for GammaLargeShape { + fn ind_sample<R: Rng>(&self, rng: &mut R) -> f64 { + loop { + let StandardNormal(x) = rng.gen::<StandardNormal>(); + let v_cbrt = 1.0 + self.c * x; + if v_cbrt <= 0.0 { // a^3 <= 0 iff a <= 0 + continue + } + + let v = v_cbrt * v_cbrt * v_cbrt; + let Open01(u) = rng.gen::<Open01<f64>>(); + + let x_sqr = x * x; + if u < 1.0 - 0.0331 * x_sqr * x_sqr || + u.ln() < 0.5 * x_sqr + self.d * (1.0 - v + v.ln()) { + return self.d * v * self.scale + } + } + } +} + +/// The chi-squared distribution `χ²(k)`, where `k` is the degrees of +/// freedom. +/// +/// For `k > 0` integral, this distribution is the sum of the squares +/// of `k` independent standard normal random variables. For other +/// `k`, this uses the equivalent characterisation +/// `χ²(k) = Gamma(k/2, 2)`. +/// +/// # Example +/// +/// ```rust +/// use rand::distributions::{ChiSquared, IndependentSample}; +/// +/// let chi = ChiSquared::new(11.0); +/// let v = chi.ind_sample(&mut rand::thread_rng()); +/// println!("{} is from a χ²(11) distribution", v) +/// ``` +#[derive(Clone, Copy, Debug)] +pub struct ChiSquared { + repr: ChiSquaredRepr, +} + +#[derive(Clone, Copy, Debug)] +enum ChiSquaredRepr { + // k == 1, Gamma(alpha, ..) is particularly slow for alpha < 1, + // e.g. when alpha = 1/2 as it would be for this case, so special- + // casing and using the definition of N(0,1)^2 is faster. + DoFExactlyOne, + DoFAnythingElse(Gamma), +} + +impl ChiSquared { + /// Create a new chi-squared distribution with degrees-of-freedom + /// `k`. Panics if `k < 0`. + pub fn new(k: f64) -> ChiSquared { + let repr = if k == 1.0 { + DoFExactlyOne + } else { + assert!(k > 0.0, "ChiSquared::new called with `k` < 0"); + DoFAnythingElse(Gamma::new(0.5 * k, 2.0)) + }; + ChiSquared { repr: repr } + } +} +impl Sample<f64> for ChiSquared { + fn sample<R: Rng>(&mut self, rng: &mut R) -> f64 { self.ind_sample(rng) } +} +impl IndependentSample<f64> for ChiSquared { + fn ind_sample<R: Rng>(&self, rng: &mut R) -> f64 { + match self.repr { + DoFExactlyOne => { + // k == 1 => N(0,1)^2 + let StandardNormal(norm) = rng.gen::<StandardNormal>(); + norm * norm + } + DoFAnythingElse(ref g) => g.ind_sample(rng) + } + } +} + +/// The Fisher F distribution `F(m, n)`. +/// +/// This distribution is equivalent to the ratio of two normalised +/// chi-squared distributions, that is, `F(m,n) = (χ²(m)/m) / +/// (χ²(n)/n)`. +/// +/// # Example +/// +/// ```rust +/// use rand::distributions::{FisherF, IndependentSample}; +/// +/// let f = FisherF::new(2.0, 32.0); +/// let v = f.ind_sample(&mut rand::thread_rng()); +/// println!("{} is from an F(2, 32) distribution", v) +/// ``` +#[derive(Clone, Copy, Debug)] +pub struct FisherF { + numer: ChiSquared, + denom: ChiSquared, + // denom_dof / numer_dof so that this can just be a straight + // multiplication, rather than a division. + dof_ratio: f64, +} + +impl FisherF { + /// Create a new `FisherF` distribution, with the given + /// parameter. Panics if either `m` or `n` are not positive. + pub fn new(m: f64, n: f64) -> FisherF { + assert!(m > 0.0, "FisherF::new called with `m < 0`"); + assert!(n > 0.0, "FisherF::new called with `n < 0`"); + + FisherF { + numer: ChiSquared::new(m), + denom: ChiSquared::new(n), + dof_ratio: n / m + } + } +} +impl Sample<f64> for FisherF { + fn sample<R: Rng>(&mut self, rng: &mut R) -> f64 { self.ind_sample(rng) } +} +impl IndependentSample<f64> for FisherF { + fn ind_sample<R: Rng>(&self, rng: &mut R) -> f64 { + self.numer.ind_sample(rng) / self.denom.ind_sample(rng) * self.dof_ratio + } +} + +/// The Student t distribution, `t(nu)`, where `nu` is the degrees of +/// freedom. +/// +/// # Example +/// +/// ```rust +/// use rand::distributions::{StudentT, IndependentSample}; +/// +/// let t = StudentT::new(11.0); +/// let v = t.ind_sample(&mut rand::thread_rng()); +/// println!("{} is from a t(11) distribution", v) +/// ``` +#[derive(Clone, Copy, Debug)] +pub struct StudentT { + chi: ChiSquared, + dof: f64 +} + +impl StudentT { + /// Create a new Student t distribution with `n` degrees of + /// freedom. Panics if `n <= 0`. + pub fn new(n: f64) -> StudentT { + assert!(n > 0.0, "StudentT::new called with `n <= 0`"); + StudentT { + chi: ChiSquared::new(n), + dof: n + } + } +} +impl Sample<f64> for StudentT { + fn sample<R: Rng>(&mut self, rng: &mut R) -> f64 { self.ind_sample(rng) } +} +impl IndependentSample<f64> for StudentT { + fn ind_sample<R: Rng>(&self, rng: &mut R) -> f64 { + let StandardNormal(norm) = rng.gen::<StandardNormal>(); + norm * (self.dof / self.chi.ind_sample(rng)).sqrt() + } +} + +#[cfg(test)] +mod test { + use distributions::{Sample, IndependentSample}; + use super::{ChiSquared, StudentT, FisherF}; + + #[test] + fn test_chi_squared_one() { + let mut chi = ChiSquared::new(1.0); + let mut rng = ::test::rng(); + for _ in 0..1000 { + chi.sample(&mut rng); + chi.ind_sample(&mut rng); + } + } + #[test] + fn test_chi_squared_small() { + let mut chi = ChiSquared::new(0.5); + let mut rng = ::test::rng(); + for _ in 0..1000 { + chi.sample(&mut rng); + chi.ind_sample(&mut rng); + } + } + #[test] + fn test_chi_squared_large() { + let mut chi = ChiSquared::new(30.0); + let mut rng = ::test::rng(); + for _ in 0..1000 { + chi.sample(&mut rng); + chi.ind_sample(&mut rng); + } + } + #[test] + #[should_panic] + fn test_chi_squared_invalid_dof() { + ChiSquared::new(-1.0); + } + + #[test] + fn test_f() { + let mut f = FisherF::new(2.0, 32.0); + let mut rng = ::test::rng(); + for _ in 0..1000 { + f.sample(&mut rng); + f.ind_sample(&mut rng); + } + } + + #[test] + fn test_t() { + let mut t = StudentT::new(11.0); + let mut rng = ::test::rng(); + for _ in 0..1000 { + t.sample(&mut rng); + t.ind_sample(&mut rng); + } + } +} diff --git a/rand/src/distributions/mod.rs b/rand/src/distributions/mod.rs new file mode 100644 index 0000000..5de8efb --- /dev/null +++ b/rand/src/distributions/mod.rs @@ -0,0 +1,409 @@ +// Copyright 2013 The Rust Project Developers. See the COPYRIGHT +// file at the top-level directory of this distribution and at +// http://rust-lang.org/COPYRIGHT. +// +// Licensed under the Apache License, Version 2.0 <LICENSE-APACHE or +// http://www.apache.org/licenses/LICENSE-2.0> or the MIT license +// <LICENSE-MIT or http://opensource.org/licenses/MIT>, at your +// option. This file may not be copied, modified, or distributed +// except according to those terms. + +//! Sampling from random distributions. +//! +//! This is a generalization of `Rand` to allow parameters to control the +//! exact properties of the generated values, e.g. the mean and standard +//! deviation of a normal distribution. The `Sample` trait is the most +//! general, and allows for generating values that change some state +//! internally. The `IndependentSample` trait is for generating values +//! that do not need to record state. + +use core::marker; + +use {Rng, Rand}; + +pub use self::range::Range; +#[cfg(feature="std")] +pub use self::gamma::{Gamma, ChiSquared, FisherF, StudentT}; +#[cfg(feature="std")] +pub use self::normal::{Normal, LogNormal}; +#[cfg(feature="std")] +pub use self::exponential::Exp; + +pub mod range; +#[cfg(feature="std")] +pub mod gamma; +#[cfg(feature="std")] +pub mod normal; +#[cfg(feature="std")] +pub mod exponential; + +#[cfg(feature="std")] +mod ziggurat_tables; + +/// Types that can be used to create a random instance of `Support`. +pub trait Sample<Support> { + /// Generate a random value of `Support`, using `rng` as the + /// source of randomness. + fn sample<R: Rng>(&mut self, rng: &mut R) -> Support; +} + +/// `Sample`s that do not require keeping track of state. +/// +/// Since no state is recorded, each sample is (statistically) +/// independent of all others, assuming the `Rng` used has this +/// property. +// FIXME maybe having this separate is overkill (the only reason is to +// take &self rather than &mut self)? or maybe this should be the +// trait called `Sample` and the other should be `DependentSample`. +pub trait IndependentSample<Support>: Sample<Support> { + /// Generate a random value. + fn ind_sample<R: Rng>(&self, &mut R) -> Support; +} + +/// A wrapper for generating types that implement `Rand` via the +/// `Sample` & `IndependentSample` traits. +#[derive(Debug)] +pub struct RandSample<Sup> { + _marker: marker::PhantomData<fn() -> Sup>, +} + +impl<Sup> Copy for RandSample<Sup> {} +impl<Sup> Clone for RandSample<Sup> { + fn clone(&self) -> Self { *self } +} + +impl<Sup: Rand> Sample<Sup> for RandSample<Sup> { + fn sample<R: Rng>(&mut self, rng: &mut R) -> Sup { self.ind_sample(rng) } +} + +impl<Sup: Rand> IndependentSample<Sup> for RandSample<Sup> { + fn ind_sample<R: Rng>(&self, rng: &mut R) -> Sup { + rng.gen() + } +} + +impl<Sup> RandSample<Sup> { + pub fn new() -> RandSample<Sup> { + RandSample { _marker: marker::PhantomData } + } +} + +/// A value with a particular weight for use with `WeightedChoice`. +#[derive(Copy, Clone, Debug)] +pub struct Weighted<T> { + /// The numerical weight of this item + pub weight: u32, + /// The actual item which is being weighted + pub item: T, +} + +/// A distribution that selects from a finite collection of weighted items. +/// +/// Each item has an associated weight that influences how likely it +/// is to be chosen: higher weight is more likely. +/// +/// The `Clone` restriction is a limitation of the `Sample` and +/// `IndependentSample` traits. Note that `&T` is (cheaply) `Clone` for +/// all `T`, as is `u32`, so one can store references or indices into +/// another vector. +/// +/// # Example +/// +/// ```rust +/// use rand::distributions::{Weighted, WeightedChoice, IndependentSample}; +/// +/// let mut items = vec!(Weighted { weight: 2, item: 'a' }, +/// Weighted { weight: 4, item: 'b' }, +/// Weighted { weight: 1, item: 'c' }); +/// let wc = WeightedChoice::new(&mut items); +/// let mut rng = rand::thread_rng(); +/// for _ in 0..16 { +/// // on average prints 'a' 4 times, 'b' 8 and 'c' twice. +/// println!("{}", wc.ind_sample(&mut rng)); +/// } +/// ``` +#[derive(Debug)] +pub struct WeightedChoice<'a, T:'a> { + items: &'a mut [Weighted<T>], + weight_range: Range<u32> +} + +impl<'a, T: Clone> WeightedChoice<'a, T> { + /// Create a new `WeightedChoice`. + /// + /// Panics if: + /// + /// - `items` is empty + /// - the total weight is 0 + /// - the total weight is larger than a `u32` can contain. + pub fn new(items: &'a mut [Weighted<T>]) -> WeightedChoice<'a, T> { + // strictly speaking, this is subsumed by the total weight == 0 case + assert!(!items.is_empty(), "WeightedChoice::new called with no items"); + + let mut running_total: u32 = 0; + + // we convert the list from individual weights to cumulative + // weights so we can binary search. This *could* drop elements + // with weight == 0 as an optimisation. + for item in items.iter_mut() { + running_total = match running_total.checked_add(item.weight) { + Some(n) => n, + None => panic!("WeightedChoice::new called with a total weight \ + larger than a u32 can contain") + }; + + item.weight = running_total; + } + assert!(running_total != 0, "WeightedChoice::new called with a total weight of 0"); + + WeightedChoice { + items: items, + // we're likely to be generating numbers in this range + // relatively often, so might as well cache it + weight_range: Range::new(0, running_total) + } + } +} + +impl<'a, T: Clone> Sample<T> for WeightedChoice<'a, T> { + fn sample<R: Rng>(&mut self, rng: &mut R) -> T { self.ind_sample(rng) } +} + +impl<'a, T: Clone> IndependentSample<T> for WeightedChoice<'a, T> { + fn ind_sample<R: Rng>(&self, rng: &mut R) -> T { + // we want to find the first element that has cumulative + // weight > sample_weight, which we do by binary since the + // cumulative weights of self.items are sorted. + + // choose a weight in [0, total_weight) + let sample_weight = self.weight_range.ind_sample(rng); + + // short circuit when it's the first item + if sample_weight < self.items[0].weight { + return self.items[0].item.clone(); + } + + let mut idx = 0; + let mut modifier = self.items.len(); + + // now we know that every possibility has an element to the + // left, so we can just search for the last element that has + // cumulative weight <= sample_weight, then the next one will + // be "it". (Note that this greatest element will never be the + // last element of the vector, since sample_weight is chosen + // in [0, total_weight) and the cumulative weight of the last + // one is exactly the total weight.) + while modifier > 1 { + let i = idx + modifier / 2; + if self.items[i].weight <= sample_weight { + // we're small, so look to the right, but allow this + // exact element still. + idx = i; + // we need the `/ 2` to round up otherwise we'll drop + // the trailing elements when `modifier` is odd. + modifier += 1; + } else { + // otherwise we're too big, so go left. (i.e. do + // nothing) + } + modifier /= 2; + } + return self.items[idx + 1].item.clone(); + } +} + +/// Sample a random number using the Ziggurat method (specifically the +/// ZIGNOR variant from Doornik 2005). Most of the arguments are +/// directly from the paper: +/// +/// * `rng`: source of randomness +/// * `symmetric`: whether this is a symmetric distribution, or one-sided with P(x < 0) = 0. +/// * `X`: the $x_i$ abscissae. +/// * `F`: precomputed values of the PDF at the $x_i$, (i.e. $f(x_i)$) +/// * `F_DIFF`: precomputed values of $f(x_i) - f(x_{i+1})$ +/// * `pdf`: the probability density function +/// * `zero_case`: manual sampling from the tail when we chose the +/// bottom box (i.e. i == 0) + +// the perf improvement (25-50%) is definitely worth the extra code +// size from force-inlining. +#[cfg(feature="std")] +#[inline(always)] +fn ziggurat<R: Rng, P, Z>( + rng: &mut R, + symmetric: bool, + x_tab: ziggurat_tables::ZigTable, + f_tab: ziggurat_tables::ZigTable, + mut pdf: P, + mut zero_case: Z) + -> f64 where P: FnMut(f64) -> f64, Z: FnMut(&mut R, f64) -> f64 { + const SCALE: f64 = (1u64 << 53) as f64; + loop { + // reimplement the f64 generation as an optimisation suggested + // by the Doornik paper: we have a lot of precision-space + // (i.e. there are 11 bits of the 64 of a u64 to use after + // creating a f64), so we might as well reuse some to save + // generating a whole extra random number. (Seems to be 15% + // faster.) + // + // This unfortunately misses out on the benefits of direct + // floating point generation if an RNG like dSMFT is + // used. (That is, such RNGs create floats directly, highly + // efficiently and overload next_f32/f64, so by not calling it + // this may be slower than it would be otherwise.) + // FIXME: investigate/optimise for the above. + let bits: u64 = rng.gen(); + let i = (bits & 0xff) as usize; + let f = (bits >> 11) as f64 / SCALE; + + // u is either U(-1, 1) or U(0, 1) depending on if this is a + // symmetric distribution or not. + let u = if symmetric {2.0 * f - 1.0} else {f}; + let x = u * x_tab[i]; + + let test_x = if symmetric { x.abs() } else {x}; + + // algebraically equivalent to |u| < x_tab[i+1]/x_tab[i] (or u < x_tab[i+1]/x_tab[i]) + if test_x < x_tab[i + 1] { + return x; + } + if i == 0 { + return zero_case(rng, u); + } + // algebraically equivalent to f1 + DRanU()*(f0 - f1) < 1 + if f_tab[i + 1] + (f_tab[i] - f_tab[i + 1]) * rng.gen::<f64>() < pdf(x) { + return x; + } + } +} + +#[cfg(test)] +mod tests { + + use {Rng, Rand}; + use super::{RandSample, WeightedChoice, Weighted, Sample, IndependentSample}; + + #[derive(PartialEq, Debug)] + struct ConstRand(usize); + impl Rand for ConstRand { + fn rand<R: Rng>(_: &mut R) -> ConstRand { + ConstRand(0) + } + } + + // 0, 1, 2, 3, ... + struct CountingRng { i: u32 } + impl Rng for CountingRng { + fn next_u32(&mut self) -> u32 { + self.i += 1; + self.i - 1 + } + fn next_u64(&mut self) -> u64 { + self.next_u32() as u64 + } + } + + #[test] + fn test_rand_sample() { + let mut rand_sample = RandSample::<ConstRand>::new(); + + assert_eq!(rand_sample.sample(&mut ::test::rng()), ConstRand(0)); + assert_eq!(rand_sample.ind_sample(&mut ::test::rng()), ConstRand(0)); + } + #[test] + fn test_weighted_choice() { + // this makes assumptions about the internal implementation of + // WeightedChoice, specifically: it doesn't reorder the items, + // it doesn't do weird things to the RNG (so 0 maps to 0, 1 to + // 1, internally; modulo a modulo operation). + + macro_rules! t { + ($items:expr, $expected:expr) => {{ + let mut items = $items; + let wc = WeightedChoice::new(&mut items); + let expected = $expected; + + let mut rng = CountingRng { i: 0 }; + + for &val in expected.iter() { + assert_eq!(wc.ind_sample(&mut rng), val) + } + }} + } + + t!(vec!(Weighted { weight: 1, item: 10}), [10]); + + // skip some + t!(vec!(Weighted { weight: 0, item: 20}, + Weighted { weight: 2, item: 21}, + Weighted { weight: 0, item: 22}, + Weighted { weight: 1, item: 23}), + [21,21, 23]); + + // different weights + t!(vec!(Weighted { weight: 4, item: 30}, + Weighted { weight: 3, item: 31}), + [30,30,30,30, 31,31,31]); + + // check that we're binary searching + // correctly with some vectors of odd + // length. + t!(vec!(Weighted { weight: 1, item: 40}, + Weighted { weight: 1, item: 41}, + Weighted { weight: 1, item: 42}, + Weighted { weight: 1, item: 43}, + Weighted { weight: 1, item: 44}), + [40, 41, 42, 43, 44]); + t!(vec!(Weighted { weight: 1, item: 50}, + Weighted { weight: 1, item: 51}, + Weighted { weight: 1, item: 52}, + Weighted { weight: 1, item: 53}, + Weighted { weight: 1, item: 54}, + Weighted { weight: 1, item: 55}, + Weighted { weight: 1, item: 56}), + [50, 51, 52, 53, 54, 55, 56]); + } + + #[test] + fn test_weighted_clone_initialization() { + let initial : Weighted<u32> = Weighted {weight: 1, item: 1}; + let clone = initial.clone(); + assert_eq!(initial.weight, clone.weight); + assert_eq!(initial.item, clone.item); + } + + #[test] #[should_panic] + fn test_weighted_clone_change_weight() { + let initial : Weighted<u32> = Weighted {weight: 1, item: 1}; + let mut clone = initial.clone(); + clone.weight = 5; + assert_eq!(initial.weight, clone.weight); + } + + #[test] #[should_panic] + fn test_weighted_clone_change_item() { + let initial : Weighted<u32> = Weighted {weight: 1, item: 1}; + let mut clone = initial.clone(); + clone.item = 5; + assert_eq!(initial.item, clone.item); + + } + + #[test] #[should_panic] + fn test_weighted_choice_no_items() { + WeightedChoice::<isize>::new(&mut []); + } + #[test] #[should_panic] + fn test_weighted_choice_zero_weight() { + WeightedChoice::new(&mut [Weighted { weight: 0, item: 0}, + Weighted { weight: 0, item: 1}]); + } + #[test] #[should_panic] + fn test_weighted_choice_weight_overflows() { + let x = ::std::u32::MAX / 2; // x + x + 2 is the overflow + WeightedChoice::new(&mut [Weighted { weight: x, item: 0 }, + Weighted { weight: 1, item: 1 }, + Weighted { weight: x, item: 2 }, + Weighted { weight: 1, item: 3 }]); + } +} diff --git a/rand/src/distributions/normal.rs b/rand/src/distributions/normal.rs new file mode 100644 index 0000000..280613d --- /dev/null +++ b/rand/src/distributions/normal.rs @@ -0,0 +1,201 @@ +// Copyright 2013 The Rust Project Developers. See the COPYRIGHT +// file at the top-level directory of this distribution and at +// http://rust-lang.org/COPYRIGHT. +// +// Licensed under the Apache License, Version 2.0 <LICENSE-APACHE or +// http://www.apache.org/licenses/LICENSE-2.0> or the MIT license +// <LICENSE-MIT or http://opensource.org/licenses/MIT>, at your +// option. This file may not be copied, modified, or distributed +// except according to those terms. + +//! The normal and derived distributions. + +use {Rng, Rand, Open01}; +use distributions::{ziggurat, ziggurat_tables, Sample, IndependentSample}; + +/// A wrapper around an `f64` to generate N(0, 1) random numbers +/// (a.k.a. a standard normal, or Gaussian). +/// +/// See `Normal` for the general normal distribution. +/// +/// Implemented via the ZIGNOR variant[1] of the Ziggurat method. +/// +/// [1]: Jurgen A. Doornik (2005). [*An Improved Ziggurat Method to +/// Generate Normal Random +/// Samples*](http://www.doornik.com/research/ziggurat.pdf). Nuffield +/// College, Oxford +/// +/// # Example +/// +/// ```rust +/// use rand::distributions::normal::StandardNormal; +/// +/// let StandardNormal(x) = rand::random(); +/// println!("{}", x); +/// ``` +#[derive(Clone, Copy, Debug)] +pub struct StandardNormal(pub f64); + +impl Rand for StandardNormal { + fn rand<R:Rng>(rng: &mut R) -> StandardNormal { + #[inline] + fn pdf(x: f64) -> f64 { + (-x*x/2.0).exp() + } + #[inline] + fn zero_case<R:Rng>(rng: &mut R, u: f64) -> f64 { + // compute a random number in the tail by hand + + // strange initial conditions, because the loop is not + // do-while, so the condition should be true on the first + // run, they get overwritten anyway (0 < 1, so these are + // good). + let mut x = 1.0f64; + let mut y = 0.0f64; + + while -2.0 * y < x * x { + let Open01(x_) = rng.gen::<Open01<f64>>(); + let Open01(y_) = rng.gen::<Open01<f64>>(); + + x = x_.ln() / ziggurat_tables::ZIG_NORM_R; + y = y_.ln(); + } + + if u < 0.0 { x - ziggurat_tables::ZIG_NORM_R } else { ziggurat_tables::ZIG_NORM_R - x } + } + + StandardNormal(ziggurat( + rng, + true, // this is symmetric + &ziggurat_tables::ZIG_NORM_X, + &ziggurat_tables::ZIG_NORM_F, + pdf, zero_case)) + } +} + +/// The normal distribution `N(mean, std_dev**2)`. +/// +/// This uses the ZIGNOR variant of the Ziggurat method, see +/// `StandardNormal` for more details. +/// +/// # Example +/// +/// ```rust +/// use rand::distributions::{Normal, IndependentSample}; +/// +/// // mean 2, standard deviation 3 +/// let normal = Normal::new(2.0, 3.0); +/// let v = normal.ind_sample(&mut rand::thread_rng()); +/// println!("{} is from a N(2, 9) distribution", v) +/// ``` +#[derive(Clone, Copy, Debug)] +pub struct Normal { + mean: f64, + std_dev: f64, +} + +impl Normal { + /// Construct a new `Normal` distribution with the given mean and + /// standard deviation. + /// + /// # Panics + /// + /// Panics if `std_dev < 0`. + #[inline] + pub fn new(mean: f64, std_dev: f64) -> Normal { + assert!(std_dev >= 0.0, "Normal::new called with `std_dev` < 0"); + Normal { + mean: mean, + std_dev: std_dev + } + } +} +impl Sample<f64> for Normal { + fn sample<R: Rng>(&mut self, rng: &mut R) -> f64 { self.ind_sample(rng) } +} +impl IndependentSample<f64> for Normal { + fn ind_sample<R: Rng>(&self, rng: &mut R) -> f64 { + let StandardNormal(n) = rng.gen::<StandardNormal>(); + self.mean + self.std_dev * n + } +} + + +/// The log-normal distribution `ln N(mean, std_dev**2)`. +/// +/// If `X` is log-normal distributed, then `ln(X)` is `N(mean, +/// std_dev**2)` distributed. +/// +/// # Example +/// +/// ```rust +/// use rand::distributions::{LogNormal, IndependentSample}; +/// +/// // mean 2, standard deviation 3 +/// let log_normal = LogNormal::new(2.0, 3.0); +/// let v = log_normal.ind_sample(&mut rand::thread_rng()); +/// println!("{} is from an ln N(2, 9) distribution", v) +/// ``` +#[derive(Clone, Copy, Debug)] +pub struct LogNormal { + norm: Normal +} + +impl LogNormal { + /// Construct a new `LogNormal` distribution with the given mean + /// and standard deviation. + /// + /// # Panics + /// + /// Panics if `std_dev < 0`. + #[inline] + pub fn new(mean: f64, std_dev: f64) -> LogNormal { + assert!(std_dev >= 0.0, "LogNormal::new called with `std_dev` < 0"); + LogNormal { norm: Normal::new(mean, std_dev) } + } +} +impl Sample<f64> for LogNormal { + fn sample<R: Rng>(&mut self, rng: &mut R) -> f64 { self.ind_sample(rng) } +} +impl IndependentSample<f64> for LogNormal { + fn ind_sample<R: Rng>(&self, rng: &mut R) -> f64 { + self.norm.ind_sample(rng).exp() + } +} + +#[cfg(test)] +mod tests { + use distributions::{Sample, IndependentSample}; + use super::{Normal, LogNormal}; + + #[test] + fn test_normal() { + let mut norm = Normal::new(10.0, 10.0); + let mut rng = ::test::rng(); + for _ in 0..1000 { + norm.sample(&mut rng); + norm.ind_sample(&mut rng); + } + } + #[test] + #[should_panic] + fn test_normal_invalid_sd() { + Normal::new(10.0, -1.0); + } + + + #[test] + fn test_log_normal() { + let mut lnorm = LogNormal::new(10.0, 10.0); + let mut rng = ::test::rng(); + for _ in 0..1000 { + lnorm.sample(&mut rng); + lnorm.ind_sample(&mut rng); + } + } + #[test] + #[should_panic] + fn test_log_normal_invalid_sd() { + LogNormal::new(10.0, -1.0); + } +} diff --git a/rand/src/distributions/range.rs b/rand/src/distributions/range.rs new file mode 100644 index 0000000..935a00a --- /dev/null +++ b/rand/src/distributions/range.rs @@ -0,0 +1,241 @@ +// Copyright 2013 The Rust Project Developers. See the COPYRIGHT +// file at the top-level directory of this distribution and at +// http://rust-lang.org/COPYRIGHT. +// +// Licensed under the Apache License, Version 2.0 <LICENSE-APACHE or +// http://www.apache.org/licenses/LICENSE-2.0> or the MIT license +// <LICENSE-MIT or http://opensource.org/licenses/MIT>, at your +// option. This file may not be copied, modified, or distributed +// except according to those terms. + +//! Generating numbers between two others. + +// this is surprisingly complicated to be both generic & correct + +use core::num::Wrapping as w; + +use Rng; +use distributions::{Sample, IndependentSample}; + +/// Sample values uniformly between two bounds. +/// +/// This gives a uniform distribution (assuming the RNG used to sample +/// it is itself uniform & the `SampleRange` implementation for the +/// given type is correct), even for edge cases like `low = 0u8`, +/// `high = 170u8`, for which a naive modulo operation would return +/// numbers less than 85 with double the probability to those greater +/// than 85. +/// +/// Types should attempt to sample in `[low, high)`, i.e., not +/// including `high`, but this may be very difficult. All the +/// primitive integer types satisfy this property, and the float types +/// normally satisfy it, but rounding may mean `high` can occur. +/// +/// # Example +/// +/// ```rust +/// use rand::distributions::{IndependentSample, Range}; +/// +/// fn main() { +/// let between = Range::new(10, 10000); +/// let mut rng = rand::thread_rng(); +/// let mut sum = 0; +/// for _ in 0..1000 { +/// sum += between.ind_sample(&mut rng); +/// } +/// println!("{}", sum); +/// } +/// ``` +#[derive(Clone, Copy, Debug)] +pub struct Range<X> { + low: X, + range: X, + accept_zone: X +} + +impl<X: SampleRange + PartialOrd> Range<X> { + /// Create a new `Range` instance that samples uniformly from + /// `[low, high)`. Panics if `low >= high`. + pub fn new(low: X, high: X) -> Range<X> { + assert!(low < high, "Range::new called with `low >= high`"); + SampleRange::construct_range(low, high) + } +} + +impl<Sup: SampleRange> Sample<Sup> for Range<Sup> { + #[inline] + fn sample<R: Rng>(&mut self, rng: &mut R) -> Sup { self.ind_sample(rng) } +} +impl<Sup: SampleRange> IndependentSample<Sup> for Range<Sup> { + fn ind_sample<R: Rng>(&self, rng: &mut R) -> Sup { + SampleRange::sample_range(self, rng) + } +} + +/// The helper trait for types that have a sensible way to sample +/// uniformly between two values. This should not be used directly, +/// and is only to facilitate `Range`. +pub trait SampleRange : Sized { + /// Construct the `Range` object that `sample_range` + /// requires. This should not ever be called directly, only via + /// `Range::new`, which will check that `low < high`, so this + /// function doesn't have to repeat the check. + fn construct_range(low: Self, high: Self) -> Range<Self>; + + /// Sample a value from the given `Range` with the given `Rng` as + /// a source of randomness. + fn sample_range<R: Rng>(r: &Range<Self>, rng: &mut R) -> Self; +} + +macro_rules! integer_impl { + ($ty:ty, $unsigned:ident) => { + impl SampleRange for $ty { + // we play free and fast with unsigned vs signed here + // (when $ty is signed), but that's fine, since the + // contract of this macro is for $ty and $unsigned to be + // "bit-equal", so casting between them is a no-op & a + // bijection. + + #[inline] + fn construct_range(low: $ty, high: $ty) -> Range<$ty> { + let range = (w(high as $unsigned) - w(low as $unsigned)).0; + let unsigned_max: $unsigned = ::core::$unsigned::MAX; + + // this is the largest number that fits into $unsigned + // that `range` divides evenly, so, if we've sampled + // `n` uniformly from this region, then `n % range` is + // uniform in [0, range) + let zone = unsigned_max - unsigned_max % range; + + Range { + low: low, + range: range as $ty, + accept_zone: zone as $ty + } + } + + #[inline] + fn sample_range<R: Rng>(r: &Range<$ty>, rng: &mut R) -> $ty { + loop { + // rejection sample + let v = rng.gen::<$unsigned>(); + // until we find something that fits into the + // region which r.range evenly divides (this will + // be uniformly distributed) + if v < r.accept_zone as $unsigned { + // and return it, with some adjustments + return (w(r.low) + w((v % r.range as $unsigned) as $ty)).0; + } + } + } + } + } +} + +integer_impl! { i8, u8 } +integer_impl! { i16, u16 } +integer_impl! { i32, u32 } +integer_impl! { i64, u64 } +#[cfg(feature = "i128_support")] +integer_impl! { i128, u128 } +integer_impl! { isize, usize } +integer_impl! { u8, u8 } +integer_impl! { u16, u16 } +integer_impl! { u32, u32 } +integer_impl! { u64, u64 } +#[cfg(feature = "i128_support")] +integer_impl! { u128, u128 } +integer_impl! { usize, usize } + +macro_rules! float_impl { + ($ty:ty) => { + impl SampleRange for $ty { + fn construct_range(low: $ty, high: $ty) -> Range<$ty> { + Range { + low: low, + range: high - low, + accept_zone: 0.0 // unused + } + } + fn sample_range<R: Rng>(r: &Range<$ty>, rng: &mut R) -> $ty { + r.low + r.range * rng.gen::<$ty>() + } + } + } +} + +float_impl! { f32 } +float_impl! { f64 } + +#[cfg(test)] +mod tests { + use distributions::{Sample, IndependentSample}; + use super::Range as Range; + + #[should_panic] + #[test] + fn test_range_bad_limits_equal() { + Range::new(10, 10); + } + #[should_panic] + #[test] + fn test_range_bad_limits_flipped() { + Range::new(10, 5); + } + + #[test] + fn test_integers() { + let mut rng = ::test::rng(); + macro_rules! t { + ($($ty:ident),*) => {{ + $( + let v: &[($ty, $ty)] = &[(0, 10), + (10, 127), + (::core::$ty::MIN, ::core::$ty::MAX)]; + for &(low, high) in v.iter() { + let mut sampler: Range<$ty> = Range::new(low, high); + for _ in 0..1000 { + let v = sampler.sample(&mut rng); + assert!(low <= v && v < high); + let v = sampler.ind_sample(&mut rng); + assert!(low <= v && v < high); + } + } + )* + }} + } + #[cfg(not(feature = "i128_support"))] + t!(i8, i16, i32, i64, isize, + u8, u16, u32, u64, usize); + #[cfg(feature = "i128_support")] + t!(i8, i16, i32, i64, i128, isize, + u8, u16, u32, u64, u128, usize); + } + + #[test] + fn test_floats() { + let mut rng = ::test::rng(); + macro_rules! t { + ($($ty:ty),*) => {{ + $( + let v: &[($ty, $ty)] = &[(0.0, 100.0), + (-1e35, -1e25), + (1e-35, 1e-25), + (-1e35, 1e35)]; + for &(low, high) in v.iter() { + let mut sampler: Range<$ty> = Range::new(low, high); + for _ in 0..1000 { + let v = sampler.sample(&mut rng); + assert!(low <= v && v < high); + let v = sampler.ind_sample(&mut rng); + assert!(low <= v && v < high); + } + } + )* + }} + } + + t!(f32, f64) + } + +} diff --git a/rand/src/distributions/ziggurat_tables.rs b/rand/src/distributions/ziggurat_tables.rs new file mode 100644 index 0000000..b6de4bf --- /dev/null +++ b/rand/src/distributions/ziggurat_tables.rs @@ -0,0 +1,280 @@ +// Copyright 2013 The Rust Project Developers. See the COPYRIGHT +// file at the top-level directory of this distribution and at +// http://rust-lang.org/COPYRIGHT. +// +// Licensed under the Apache License, Version 2.0 <LICENSE-APACHE or +// http://www.apache.org/licenses/LICENSE-2.0> or the MIT license +// <LICENSE-MIT or http://opensource.org/licenses/MIT>, at your +// option. This file may not be copied, modified, or distributed +// except according to those terms. + +// Tables for distributions which are sampled using the ziggurat +// algorithm. Autogenerated by `ziggurat_tables.py`. + +pub type ZigTable = &'static [f64; 257]; +pub const ZIG_NORM_R: f64 = 3.654152885361008796; +pub static ZIG_NORM_X: [f64; 257] = + [3.910757959537090045, 3.654152885361008796, 3.449278298560964462, 3.320244733839166074, + 3.224575052047029100, 3.147889289517149969, 3.083526132001233044, 3.027837791768635434, + 2.978603279880844834, 2.934366867207854224, 2.894121053612348060, 2.857138730872132548, + 2.822877396825325125, 2.790921174000785765, 2.760944005278822555, 2.732685359042827056, + 2.705933656121858100, 2.680514643284522158, 2.656283037575502437, 2.633116393630324570, + 2.610910518487548515, 2.589575986706995181, 2.569035452680536569, 2.549221550323460761, + 2.530075232158516929, 2.511544441625342294, 2.493583041269680667, 2.476149939669143318, + 2.459208374333311298, 2.442725318198956774, 2.426670984935725972, 2.411018413899685520, + 2.395743119780480601, 2.380822795170626005, 2.366237056715818632, 2.351967227377659952, + 2.337996148795031370, 2.324308018869623016, 2.310888250599850036, 2.297723348901329565, + 2.284800802722946056, 2.272108990226823888, 2.259637095172217780, 2.247375032945807760, + 2.235313384928327984, 2.223443340090905718, 2.211756642882544366, 2.200245546609647995, + 2.188902771624720689, 2.177721467738641614, 2.166695180352645966, 2.155817819875063268, + 2.145083634046203613, 2.134487182844320152, 2.124023315687815661, 2.113687150684933957, + 2.103474055713146829, 2.093379631137050279, 2.083399693996551783, 2.073530263516978778, + 2.063767547809956415, 2.054107931648864849, 2.044547965215732788, 2.035084353727808715, + 2.025713947862032960, 2.016433734904371722, 2.007240830558684852, 1.998132471356564244, + 1.989106007615571325, 1.980158896898598364, 1.971288697931769640, 1.962493064942461896, + 1.953769742382734043, 1.945116560006753925, 1.936531428273758904, 1.928012334050718257, + 1.919557336591228847, 1.911164563769282232, 1.902832208548446369, 1.894558525668710081, + 1.886341828534776388, 1.878180486290977669, 1.870072921069236838, 1.862017605397632281, + 1.854013059758148119, 1.846057850283119750, 1.838150586580728607, 1.830289919680666566, + 1.822474540091783224, 1.814703175964167636, 1.806974591348693426, 1.799287584547580199, + 1.791640986550010028, 1.784033659547276329, 1.776464495522344977, 1.768932414909077933, + 1.761436365316706665, 1.753975320315455111, 1.746548278279492994, 1.739154261283669012, + 1.731792314050707216, 1.724461502945775715, 1.717160915015540690, 1.709889657069006086, + 1.702646854797613907, 1.695431651932238548, 1.688243209434858727, 1.681080704722823338, + 1.673943330923760353, 1.666830296159286684, 1.659740822855789499, 1.652674147080648526, + 1.645629517902360339, 1.638606196773111146, 1.631603456932422036, 1.624620582830568427, + 1.617656869570534228, 1.610711622367333673, 1.603784156023583041, 1.596873794420261339, + 1.589979870021648534, 1.583101723393471438, 1.576238702733332886, 1.569390163412534456, + 1.562555467528439657, 1.555733983466554893, 1.548925085471535512, 1.542128153226347553, + 1.535342571438843118, 1.528567729435024614, 1.521803020758293101, 1.515047842773992404, + 1.508301596278571965, 1.501563685112706548, 1.494833515777718391, 1.488110497054654369, + 1.481394039625375747, 1.474683555695025516, 1.467978458615230908, 1.461278162507407830, + 1.454582081885523293, 1.447889631277669675, 1.441200224845798017, 1.434513276002946425, + 1.427828197027290358, 1.421144398672323117, 1.414461289772464658, 1.407778276843371534, + 1.401094763676202559, 1.394410150925071257, 1.387723835686884621, 1.381035211072741964, + 1.374343665770030531, 1.367648583594317957, 1.360949343030101844, 1.354245316759430606, + 1.347535871177359290, 1.340820365893152122, 1.334098153216083604, 1.327368577624624679, + 1.320630975217730096, 1.313884673146868964, 1.307128989027353860, 1.300363230327433728, + 1.293586693733517645, 1.286798664489786415, 1.279998415710333237, 1.273185207661843732, + 1.266358287014688333, 1.259516886060144225, 1.252660221891297887, 1.245787495544997903, + 1.238897891102027415, 1.231990574742445110, 1.225064693752808020, 1.218119375481726552, + 1.211153726239911244, 1.204166830140560140, 1.197157747875585931, 1.190125515422801650, + 1.183069142678760732, 1.175987612011489825, 1.168879876726833800, 1.161744859441574240, + 1.154581450355851802, 1.147388505416733873, 1.140164844363995789, 1.132909248648336975, + 1.125620459211294389, 1.118297174115062909, 1.110938046009249502, 1.103541679420268151, + 1.096106627847603487, 1.088631390649514197, 1.081114409698889389, 1.073554065787871714, + 1.065948674757506653, 1.058296483326006454, 1.050595664586207123, 1.042844313139370538, + 1.035040439828605274, 1.027181966030751292, 1.019266717460529215, 1.011292417434978441, + 1.003256679539591412, 0.995156999629943084, 0.986990747093846266, 0.978755155288937750, + 0.970447311058864615, 0.962064143217605250, 0.953602409875572654, 0.945058684462571130, + 0.936429340280896860, 0.927710533396234771, 0.918898183643734989, 0.909987953490768997, + 0.900975224455174528, 0.891855070726792376, 0.882622229578910122, 0.873271068082494550, + 0.863795545546826915, 0.854189171001560554, 0.844444954902423661, 0.834555354079518752, + 0.824512208745288633, 0.814306670128064347, 0.803929116982664893, 0.793369058833152785, + 0.782615023299588763, 0.771654424216739354, 0.760473406422083165, 0.749056662009581653, + 0.737387211425838629, 0.725446140901303549, 0.713212285182022732, 0.700661841097584448, + 0.687767892786257717, 0.674499822827436479, 0.660822574234205984, 0.646695714884388928, + 0.632072236375024632, 0.616896989996235545, 0.601104617743940417, 0.584616766093722262, + 0.567338257040473026, 0.549151702313026790, 0.529909720646495108, 0.509423329585933393, + 0.487443966121754335, 0.463634336771763245, 0.437518402186662658, 0.408389134588000746, + 0.375121332850465727, 0.335737519180459465, 0.286174591747260509, 0.215241895913273806, + 0.000000000000000000]; +pub static ZIG_NORM_F: [f64; 257] = + [0.000477467764586655, 0.001260285930498598, 0.002609072746106363, 0.004037972593371872, + 0.005522403299264754, 0.007050875471392110, 0.008616582769422917, 0.010214971439731100, + 0.011842757857943104, 0.013497450601780807, 0.015177088307982072, 0.016880083152595839, + 0.018605121275783350, 0.020351096230109354, 0.022117062707379922, 0.023902203305873237, + 0.025705804008632656, 0.027527235669693315, 0.029365939758230111, 0.031221417192023690, + 0.033093219458688698, 0.034980941461833073, 0.036884215688691151, 0.038802707404656918, + 0.040736110656078753, 0.042684144916619378, 0.044646552251446536, 0.046623094902089664, + 0.048613553216035145, 0.050617723861121788, 0.052635418276973649, 0.054666461325077916, + 0.056710690106399467, 0.058767952921137984, 0.060838108349751806, 0.062921024437977854, + 0.065016577971470438, 0.067124653828023989, 0.069245144397250269, 0.071377949059141965, + 0.073522973714240991, 0.075680130359194964, 0.077849336702372207, 0.080030515814947509, + 0.082223595813495684, 0.084428509570654661, 0.086645194450867782, 0.088873592068594229, + 0.091113648066700734, 0.093365311913026619, 0.095628536713353335, 0.097903279039215627, + 0.100189498769172020, 0.102487158942306270, 0.104796225622867056, 0.107116667775072880, + 0.109448457147210021, 0.111791568164245583, 0.114145977828255210, 0.116511665626037014, + 0.118888613443345698, 0.121276805485235437, 0.123676228202051403, 0.126086870220650349, + 0.128508722280473636, 0.130941777174128166, 0.133386029692162844, 0.135841476571757352, + 0.138308116449064322, 0.140785949814968309, 0.143274978974047118, 0.145775208006537926, + 0.148286642733128721, 0.150809290682410169, 0.153343161060837674, 0.155888264725064563, + 0.158444614156520225, 0.161012223438117663, 0.163591108232982951, 0.166181285765110071, + 0.168782774801850333, 0.171395595638155623, 0.174019770082499359, 0.176655321444406654, + 0.179302274523530397, 0.181960655600216487, 0.184630492427504539, 0.187311814224516926, + 0.190004651671193070, 0.192709036904328807, 0.195425003514885592, 0.198152586546538112, + 0.200891822495431333, 0.203642749311121501, 0.206405406398679298, 0.209179834621935651, + 0.211966076307852941, 0.214764175252008499, 0.217574176725178370, 0.220396127481011589, + 0.223230075764789593, 0.226076071323264877, 0.228934165415577484, 0.231804410825248525, + 0.234686861873252689, 0.237581574432173676, 0.240488605941449107, 0.243408015423711988, + 0.246339863502238771, 0.249284212419516704, 0.252241126056943765, 0.255210669955677150, + 0.258192911338648023, 0.261187919133763713, 0.264195763998317568, 0.267216518344631837, + 0.270250256366959984, 0.273297054069675804, 0.276356989296781264, 0.279430141762765316, + 0.282516593084849388, 0.285616426816658109, 0.288729728483353931, 0.291856585618280984, + 0.294997087801162572, 0.298151326697901342, 0.301319396102034120, 0.304501391977896274, + 0.307697412505553769, 0.310907558127563710, 0.314131931597630143, 0.317370638031222396, + 0.320623784958230129, 0.323891482377732021, 0.327173842814958593, 0.330470981380537099, + 0.333783015832108509, 0.337110066638412809, 0.340452257045945450, 0.343809713148291340, + 0.347182563958251478, 0.350570941482881204, 0.353974980801569250, 0.357394820147290515, + 0.360830600991175754, 0.364282468130549597, 0.367750569780596226, 0.371235057669821344, + 0.374736087139491414, 0.378253817247238111, 0.381788410875031348, 0.385340034841733958, + 0.388908860020464597, 0.392495061461010764, 0.396098818517547080, 0.399720314981931668, + 0.403359739222868885, 0.407017284331247953, 0.410693148271983222, 0.414387534042706784, + 0.418100649839684591, 0.421832709231353298, 0.425583931339900579, 0.429354541031341519, + 0.433144769114574058, 0.436954852549929273, 0.440785034667769915, 0.444635565397727750, + 0.448506701509214067, 0.452398706863882505, 0.456311852680773566, 0.460246417814923481, + 0.464202689050278838, 0.468180961407822172, 0.472181538469883255, 0.476204732721683788, + 0.480250865911249714, 0.484320269428911598, 0.488413284707712059, 0.492530263646148658, + 0.496671569054796314, 0.500837575128482149, 0.505028667945828791, 0.509245245998136142, + 0.513487720749743026, 0.517756517232200619, 0.522052074674794864, 0.526374847174186700, + 0.530725304406193921, 0.535103932383019565, 0.539511234259544614, 0.543947731192649941, + 0.548413963257921133, 0.552910490428519918, 0.557437893621486324, 0.561996775817277916, + 0.566587763258951771, 0.571211506738074970, 0.575868682975210544, 0.580559996103683473, + 0.585286179266300333, 0.590047996335791969, 0.594846243770991268, 0.599681752622167719, + 0.604555390700549533, 0.609468064928895381, 0.614420723892076803, 0.619414360609039205, + 0.624450015550274240, 0.629528779928128279, 0.634651799290960050, 0.639820277456438991, + 0.645035480824251883, 0.650298743114294586, 0.655611470583224665, 0.660975147780241357, + 0.666391343912380640, 0.671861719900766374, 0.677388036222513090, 0.682972161648791376, + 0.688616083008527058, 0.694321916130032579, 0.700091918140490099, 0.705928501336797409, + 0.711834248882358467, 0.717811932634901395, 0.723864533472881599, 0.729995264565802437, + 0.736207598131266683, 0.742505296344636245, 0.748892447223726720, 0.755373506511754500, + 0.761953346841546475, 0.768637315803334831, 0.775431304986138326, 0.782341832659861902, + 0.789376143571198563, 0.796542330428254619, 0.803849483176389490, 0.811307874318219935, + 0.818929191609414797, 0.826726833952094231, 0.834716292992930375, 0.842915653118441077, + 0.851346258465123684, 0.860033621203008636, 0.869008688043793165, 0.878309655816146839, + 0.887984660763399880, 0.898095921906304051, 0.908726440060562912, 0.919991505048360247, + 0.932060075968990209, 0.945198953453078028, 0.959879091812415930, 0.977101701282731328, + 1.000000000000000000]; +pub const ZIG_EXP_R: f64 = 7.697117470131050077; +pub static ZIG_EXP_X: [f64; 257] = + [8.697117470131052741, 7.697117470131050077, 6.941033629377212577, 6.478378493832569696, + 6.144164665772472667, 5.882144315795399869, 5.666410167454033697, 5.482890627526062488, + 5.323090505754398016, 5.181487281301500047, 5.054288489981304089, 4.938777085901250530, + 4.832939741025112035, 4.735242996601741083, 4.644491885420085175, 4.559737061707351380, + 4.480211746528421912, 4.405287693473573185, 4.334443680317273007, 4.267242480277365857, + 4.203313713735184365, 4.142340865664051464, 4.084051310408297830, 4.028208544647936762, + 3.974606066673788796, 3.923062500135489739, 3.873417670399509127, 3.825529418522336744, + 3.779270992411667862, 3.734528894039797375, 3.691201090237418825, 3.649195515760853770, + 3.608428813128909507, 3.568825265648337020, 3.530315889129343354, 3.492837654774059608, + 3.456332821132760191, 3.420748357251119920, 3.386035442460300970, 3.352149030900109405, + 3.319047470970748037, 3.286692171599068679, 3.255047308570449882, 3.224079565286264160, + 3.193757903212240290, 3.164053358025972873, 3.134938858084440394, 3.106389062339824481, + 3.078380215254090224, 3.050890016615455114, 3.023897504455676621, 2.997382949516130601, + 2.971327759921089662, 2.945714394895045718, 2.920526286512740821, 2.895747768600141825, + 2.871364012015536371, 2.847360965635188812, 2.823725302450035279, 2.800444370250737780, + 2.777506146439756574, 2.754899196562344610, 2.732612636194700073, 2.710636095867928752, + 2.688959688741803689, 2.667573980773266573, 2.646469963151809157, 2.625639026797788489, + 2.605072938740835564, 2.584763820214140750, 2.564704126316905253, 2.544886627111869970, + 2.525304390037828028, 2.505950763528594027, 2.486819361740209455, 2.467904050297364815, + 2.449198932978249754, 2.430698339264419694, 2.412396812688870629, 2.394289099921457886, + 2.376370140536140596, 2.358635057409337321, 2.341079147703034380, 2.323697874390196372, + 2.306486858283579799, 2.289441870532269441, 2.272558825553154804, 2.255833774367219213, + 2.239262898312909034, 2.222842503111036816, 2.206569013257663858, 2.190438966723220027, + 2.174449009937774679, 2.158595893043885994, 2.142876465399842001, 2.127287671317368289, + 2.111826546019042183, 2.096490211801715020, 2.081275874393225145, 2.066180819490575526, + 2.051202409468584786, 2.036338080248769611, 2.021585338318926173, 2.006941757894518563, + 1.992404978213576650, 1.977972700957360441, 1.963642687789548313, 1.949412758007184943, + 1.935280786297051359, 1.921244700591528076, 1.907302480018387536, 1.893452152939308242, + 1.879691795072211180, 1.866019527692827973, 1.852433515911175554, 1.838931967018879954, + 1.825513128903519799, 1.812175288526390649, 1.798916770460290859, 1.785735935484126014, + 1.772631179231305643, 1.759600930889074766, 1.746643651946074405, 1.733757834985571566, + 1.720942002521935299, 1.708194705878057773, 1.695514524101537912, 1.682900062917553896, + 1.670349953716452118, 1.657862852574172763, 1.645437439303723659, 1.633072416535991334, + 1.620766508828257901, 1.608518461798858379, 1.596327041286483395, 1.584191032532688892, + 1.572109239386229707, 1.560080483527888084, 1.548103603714513499, 1.536177455041032092, + 1.524300908219226258, 1.512472848872117082, 1.500692176842816750, 1.488957805516746058, + 1.477268661156133867, 1.465623682245745352, 1.454021818848793446, 1.442462031972012504, + 1.430943292938879674, 1.419464582769983219, 1.408024891569535697, 1.396623217917042137, + 1.385258568263121992, 1.373929956328490576, 1.362636402505086775, 1.351376933258335189, + 1.340150580529504643, 1.328956381137116560, 1.317793376176324749, 1.306660610415174117, + 1.295557131686601027, 1.284481990275012642, 1.273434238296241139, 1.262412929069615330, + 1.251417116480852521, 1.240445854334406572, 1.229498195693849105, 1.218573192208790124, + 1.207669893426761121, 1.196787346088403092, 1.185924593404202199, 1.175080674310911677, + 1.164254622705678921, 1.153445466655774743, 1.142652227581672841, 1.131873919411078511, + 1.121109547701330200, 1.110358108727411031, 1.099618588532597308, 1.088889961938546813, + 1.078171191511372307, 1.067461226479967662, 1.056759001602551429, 1.046063435977044209, + 1.035373431790528542, 1.024687873002617211, 1.014005623957096480, 1.003325527915696735, + 0.992646405507275897, 0.981967053085062602, 0.971286240983903260, 0.960602711668666509, + 0.949915177764075969, 0.939222319955262286, 0.928522784747210395, 0.917815182070044311, + 0.907098082715690257, 0.896370015589889935, 0.885629464761751528, 0.874874866291025066, + 0.864104604811004484, 0.853317009842373353, 0.842510351810368485, 0.831682837734273206, + 0.820832606554411814, 0.809957724057418282, 0.799056177355487174, 0.788125868869492430, + 0.777164609759129710, 0.766170112735434672, 0.755139984181982249, 0.744071715500508102, + 0.732962673584365398, 0.721810090308756203, 0.710611050909655040, 0.699362481103231959, + 0.688061132773747808, 0.676703568029522584, 0.665286141392677943, 0.653804979847664947, + 0.642255960424536365, 0.630634684933490286, 0.618936451394876075, 0.607156221620300030, + 0.595288584291502887, 0.583327712748769489, 0.571267316532588332, 0.559100585511540626, + 0.546820125163310577, 0.534417881237165604, 0.521885051592135052, 0.509211982443654398, + 0.496388045518671162, 0.483401491653461857, 0.470239275082169006, 0.456886840931420235, + 0.443327866073552401, 0.429543940225410703, 0.415514169600356364, 0.401214678896277765, + 0.386617977941119573, 0.371692145329917234, 0.356399760258393816, 0.340696481064849122, + 0.324529117016909452, 0.307832954674932158, 0.290527955491230394, 0.272513185478464703, + 0.253658363385912022, 0.233790483059674731, 0.212671510630966620, 0.189958689622431842, + 0.165127622564187282, 0.137304980940012589, 0.104838507565818778, 0.063852163815001570, + 0.000000000000000000]; +pub static ZIG_EXP_F: [f64; 257] = + [0.000167066692307963, 0.000454134353841497, 0.000967269282327174, 0.001536299780301573, + 0.002145967743718907, 0.002788798793574076, 0.003460264777836904, 0.004157295120833797, + 0.004877655983542396, 0.005619642207205489, 0.006381905937319183, 0.007163353183634991, + 0.007963077438017043, 0.008780314985808977, 0.009614413642502212, 0.010464810181029981, + 0.011331013597834600, 0.012212592426255378, 0.013109164931254991, 0.014020391403181943, + 0.014945968011691148, 0.015885621839973156, 0.016839106826039941, 0.017806200410911355, + 0.018786700744696024, 0.019780424338009740, 0.020787204072578114, 0.021806887504283581, + 0.022839335406385240, 0.023884420511558174, 0.024942026419731787, 0.026012046645134221, + 0.027094383780955803, 0.028188948763978646, 0.029295660224637411, 0.030414443910466622, + 0.031545232172893622, 0.032687963508959555, 0.033842582150874358, 0.035009037697397431, + 0.036187284781931443, 0.037377282772959382, 0.038578995503074871, 0.039792391023374139, + 0.041017441380414840, 0.042254122413316254, 0.043502413568888197, 0.044762297732943289, + 0.046033761076175184, 0.047316792913181561, 0.048611385573379504, 0.049917534282706379, + 0.051235237055126281, 0.052564494593071685, 0.053905310196046080, 0.055257689676697030, + 0.056621641283742870, 0.057997175631200659, 0.059384305633420280, 0.060783046445479660, + 0.062193415408541036, 0.063615431999807376, 0.065049117786753805, 0.066494496385339816, + 0.067951593421936643, 0.069420436498728783, 0.070901055162371843, 0.072393480875708752, + 0.073897746992364746, 0.075413888734058410, 0.076941943170480517, 0.078481949201606435, + 0.080033947542319905, 0.081597980709237419, 0.083174093009632397, 0.084762330532368146, + 0.086362741140756927, 0.087975374467270231, 0.089600281910032886, 0.091237516631040197, + 0.092887133556043569, 0.094549189376055873, 0.096223742550432825, 0.097910853311492213, + 0.099610583670637132, 0.101322997425953631, 0.103048160171257702, 0.104786139306570145, + 0.106537004050001632, 0.108300825451033755, 0.110077676405185357, 0.111867631670056283, + 0.113670767882744286, 0.115487163578633506, 0.117316899211555525, 0.119160057175327641, + 0.121016721826674792, 0.122886979509545108, 0.124770918580830933, 0.126668629437510671, + 0.128580204545228199, 0.130505738468330773, 0.132445327901387494, 0.134399071702213602, + 0.136367070926428829, 0.138349428863580176, 0.140346251074862399, 0.142357645432472146, + 0.144383722160634720, 0.146424593878344889, 0.148480375643866735, 0.150551185001039839, + 0.152637142027442801, 0.154738369384468027, 0.156854992369365148, 0.158987138969314129, + 0.161134939917591952, 0.163298528751901734, 0.165478041874935922, 0.167673618617250081, + 0.169885401302527550, 0.172113535315319977, 0.174358169171353411, 0.176619454590494829, + 0.178897546572478278, 0.181192603475496261, 0.183504787097767436, 0.185834262762197083, + 0.188181199404254262, 0.190545769663195363, 0.192928149976771296, 0.195328520679563189, + 0.197747066105098818, 0.200183974691911210, 0.202639439093708962, 0.205113656293837654, + 0.207606827724221982, 0.210119159388988230, 0.212650861992978224, 0.215202151075378628, + 0.217773247148700472, 0.220364375843359439, 0.222975768058120111, 0.225607660116683956, + 0.228260293930716618, 0.230933917169627356, 0.233628783437433291, 0.236345152457059560, + 0.239083290262449094, 0.241843469398877131, 0.244625969131892024, 0.247431075665327543, + 0.250259082368862240, 0.253110290015629402, 0.255985007030415324, 0.258883549749016173, + 0.261806242689362922, 0.264753418835062149, 0.267725419932044739, 0.270722596799059967, + 0.273745309652802915, 0.276793928448517301, 0.279868833236972869, 0.282970414538780746, + 0.286099073737076826, 0.289255223489677693, 0.292439288161892630, 0.295651704281261252, + 0.298892921015581847, 0.302163400675693528, 0.305463619244590256, 0.308794066934560185, + 0.312155248774179606, 0.315547685227128949, 0.318971912844957239, 0.322428484956089223, + 0.325917972393556354, 0.329440964264136438, 0.332998068761809096, 0.336589914028677717, + 0.340217149066780189, 0.343880444704502575, 0.347580494621637148, 0.351318016437483449, + 0.355093752866787626, 0.358908472948750001, 0.362762973354817997, 0.366658079781514379, + 0.370594648435146223, 0.374573567615902381, 0.378595759409581067, 0.382662181496010056, + 0.386773829084137932, 0.390931736984797384, 0.395136981833290435, 0.399390684475231350, + 0.403694012530530555, 0.408048183152032673, 0.412454465997161457, 0.416914186433003209, + 0.421428728997616908, 0.425999541143034677, 0.430628137288459167, 0.435316103215636907, + 0.440065100842354173, 0.444876873414548846, 0.449753251162755330, 0.454696157474615836, + 0.459707615642138023, 0.464789756250426511, 0.469944825283960310, 0.475175193037377708, + 0.480483363930454543, 0.485871987341885248, 0.491343869594032867, 0.496901987241549881, + 0.502549501841348056, 0.508289776410643213, 0.514126393814748894, 0.520063177368233931, + 0.526104213983620062, 0.532253880263043655, 0.538516872002862246, 0.544898237672440056, + 0.551403416540641733, 0.558038282262587892, 0.564809192912400615, 0.571723048664826150, + 0.578787358602845359, 0.586010318477268366, 0.593400901691733762, 0.600968966365232560, + 0.608725382079622346, 0.616682180915207878, 0.624852738703666200, 0.633251994214366398, + 0.641896716427266423, 0.650805833414571433, 0.660000841079000145, 0.669506316731925177, + 0.679350572264765806, 0.689566496117078431, 0.700192655082788606, 0.711274760805076456, + 0.722867659593572465, 0.735038092431424039, 0.747868621985195658, 0.761463388849896838, + 0.775956852040116218, 0.791527636972496285, 0.808421651523009044, 0.826993296643051101, + 0.847785500623990496, 0.871704332381204705, 0.900469929925747703, 0.938143680862176477, + 1.000000000000000000]; diff --git a/rand/src/jitter.rs b/rand/src/jitter.rs new file mode 100644 index 0000000..3693481 --- /dev/null +++ b/rand/src/jitter.rs @@ -0,0 +1,754 @@ +// Copyright 2017 The Rust Project Developers. See the COPYRIGHT +// file at the top-level directory of this distribution and at +// http://rust-lang.org/COPYRIGHT. +// +// Licensed under the Apache License, Version 2.0 <LICENSE-APACHE or +// http://www.apache.org/licenses/LICENSE-2.0> or the MIT license +// <LICENSE-MIT or http://opensource.org/licenses/MIT>, at your +// option. This file may not be copied, modified, or distributed +// except according to those terms. +// +// Based on jitterentropy-library, http://www.chronox.de/jent.html. +// Copyright Stephan Mueller <smueller@chronox.de>, 2014 - 2017. +// +// With permission from Stephan Mueller to relicense the Rust translation under +// the MIT license. + +//! Non-physical true random number generator based on timing jitter. + +use Rng; + +use core::{fmt, mem, ptr}; +#[cfg(feature="std")] +use std::sync::atomic::{AtomicUsize, ATOMIC_USIZE_INIT, Ordering}; + +const MEMORY_BLOCKS: usize = 64; +const MEMORY_BLOCKSIZE: usize = 32; +const MEMORY_SIZE: usize = MEMORY_BLOCKS * MEMORY_BLOCKSIZE; + +/// A true random number generator based on jitter in the CPU execution time, +/// and jitter in memory access time. +/// +/// This is a true random number generator, as opposed to pseudo-random +/// generators. Random numbers generated by `JitterRng` can be seen as fresh +/// entropy. A consequence is that is orders of magnitude slower than `OsRng` +/// and PRNGs (about 10^3 .. 10^6 slower). +/// +/// There are very few situations where using this RNG is appropriate. Only very +/// few applications require true entropy. A normal PRNG can be statistically +/// indistinguishable, and a cryptographic PRNG should also be as impossible to +/// predict. +/// +/// Use of `JitterRng` is recommended for initializing cryptographic PRNGs when +/// `OsRng` is not available. +/// +/// This implementation is based on +/// [Jitterentropy](http://www.chronox.de/jent.html) version 2.1.0. +// +// Note: the C implementation relies on being compiled without optimizations. +// This implementation goes through lengths to make the compiler not optimise +// out what is technically dead code, but that does influence timing jitter. +pub struct JitterRng { + data: u64, // Actual random number + // Number of rounds to run the entropy collector per 64 bits + rounds: u32, + // Timer and previous time stamp, used by `measure_jitter` + timer: fn() -> u64, + prev_time: u64, + // Deltas used for the stuck test + last_delta: i64, + last_delta2: i64, + // Memory for the Memory Access noise source + mem_prev_index: usize, + mem: [u8; MEMORY_SIZE], + // Make `next_u32` not waste 32 bits + data_remaining: Option<u32>, +} + +// Custom Debug implementation that does not expose the internal state +impl fmt::Debug for JitterRng { + fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result { + write!(f, "JitterRng {{}}") + } +} + +/// An error that can occur when `test_timer` fails. +#[derive(Debug, Clone, PartialEq, Eq)] +pub enum TimerError { + /// No timer available. + NoTimer, + /// Timer too coarse to use as an entropy source. + CoarseTimer, + /// Timer is not monotonically increasing. + NotMonotonic, + /// Variations of deltas of time too small. + TinyVariantions, + /// Too many stuck results (indicating no added entropy). + TooManyStuck, + #[doc(hidden)] + __Nonexhaustive, +} + +impl TimerError { + fn description(&self) -> &'static str { + match *self { + TimerError::NoTimer => "no timer available", + TimerError::CoarseTimer => "coarse timer", + TimerError::NotMonotonic => "timer not monotonic", + TimerError::TinyVariantions => "time delta variations too small", + TimerError::TooManyStuck => "too many stuck results", + TimerError::__Nonexhaustive => unreachable!(), + } + } +} + +impl fmt::Display for TimerError { + fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result { + write!(f, "{}", self.description()) + } +} + +#[cfg(feature="std")] +impl ::std::error::Error for TimerError { + fn description(&self) -> &str { + self.description() + } +} + +// Initialise to zero; must be positive +#[cfg(feature="std")] +static JITTER_ROUNDS: AtomicUsize = ATOMIC_USIZE_INIT; + +impl JitterRng { + /// Create a new `JitterRng`. + /// Makes use of `std::time` for a timer. + /// + /// During initialization CPU execution timing jitter is measured a few + /// hundred times. If this does not pass basic quality tests, an error is + /// returned. The test result is cached to make subsequent calls faster. + #[cfg(feature="std")] + pub fn new() -> Result<JitterRng, TimerError> { + let mut ec = JitterRng::new_with_timer(platform::get_nstime); + let mut rounds = JITTER_ROUNDS.load(Ordering::Relaxed) as u32; + if rounds == 0 { + // No result yet: run test. + // This allows the timer test to run multiple times; we don't care. + rounds = ec.test_timer()?; + JITTER_ROUNDS.store(rounds as usize, Ordering::Relaxed); + } + ec.set_rounds(rounds); + Ok(ec) + } + + /// Create a new `JitterRng`. + /// A custom timer can be supplied, making it possible to use `JitterRng` in + /// `no_std` environments. + /// + /// The timer must have nanosecond precision. + /// + /// This method is more low-level than `new()`. It is the responsibility of + /// the caller to run `test_timer` before using any numbers generated with + /// `JitterRng`, and optionally call `set_rounds()`. + pub fn new_with_timer(timer: fn() -> u64) -> JitterRng { + let mut ec = JitterRng { + data: 0, + rounds: 64, + timer: timer, + prev_time: 0, + last_delta: 0, + last_delta2: 0, + mem_prev_index: 0, + mem: [0; MEMORY_SIZE], + data_remaining: None, + }; + + // Fill `data`, `prev_time`, `last_delta` and `last_delta2` with + // non-zero values. + ec.prev_time = timer(); + ec.gen_entropy(); + + // Do a single read from `self.mem` to make sure the Memory Access noise + // source is not optimised out. + // Note: this read is important, it effects optimisations for the entire + // module! + black_box(ec.mem[0]); + + ec + } + + /// Configures how many rounds are used to generate each 64-bit value. + /// This must be greater than zero, and has a big impact on performance + /// and output quality. + /// + /// `new_with_timer` conservatively uses 64 rounds, but often less rounds + /// can be used. The `test_timer()` function returns the minimum number of + /// rounds required for full strength (platform dependent), so one may use + /// `rng.set_rounds(rng.test_timer()?);` or cache the value. + pub fn set_rounds(&mut self, rounds: u32) { + assert!(rounds > 0); + self.rounds = rounds; + } + + // Calculate a random loop count used for the next round of an entropy + // collection, based on bits from a fresh value from the timer. + // + // The timer is folded to produce a number that contains at most `n_bits` + // bits. + // + // Note: A constant should be added to the resulting random loop count to + // prevent loops that run 0 times. + #[inline(never)] + fn random_loop_cnt(&mut self, n_bits: u32) -> u32 { + let mut rounds = 0; + + let mut time = (self.timer)(); + // Mix with the current state of the random number balance the random + // loop counter a bit more. + time ^= self.data; + + // We fold the time value as much as possible to ensure that as many + // bits of the time stamp are included as possible. + let folds = (64 + n_bits - 1) / n_bits; + let mask = (1 << n_bits) - 1; + for _ in 0..folds { + rounds ^= time & mask; + time = time >> n_bits; + } + + rounds as u32 + } + + // CPU jitter noise source + // Noise source based on the CPU execution time jitter + // + // This function injects the individual bits of the time value into the + // entropy pool using an LFSR. + // + // The code is deliberately inefficient with respect to the bit shifting. + // This function not only acts as folding operation, but this function's + // execution is used to measure the CPU execution time jitter. Any change to + // the loop in this function implies that careful retesting must be done. + #[inline(never)] + fn lfsr_time(&mut self, time: u64, var_rounds: bool) { + fn lfsr(mut data: u64, time: u64) -> u64{ + for i in 1..65 { + let mut tmp = time << (64 - i); + tmp = tmp >> (64 - 1); + + // Fibonacci LSFR with polynomial of + // x^64 + x^61 + x^56 + x^31 + x^28 + x^23 + 1 which is + // primitive according to + // http://poincare.matf.bg.ac.rs/~ezivkovm/publications/primpol1.pdf + // (the shift values are the polynomial values minus one + // due to counting bits from 0 to 63). As the current + // position is always the LSB, the polynomial only needs + // to shift data in from the left without wrap. + data ^= tmp; + data ^= (data >> 63) & 1; + data ^= (data >> 60) & 1; + data ^= (data >> 55) & 1; + data ^= (data >> 30) & 1; + data ^= (data >> 27) & 1; + data ^= (data >> 22) & 1; + data = data.rotate_left(1); + } + data + } + + // Note: in the reference implementation only the last round effects + // `self.data`, all the other results are ignored. To make sure the + // other rounds are not optimised out, we first run all but the last + // round on a throw-away value instead of the real `self.data`. + let mut lfsr_loop_cnt = 0; + if var_rounds { lfsr_loop_cnt = self.random_loop_cnt(4) }; + + let mut throw_away: u64 = 0; + for _ in 0..lfsr_loop_cnt { + throw_away = lfsr(throw_away, time); + } + black_box(throw_away); + + self.data = lfsr(self.data, time); + } + + // Memory Access noise source + // This is a noise source based on variations in memory access times + // + // This function performs memory accesses which will add to the timing + // variations due to an unknown amount of CPU wait states that need to be + // added when accessing memory. The memory size should be larger than the L1 + // caches as outlined in the documentation and the associated testing. + // + // The L1 cache has a very high bandwidth, albeit its access rate is usually + // slower than accessing CPU registers. Therefore, L1 accesses only add + // minimal variations as the CPU has hardly to wait. Starting with L2, + // significant variations are added because L2 typically does not belong to + // the CPU any more and therefore a wider range of CPU wait states is + // necessary for accesses. L3 and real memory accesses have even a wider + // range of wait states. However, to reliably access either L3 or memory, + // the `self.mem` memory must be quite large which is usually not desirable. + #[inline(never)] + fn memaccess(&mut self, var_rounds: bool) { + let mut acc_loop_cnt = 128; + if var_rounds { acc_loop_cnt += self.random_loop_cnt(4) }; + + let mut index = self.mem_prev_index; + for _ in 0..acc_loop_cnt { + // Addition of memblocksize - 1 to index with wrap around logic to + // ensure that every memory location is hit evenly. + // The modulus also allows the compiler to remove the indexing + // bounds check. + index = (index + MEMORY_BLOCKSIZE - 1) % MEMORY_SIZE; + + // memory access: just add 1 to one byte + // memory access implies read from and write to memory location + let tmp = self.mem[index]; + self.mem[index] = tmp.wrapping_add(1); + } + self.mem_prev_index = index; + } + + + // Stuck test by checking the: + // - 1st derivation of the jitter measurement (time delta) + // - 2nd derivation of the jitter measurement (delta of time deltas) + // - 3rd derivation of the jitter measurement (delta of delta of time + // deltas) + // + // All values must always be non-zero. + // This test is a heuristic to see whether the last measurement holds + // entropy. + fn stuck(&mut self, current_delta: i64) -> bool { + let delta2 = self.last_delta - current_delta; + let delta3 = delta2 - self.last_delta2; + + self.last_delta = current_delta; + self.last_delta2 = delta2; + + current_delta == 0 || delta2 == 0 || delta3 == 0 + } + + // This is the heart of the entropy generation: calculate time deltas and + // use the CPU jitter in the time deltas. The jitter is injected into the + // entropy pool. + // + // Ensure that `self.prev_time` is primed before using the output of this + // function. This can be done by calling this function and not using its + // result. + fn measure_jitter(&mut self) -> Option<()> { + // Invoke one noise source before time measurement to add variations + self.memaccess(true); + + // Get time stamp and calculate time delta to previous + // invocation to measure the timing variations + let time = (self.timer)(); + // Note: wrapping_sub combined with a cast to `i64` generates a correct + // delta, even in the unlikely case this is a timer that is not strictly + // monotonic. + let current_delta = time.wrapping_sub(self.prev_time) as i64; + self.prev_time = time; + + // Call the next noise source which also injects the data + self.lfsr_time(current_delta as u64, true); + + // Check whether we have a stuck measurement (i.e. does the last + // measurement holds entropy?). + if self.stuck(current_delta) { return None }; + + // Rotate the data buffer by a prime number (any odd number would + // do) to ensure that every bit position of the input time stamp + // has an even chance of being merged with a bit position in the + // entropy pool. We do not use one here as the adjacent bits in + // successive time deltas may have some form of dependency. The + // chosen value of 7 implies that the low 7 bits of the next + // time delta value is concatenated with the current time delta. + self.data = self.data.rotate_left(7); + + Some(()) + } + + // Shuffle the pool a bit by mixing some value with a bijective function + // (XOR) into the pool. + // + // The function generates a mixer value that depends on the bits set and + // the location of the set bits in the random number generated by the + // entropy source. Therefore, based on the generated random number, this + // mixer value can have 2^64 different values. That mixer value is + // initialized with the first two SHA-1 constants. After obtaining the + // mixer value, it is XORed into the random number. + // + // The mixer value is not assumed to contain any entropy. But due to the + // XOR operation, it can also not destroy any entropy present in the + // entropy pool. + #[inline(never)] + fn stir_pool(&mut self) { + // This constant is derived from the first two 32 bit initialization + // vectors of SHA-1 as defined in FIPS 180-4 section 5.3.1 + // The order does not really matter as we do not rely on the specific + // numbers. We just pick the SHA-1 constants as they have a good mix of + // bit set and unset. + const CONSTANT: u64 = 0x67452301efcdab89; + + // The start value of the mixer variable is derived from the third + // and fourth 32 bit initialization vector of SHA-1 as defined in + // FIPS 180-4 section 5.3.1 + let mut mixer = 0x98badcfe10325476; + + // This is a constant time function to prevent leaking timing + // information about the random number. + // The normal code is: + // ``` + // for i in 0..64 { + // if ((self.data >> i) & 1) == 1 { mixer ^= CONSTANT; } + // } + // ``` + // This is a bit fragile, as LLVM really wants to use branches here, and + // we rely on it to not recognise the opportunity. + for i in 0..64 { + let apply = (self.data >> i) & 1; + let mask = !apply.wrapping_sub(1); + mixer ^= CONSTANT & mask; + mixer = mixer.rotate_left(1); + } + + self.data ^= mixer; + } + + fn gen_entropy(&mut self) -> u64 { + // Prime `self.prev_time`, and run the noice sources to make sure the + // first loop round collects the expected entropy. + let _ = self.measure_jitter(); + + for _ in 0..self.rounds { + // If a stuck measurement is received, repeat measurement + // Note: we do not guard against an infinite loop, that would mean + // the timer suddenly became broken. + while self.measure_jitter().is_none() {} + } + + self.stir_pool(); + self.data + } + + /// Basic quality tests on the timer, by measuring CPU timing jitter a few + /// hundred times. + /// + /// If succesful, this will return the estimated number of rounds necessary + /// to collect 64 bits of entropy. Otherwise a `TimerError` with the cause + /// of the failure will be returned. + pub fn test_timer(&mut self) -> Result<u32, TimerError> { + // We could add a check for system capabilities such as `clock_getres` + // or check for `CONFIG_X86_TSC`, but it does not make much sense as the + // following sanity checks verify that we have a high-resolution timer. + + #[cfg(all(target_arch = "wasm32", not(target_os = "emscripten")))] + return Err(TimerError::NoTimer); + + let mut delta_sum = 0; + let mut old_delta = 0; + + let mut time_backwards = 0; + let mut count_mod = 0; + let mut count_stuck = 0; + + // TESTLOOPCOUNT needs some loops to identify edge systems. + // 100 is definitely too little. + const TESTLOOPCOUNT: u64 = 300; + const CLEARCACHE: u64 = 100; + + for i in 0..(CLEARCACHE + TESTLOOPCOUNT) { + // Measure time delta of core entropy collection logic + let time = (self.timer)(); + self.memaccess(true); + self.lfsr_time(time, true); + let time2 = (self.timer)(); + + // Test whether timer works + if time == 0 || time2 == 0 { + return Err(TimerError::NoTimer); + } + let delta = time2.wrapping_sub(time) as i64; + + // Test whether timer is fine grained enough to provide delta even + // when called shortly after each other -- this implies that we also + // have a high resolution timer + if delta == 0 { + return Err(TimerError::CoarseTimer); + } + + // Up to here we did not modify any variable that will be + // evaluated later, but we already performed some work. Thus we + // already have had an impact on the caches, branch prediction, + // etc. with the goal to clear it to get the worst case + // measurements. + if i < CLEARCACHE { continue; } + + if self.stuck(delta) { count_stuck += 1; } + + // Test whether we have an increasing timer. + if !(time2 > time) { time_backwards += 1; } + + // Count the number of times the counter increases in steps of 100ns + // or greater. + if (delta % 100) == 0 { count_mod += 1; } + + // Ensure that we have a varying delta timer which is necessary for + // the calculation of entropy -- perform this check only after the + // first loop is executed as we need to prime the old_delta value + delta_sum += (delta - old_delta).abs() as u64; + old_delta = delta; + } + + // We allow the time to run backwards for up to three times. + // This can happen if the clock is being adjusted by NTP operations. + // If such an operation just happens to interfere with our test, it + // should not fail. The value of 3 should cover the NTP case being + // performed during our test run. + if time_backwards > 3 { + return Err(TimerError::NotMonotonic); + } + + // Test that the available amount of entropy per round does not get to + // low. We expect 1 bit of entropy per round as a reasonable minimum + // (although less is possible, it means the collector loop has to run + // much more often). + // `assert!(delta_average >= log2(1))` + // `assert!(delta_sum / TESTLOOPCOUNT >= 1)` + // `assert!(delta_sum >= TESTLOOPCOUNT)` + if delta_sum < TESTLOOPCOUNT { + return Err(TimerError::TinyVariantions); + } + + // Ensure that we have variations in the time stamp below 100 for at + // least 10% of all checks -- on some platforms, the counter increments + // in multiples of 100, but not always + if count_mod > (TESTLOOPCOUNT * 9 / 10) { + return Err(TimerError::CoarseTimer); + } + + // If we have more than 90% stuck results, then this Jitter RNG is + // likely to not work well. + if count_stuck > (TESTLOOPCOUNT * 9 / 10) { + return Err(TimerError::TooManyStuck); + } + + // Estimate the number of `measure_jitter` rounds necessary for 64 bits + // of entropy. + // + // We don't try very hard to come up with a good estimate of the + // available bits of entropy per round here for two reasons: + // 1. Simple estimates of the available bits (like Shannon entropy) are + // too optimistic. + // 2) Unless we want to waste a lot of time during intialization, there + // only a small number of samples are available. + // + // Therefore we use a very simple and conservative estimate: + // `let bits_of_entropy = log2(delta_average) / 2`. + // + // The number of rounds `measure_jitter` should run to collect 64 bits + // of entropy is `64 / bits_of_entropy`. + // + // To have smaller rounding errors, intermediate values are multiplied + // by `FACTOR`. To compensate for `log2` and division rounding down, + // add 1. + let delta_average = delta_sum / TESTLOOPCOUNT; + // println!("delta_average: {}", delta_average); + + const FACTOR: u32 = 3; + fn log2(x: u64) -> u32 { 64 - x.leading_zeros() } + + // pow(δ, FACTOR) must be representable; if you have overflow reduce FACTOR + Ok(64 * 2 * FACTOR / (log2(delta_average.pow(FACTOR)) + 1)) + } + + /// Statistical test: return the timer delta of one normal run of the + /// `JitterEntropy` entropy collector. + /// + /// Setting `var_rounds` to `true` will execute the memory access and the + /// CPU jitter noice sources a variable amount of times (just like a real + /// `JitterEntropy` round). + /// + /// Setting `var_rounds` to `false` will execute the noice sources the + /// minimal number of times. This can be used to measure the minimum amount + /// of entropy one round of entropy collector can collect in the worst case. + /// + /// # Example + /// + /// Use `timer_stats` to run the [NIST SP 800-90B Entropy Estimation Suite] + /// (https://github.com/usnistgov/SP800-90B_EntropyAssessment). + /// + /// This is the recommended way to test the quality of `JitterRng`. It + /// should be run before using the RNG on untested hardware, after changes + /// that could effect how the code is optimised, and after major compiler + /// compiler changes, like a new LLVM version. + /// + /// First generate two files `jitter_rng_var.bin` and `jitter_rng_var.min`. + /// + /// Execute `python noniid_main.py -v jitter_rng_var.bin 8`, and validate it + /// with `restart.py -v jitter_rng_var.bin 8 <min-entropy>`. + /// This number is the expected amount of entropy that is at least available + /// for each round of the entropy collector. This number should be greater + /// than the amount estimated with `64 / test_timer()`. + /// + /// Execute `python noniid_main.py -v -u 4 jitter_rng_var.bin 4`, and + /// validate it with `restart.py -v -u 4 jitter_rng_var.bin 4 <min-entropy>`. + /// This number is the expected amount of entropy that is available in the + /// last 4 bits of the timer delta after running noice sources. Note that + /// a value of 3.70 is the minimum estimated entropy for true randomness. + /// + /// Execute `python noniid_main.py -v -u 4 jitter_rng_var.bin 4`, and + /// validate it with `restart.py -v -u 4 jitter_rng_var.bin 4 <min-entropy>`. + /// This number is the expected amount of entropy that is available to the + /// entropy collecter if both noice sources only run their minimal number of + /// times. This measures the absolute worst-case, and gives a lower bound + /// for the available entropy. + /// + /// ```rust,no_run + /// use rand::JitterRng; + /// + /// # use std::error::Error; + /// # use std::fs::File; + /// # use std::io::Write; + /// # + /// # fn try_main() -> Result<(), Box<Error>> { + /// fn get_nstime() -> u64 { + /// use std::time::{SystemTime, UNIX_EPOCH}; + /// + /// let dur = SystemTime::now().duration_since(UNIX_EPOCH).unwrap(); + /// // The correct way to calculate the current time is + /// // `dur.as_secs() * 1_000_000_000 + dur.subsec_nanos() as u64` + /// // But this is faster, and the difference in terms of entropy is + /// // negligible (log2(10^9) == 29.9). + /// dur.as_secs() << 30 | dur.subsec_nanos() as u64 + /// } + /// + /// // Do not initialize with `JitterRng::new`, but with `new_with_timer`. + /// // 'new' always runst `test_timer`, and can therefore fail to + /// // initialize. We want to be able to get the statistics even when the + /// // timer test fails. + /// let mut rng = JitterRng::new_with_timer(get_nstime); + /// + /// // 1_000_000 results are required for the NIST SP 800-90B Entropy + /// // Estimation Suite + /// // FIXME: this number is smaller here, otherwise the Doc-test is too slow + /// const ROUNDS: usize = 10_000; + /// let mut deltas_variable: Vec<u8> = Vec::with_capacity(ROUNDS); + /// let mut deltas_minimal: Vec<u8> = Vec::with_capacity(ROUNDS); + /// + /// for _ in 0..ROUNDS { + /// deltas_variable.push(rng.timer_stats(true) as u8); + /// deltas_minimal.push(rng.timer_stats(false) as u8); + /// } + /// + /// // Write out after the statistics collection loop, to not disturb the + /// // test results. + /// File::create("jitter_rng_var.bin")?.write(&deltas_variable)?; + /// File::create("jitter_rng_min.bin")?.write(&deltas_minimal)?; + /// # + /// # Ok(()) + /// # } + /// # + /// # fn main() { + /// # try_main().unwrap(); + /// # } + /// ``` + #[cfg(feature="std")] + pub fn timer_stats(&mut self, var_rounds: bool) -> i64 { + let time = platform::get_nstime(); + self.memaccess(var_rounds); + self.lfsr_time(time, var_rounds); + let time2 = platform::get_nstime(); + time2.wrapping_sub(time) as i64 + } +} + +#[cfg(feature="std")] +mod platform { + #[cfg(not(any(target_os = "macos", target_os = "ios", target_os = "windows", all(target_arch = "wasm32", not(target_os = "emscripten")))))] + pub fn get_nstime() -> u64 { + use std::time::{SystemTime, UNIX_EPOCH}; + + let dur = SystemTime::now().duration_since(UNIX_EPOCH).unwrap(); + // The correct way to calculate the current time is + // `dur.as_secs() * 1_000_000_000 + dur.subsec_nanos() as u64` + // But this is faster, and the difference in terms of entropy is negligible + // (log2(10^9) == 29.9). + dur.as_secs() << 30 | dur.subsec_nanos() as u64 + } + + #[cfg(any(target_os = "macos", target_os = "ios"))] + pub fn get_nstime() -> u64 { + extern crate libc; + // On Mac OS and iOS std::time::SystemTime only has 1000ns resolution. + // We use `mach_absolute_time` instead. This provides a CPU dependent unit, + // to get real nanoseconds the result should by multiplied by numer/denom + // from `mach_timebase_info`. + // But we are not interested in the exact nanoseconds, just entropy. So we + // use the raw result. + unsafe { libc::mach_absolute_time() } + } + + #[cfg(target_os = "windows")] + pub fn get_nstime() -> u64 { + extern crate winapi; + unsafe { + let mut t = super::mem::zeroed(); + winapi::um::profileapi::QueryPerformanceCounter(&mut t); + *t.QuadPart() as u64 + } + } + + #[cfg(all(target_arch = "wasm32", not(target_os = "emscripten")))] + pub fn get_nstime() -> u64 { + unreachable!() + } +} + +// A function that is opaque to the optimizer to assist in avoiding dead-code +// elimination. Taken from `bencher`. +fn black_box<T>(dummy: T) -> T { + unsafe { + let ret = ptr::read_volatile(&dummy); + mem::forget(dummy); + ret + } +} + +impl Rng for JitterRng { + fn next_u32(&mut self) -> u32 { + // We want to use both parts of the generated entropy + if let Some(high) = self.data_remaining.take() { + high + } else { + let data = self.next_u64(); + self.data_remaining = Some((data >> 32) as u32); + data as u32 + } + } + + fn next_u64(&mut self) -> u64 { + self.gen_entropy() + } + + fn fill_bytes(&mut self, dest: &mut [u8]) { + let mut left = dest; + while left.len() >= 8 { + let (l, r) = {left}.split_at_mut(8); + left = r; + let chunk: [u8; 8] = unsafe { + mem::transmute(self.next_u64().to_le()) + }; + l.copy_from_slice(&chunk); + } + let n = left.len(); + if n > 0 { + let chunk: [u8; 8] = unsafe { + mem::transmute(self.next_u64().to_le()) + }; + left.copy_from_slice(&chunk[..n]); + } + } +} + +// There are no tests included because (1) this is an "external" RNG, so output +// is not reproducible and (2) `test_timer` *will* fail on some platforms. diff --git a/rand/src/lib.rs b/rand/src/lib.rs new file mode 100644 index 0000000..7b22dd4 --- /dev/null +++ b/rand/src/lib.rs @@ -0,0 +1,1214 @@ +// Copyright 2013-2017 The Rust Project Developers. See the COPYRIGHT +// file at the top-level directory of this distribution and at +// http://rust-lang.org/COPYRIGHT. +// +// Licensed under the Apache License, Version 2.0 <LICENSE-APACHE or +// http://www.apache.org/licenses/LICENSE-2.0> or the MIT license +// <LICENSE-MIT or http://opensource.org/licenses/MIT>, at your +// option. This file may not be copied, modified, or distributed +// except according to those terms. + +//! Utilities for random number generation +//! +//! The key functions are `random()` and `Rng::gen()`. These are polymorphic and +//! so can be used to generate any type that implements `Rand`. Type inference +//! means that often a simple call to `rand::random()` or `rng.gen()` will +//! suffice, but sometimes an annotation is required, e.g. +//! `rand::random::<f64>()`. +//! +//! See the `distributions` submodule for sampling random numbers from +//! distributions like normal and exponential. +//! +//! # Usage +//! +//! This crate is [on crates.io](https://crates.io/crates/rand) and can be +//! used by adding `rand` to the dependencies in your project's `Cargo.toml`. +//! +//! ```toml +//! [dependencies] +//! rand = "0.4" +//! ``` +//! +//! and this to your crate root: +//! +//! ```rust +//! extern crate rand; +//! ``` +//! +//! # Thread-local RNG +//! +//! There is built-in support for a RNG associated with each thread stored +//! in thread-local storage. This RNG can be accessed via `thread_rng`, or +//! used implicitly via `random`. This RNG is normally randomly seeded +//! from an operating-system source of randomness, e.g. `/dev/urandom` on +//! Unix systems, and will automatically reseed itself from this source +//! after generating 32 KiB of random data. +//! +//! # Cryptographic security +//! +//! An application that requires an entropy source for cryptographic purposes +//! must use `OsRng`, which reads randomness from the source that the operating +//! system provides (e.g. `/dev/urandom` on Unixes or `CryptGenRandom()` on +//! Windows). +//! The other random number generators provided by this module are not suitable +//! for such purposes. +//! +//! *Note*: many Unix systems provide `/dev/random` as well as `/dev/urandom`. +//! This module uses `/dev/urandom` for the following reasons: +//! +//! - On Linux, `/dev/random` may block if entropy pool is empty; +//! `/dev/urandom` will not block. This does not mean that `/dev/random` +//! provides better output than `/dev/urandom`; the kernel internally runs a +//! cryptographically secure pseudorandom number generator (CSPRNG) based on +//! entropy pool for random number generation, so the "quality" of +//! `/dev/random` is not better than `/dev/urandom` in most cases. However, +//! this means that `/dev/urandom` can yield somewhat predictable randomness +//! if the entropy pool is very small, such as immediately after first +//! booting. Linux 3.17 added the `getrandom(2)` system call which solves +//! the issue: it blocks if entropy pool is not initialized yet, but it does +//! not block once initialized. `OsRng` tries to use `getrandom(2)` if +//! available, and use `/dev/urandom` fallback if not. If an application +//! does not have `getrandom` and likely to be run soon after first booting, +//! or on a system with very few entropy sources, one should consider using +//! `/dev/random` via `ReadRng`. +//! - On some systems (e.g. FreeBSD, OpenBSD and Mac OS X) there is no +//! difference between the two sources. (Also note that, on some systems +//! e.g. FreeBSD, both `/dev/random` and `/dev/urandom` may block once if +//! the CSPRNG has not seeded yet.) +//! +//! # Examples +//! +//! ```rust +//! use rand::Rng; +//! +//! let mut rng = rand::thread_rng(); +//! if rng.gen() { // random bool +//! println!("i32: {}, u32: {}", rng.gen::<i32>(), rng.gen::<u32>()) +//! } +//! ``` +//! +//! ```rust +//! let tuple = rand::random::<(f64, char)>(); +//! println!("{:?}", tuple) +//! ``` +//! +//! ## Monte Carlo estimation of π +//! +//! For this example, imagine we have a square with sides of length 2 and a unit +//! circle, both centered at the origin. Since the area of a unit circle is π, +//! we have: +//! +//! ```text +//! (area of unit circle) / (area of square) = π / 4 +//! ``` +//! +//! So if we sample many points randomly from the square, roughly π / 4 of them +//! should be inside the circle. +//! +//! We can use the above fact to estimate the value of π: pick many points in +//! the square at random, calculate the fraction that fall within the circle, +//! and multiply this fraction by 4. +//! +//! ``` +//! use rand::distributions::{IndependentSample, Range}; +//! +//! fn main() { +//! let between = Range::new(-1f64, 1.); +//! let mut rng = rand::thread_rng(); +//! +//! let total = 1_000_000; +//! let mut in_circle = 0; +//! +//! for _ in 0..total { +//! let a = between.ind_sample(&mut rng); +//! let b = between.ind_sample(&mut rng); +//! if a*a + b*b <= 1. { +//! in_circle += 1; +//! } +//! } +//! +//! // prints something close to 3.14159... +//! println!("{}", 4. * (in_circle as f64) / (total as f64)); +//! } +//! ``` +//! +//! ## Monty Hall Problem +//! +//! This is a simulation of the [Monty Hall Problem][]: +//! +//! > Suppose you're on a game show, and you're given the choice of three doors: +//! > Behind one door is a car; behind the others, goats. You pick a door, say +//! > No. 1, and the host, who knows what's behind the doors, opens another +//! > door, say No. 3, which has a goat. He then says to you, "Do you want to +//! > pick door No. 2?" Is it to your advantage to switch your choice? +//! +//! The rather unintuitive answer is that you will have a 2/3 chance of winning +//! if you switch and a 1/3 chance of winning if you don't, so it's better to +//! switch. +//! +//! This program will simulate the game show and with large enough simulation +//! steps it will indeed confirm that it is better to switch. +//! +//! [Monty Hall Problem]: http://en.wikipedia.org/wiki/Monty_Hall_problem +//! +//! ``` +//! use rand::Rng; +//! use rand::distributions::{IndependentSample, Range}; +//! +//! struct SimulationResult { +//! win: bool, +//! switch: bool, +//! } +//! +//! // Run a single simulation of the Monty Hall problem. +//! fn simulate<R: Rng>(random_door: &Range<u32>, rng: &mut R) +//! -> SimulationResult { +//! let car = random_door.ind_sample(rng); +//! +//! // This is our initial choice +//! let mut choice = random_door.ind_sample(rng); +//! +//! // The game host opens a door +//! let open = game_host_open(car, choice, rng); +//! +//! // Shall we switch? +//! let switch = rng.gen(); +//! if switch { +//! choice = switch_door(choice, open); +//! } +//! +//! SimulationResult { win: choice == car, switch: switch } +//! } +//! +//! // Returns the door the game host opens given our choice and knowledge of +//! // where the car is. The game host will never open the door with the car. +//! fn game_host_open<R: Rng>(car: u32, choice: u32, rng: &mut R) -> u32 { +//! let choices = free_doors(&[car, choice]); +//! rand::seq::sample_slice(rng, &choices, 1)[0] +//! } +//! +//! // Returns the door we switch to, given our current choice and +//! // the open door. There will only be one valid door. +//! fn switch_door(choice: u32, open: u32) -> u32 { +//! free_doors(&[choice, open])[0] +//! } +//! +//! fn free_doors(blocked: &[u32]) -> Vec<u32> { +//! (0..3).filter(|x| !blocked.contains(x)).collect() +//! } +//! +//! fn main() { +//! // The estimation will be more accurate with more simulations +//! let num_simulations = 10000; +//! +//! let mut rng = rand::thread_rng(); +//! let random_door = Range::new(0, 3); +//! +//! let (mut switch_wins, mut switch_losses) = (0, 0); +//! let (mut keep_wins, mut keep_losses) = (0, 0); +//! +//! println!("Running {} simulations...", num_simulations); +//! for _ in 0..num_simulations { +//! let result = simulate(&random_door, &mut rng); +//! +//! match (result.win, result.switch) { +//! (true, true) => switch_wins += 1, +//! (true, false) => keep_wins += 1, +//! (false, true) => switch_losses += 1, +//! (false, false) => keep_losses += 1, +//! } +//! } +//! +//! let total_switches = switch_wins + switch_losses; +//! let total_keeps = keep_wins + keep_losses; +//! +//! println!("Switched door {} times with {} wins and {} losses", +//! total_switches, switch_wins, switch_losses); +//! +//! println!("Kept our choice {} times with {} wins and {} losses", +//! total_keeps, keep_wins, keep_losses); +//! +//! // With a large number of simulations, the values should converge to +//! // 0.667 and 0.333 respectively. +//! println!("Estimated chance to win if we switch: {}", +//! switch_wins as f32 / total_switches as f32); +//! println!("Estimated chance to win if we don't: {}", +//! keep_wins as f32 / total_keeps as f32); +//! } +//! ``` + +#![doc(html_logo_url = "https://www.rust-lang.org/logos/rust-logo-128x128-blk.png", + html_favicon_url = "https://www.rust-lang.org/favicon.ico", + html_root_url = "https://docs.rs/rand/0.4")] + +#![deny(missing_debug_implementations)] + +#![cfg_attr(not(feature="std"), no_std)] +#![cfg_attr(all(feature="alloc", not(feature="std")), feature(alloc))] +#![cfg_attr(feature = "i128_support", feature(i128_type, i128))] + +#[cfg(feature="std")] extern crate std as core; +#[cfg(all(feature = "alloc", not(feature="std")))] extern crate alloc; + +use core::marker; +use core::mem; +#[cfg(feature="std")] use std::cell::RefCell; +#[cfg(feature="std")] use std::io; +#[cfg(feature="std")] use std::rc::Rc; + +// external rngs +pub use jitter::JitterRng; +#[cfg(feature="std")] pub use os::OsRng; + +// pseudo rngs +pub use isaac::{IsaacRng, Isaac64Rng}; +pub use chacha::ChaChaRng; +pub use prng::XorShiftRng; + +// local use declarations +#[cfg(target_pointer_width = "32")] +use prng::IsaacRng as IsaacWordRng; +#[cfg(target_pointer_width = "64")] +use prng::Isaac64Rng as IsaacWordRng; + +use distributions::{Range, IndependentSample}; +use distributions::range::SampleRange; + +// public modules +pub mod distributions; +pub mod jitter; +#[cfg(feature="std")] pub mod os; +#[cfg(feature="std")] pub mod read; +pub mod reseeding; +#[cfg(any(feature="std", feature = "alloc"))] pub mod seq; + +// These tiny modules are here to avoid API breakage, probably only temporarily +pub mod chacha { + //! The ChaCha random number generator. + pub use prng::ChaChaRng; +} +pub mod isaac { + //! The ISAAC random number generator. + pub use prng::{IsaacRng, Isaac64Rng}; +} + +// private modules +mod rand_impls; +mod prng; + + +/// A type that can be randomly generated using an `Rng`. +/// +/// ## Built-in Implementations +/// +/// This crate implements `Rand` for various primitive types. Assuming the +/// provided `Rng` is well-behaved, these implementations generate values with +/// the following ranges and distributions: +/// +/// * Integers (`i32`, `u32`, `isize`, `usize`, etc.): Uniformly distributed +/// over all values of the type. +/// * `char`: Uniformly distributed over all Unicode scalar values, i.e. all +/// code points in the range `0...0x10_FFFF`, except for the range +/// `0xD800...0xDFFF` (the surrogate code points). This includes +/// unassigned/reserved code points. +/// * `bool`: Generates `false` or `true`, each with probability 0.5. +/// * Floating point types (`f32` and `f64`): Uniformly distributed in the +/// half-open range `[0, 1)`. (The [`Open01`], [`Closed01`], [`Exp1`], and +/// [`StandardNormal`] wrapper types produce floating point numbers with +/// alternative ranges or distributions.) +/// +/// [`Open01`]: struct.Open01.html +/// [`Closed01`]: struct.Closed01.html +/// [`Exp1`]: distributions/exponential/struct.Exp1.html +/// [`StandardNormal`]: distributions/normal/struct.StandardNormal.html +/// +/// The following aggregate types also implement `Rand` as long as their +/// component types implement it: +/// +/// * Tuples and arrays: Each element of the tuple or array is generated +/// independently, using its own `Rand` implementation. +/// * `Option<T>`: Returns `None` with probability 0.5; otherwise generates a +/// random `T` and returns `Some(T)`. +pub trait Rand : Sized { + /// Generates a random instance of this type using the specified source of + /// randomness. + fn rand<R: Rng>(rng: &mut R) -> Self; +} + +/// A random number generator. +pub trait Rng { + /// Return the next random u32. + /// + /// This rarely needs to be called directly, prefer `r.gen()` to + /// `r.next_u32()`. + // FIXME #rust-lang/rfcs#628: Should be implemented in terms of next_u64 + fn next_u32(&mut self) -> u32; + + /// Return the next random u64. + /// + /// By default this is implemented in terms of `next_u32`. An + /// implementation of this trait must provide at least one of + /// these two methods. Similarly to `next_u32`, this rarely needs + /// to be called directly, prefer `r.gen()` to `r.next_u64()`. + fn next_u64(&mut self) -> u64 { + ((self.next_u32() as u64) << 32) | (self.next_u32() as u64) + } + + /// Return the next random f32 selected from the half-open + /// interval `[0, 1)`. + /// + /// This uses a technique described by Saito and Matsumoto at + /// MCQMC'08. Given that the IEEE floating point numbers are + /// uniformly distributed over [1,2), we generate a number in + /// this range and then offset it onto the range [0,1). Our + /// choice of bits (masking v. shifting) is arbitrary and + /// should be immaterial for high quality generators. For low + /// quality generators (ex. LCG), prefer bitshifting due to + /// correlation between sequential low order bits. + /// + /// See: + /// A PRNG specialized in double precision floating point numbers using + /// an affine transition + /// + /// * <http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/ARTICLES/dSFMT.pdf> + /// * <http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/SFMT/dSFMT-slide-e.pdf> + /// + /// By default this is implemented in terms of `next_u32`, but a + /// random number generator which can generate numbers satisfying + /// the requirements directly can overload this for performance. + /// It is required that the return value lies in `[0, 1)`. + /// + /// See `Closed01` for the closed interval `[0,1]`, and + /// `Open01` for the open interval `(0,1)`. + fn next_f32(&mut self) -> f32 { + const UPPER_MASK: u32 = 0x3F800000; + const LOWER_MASK: u32 = 0x7FFFFF; + let tmp = UPPER_MASK | (self.next_u32() & LOWER_MASK); + let result: f32 = unsafe { mem::transmute(tmp) }; + result - 1.0 + } + + /// Return the next random f64 selected from the half-open + /// interval `[0, 1)`. + /// + /// By default this is implemented in terms of `next_u64`, but a + /// random number generator which can generate numbers satisfying + /// the requirements directly can overload this for performance. + /// It is required that the return value lies in `[0, 1)`. + /// + /// See `Closed01` for the closed interval `[0,1]`, and + /// `Open01` for the open interval `(0,1)`. + fn next_f64(&mut self) -> f64 { + const UPPER_MASK: u64 = 0x3FF0000000000000; + const LOWER_MASK: u64 = 0xFFFFFFFFFFFFF; + let tmp = UPPER_MASK | (self.next_u64() & LOWER_MASK); + let result: f64 = unsafe { mem::transmute(tmp) }; + result - 1.0 + } + + /// Fill `dest` with random data. + /// + /// This has a default implementation in terms of `next_u64` and + /// `next_u32`, but should be overridden by implementations that + /// offer a more efficient solution than just calling those + /// methods repeatedly. + /// + /// This method does *not* have a requirement to bear any fixed + /// relationship to the other methods, for example, it does *not* + /// have to result in the same output as progressively filling + /// `dest` with `self.gen::<u8>()`, and any such behaviour should + /// not be relied upon. + /// + /// This method should guarantee that `dest` is entirely filled + /// with new data, and may panic if this is impossible + /// (e.g. reading past the end of a file that is being used as the + /// source of randomness). + /// + /// # Example + /// + /// ```rust + /// use rand::{thread_rng, Rng}; + /// + /// let mut v = [0u8; 13579]; + /// thread_rng().fill_bytes(&mut v); + /// println!("{:?}", &v[..]); + /// ``` + fn fill_bytes(&mut self, dest: &mut [u8]) { + // this could, in theory, be done by transmuting dest to a + // [u64], but this is (1) likely to be undefined behaviour for + // LLVM, (2) has to be very careful about alignment concerns, + // (3) adds more `unsafe` that needs to be checked, (4) + // probably doesn't give much performance gain if + // optimisations are on. + let mut count = 0; + let mut num = 0; + for byte in dest.iter_mut() { + if count == 0 { + // we could micro-optimise here by generating a u32 if + // we only need a few more bytes to fill the vector + // (i.e. at most 4). + num = self.next_u64(); + count = 8; + } + + *byte = (num & 0xff) as u8; + num >>= 8; + count -= 1; + } + } + + /// Return a random value of a `Rand` type. + /// + /// # Example + /// + /// ```rust + /// use rand::{thread_rng, Rng}; + /// + /// let mut rng = thread_rng(); + /// let x: u32 = rng.gen(); + /// println!("{}", x); + /// println!("{:?}", rng.gen::<(f64, bool)>()); + /// ``` + #[inline(always)] + fn gen<T: Rand>(&mut self) -> T where Self: Sized { + Rand::rand(self) + } + + /// Return an iterator that will yield an infinite number of randomly + /// generated items. + /// + /// # Example + /// + /// ``` + /// use rand::{thread_rng, Rng}; + /// + /// let mut rng = thread_rng(); + /// let x = rng.gen_iter::<u32>().take(10).collect::<Vec<u32>>(); + /// println!("{:?}", x); + /// println!("{:?}", rng.gen_iter::<(f64, bool)>().take(5) + /// .collect::<Vec<(f64, bool)>>()); + /// ``` + fn gen_iter<'a, T: Rand>(&'a mut self) -> Generator<'a, T, Self> where Self: Sized { + Generator { rng: self, _marker: marker::PhantomData } + } + + /// Generate a random value in the range [`low`, `high`). + /// + /// This is a convenience wrapper around + /// `distributions::Range`. If this function will be called + /// repeatedly with the same arguments, one should use `Range`, as + /// that will amortize the computations that allow for perfect + /// uniformity, as they only happen on initialization. + /// + /// # Panics + /// + /// Panics if `low >= high`. + /// + /// # Example + /// + /// ```rust + /// use rand::{thread_rng, Rng}; + /// + /// let mut rng = thread_rng(); + /// let n: u32 = rng.gen_range(0, 10); + /// println!("{}", n); + /// let m: f64 = rng.gen_range(-40.0f64, 1.3e5f64); + /// println!("{}", m); + /// ``` + fn gen_range<T: PartialOrd + SampleRange>(&mut self, low: T, high: T) -> T where Self: Sized { + assert!(low < high, "Rng.gen_range called with low >= high"); + Range::new(low, high).ind_sample(self) + } + + /// Return a bool with a 1 in n chance of true + /// + /// # Example + /// + /// ```rust + /// use rand::{thread_rng, Rng}; + /// + /// let mut rng = thread_rng(); + /// println!("{}", rng.gen_weighted_bool(3)); + /// ``` + fn gen_weighted_bool(&mut self, n: u32) -> bool where Self: Sized { + n <= 1 || self.gen_range(0, n) == 0 + } + + /// Return an iterator of random characters from the set A-Z,a-z,0-9. + /// + /// # Example + /// + /// ```rust + /// use rand::{thread_rng, Rng}; + /// + /// let s: String = thread_rng().gen_ascii_chars().take(10).collect(); + /// println!("{}", s); + /// ``` + fn gen_ascii_chars<'a>(&'a mut self) -> AsciiGenerator<'a, Self> where Self: Sized { + AsciiGenerator { rng: self } + } + + /// Return a random element from `values`. + /// + /// Return `None` if `values` is empty. + /// + /// # Example + /// + /// ``` + /// use rand::{thread_rng, Rng}; + /// + /// let choices = [1, 2, 4, 8, 16, 32]; + /// let mut rng = thread_rng(); + /// println!("{:?}", rng.choose(&choices)); + /// assert_eq!(rng.choose(&choices[..0]), None); + /// ``` + fn choose<'a, T>(&mut self, values: &'a [T]) -> Option<&'a T> where Self: Sized { + if values.is_empty() { + None + } else { + Some(&values[self.gen_range(0, values.len())]) + } + } + + /// Return a mutable pointer to a random element from `values`. + /// + /// Return `None` if `values` is empty. + fn choose_mut<'a, T>(&mut self, values: &'a mut [T]) -> Option<&'a mut T> where Self: Sized { + if values.is_empty() { + None + } else { + let len = values.len(); + Some(&mut values[self.gen_range(0, len)]) + } + } + + /// Shuffle a mutable slice in place. + /// + /// This applies Durstenfeld's algorithm for the [Fisher–Yates shuffle](https://en.wikipedia.org/wiki/Fisher%E2%80%93Yates_shuffle#The_modern_algorithm) + /// which produces an unbiased permutation. + /// + /// # Example + /// + /// ```rust + /// use rand::{thread_rng, Rng}; + /// + /// let mut rng = thread_rng(); + /// let mut y = [1, 2, 3]; + /// rng.shuffle(&mut y); + /// println!("{:?}", y); + /// rng.shuffle(&mut y); + /// println!("{:?}", y); + /// ``` + fn shuffle<T>(&mut self, values: &mut [T]) where Self: Sized { + let mut i = values.len(); + while i >= 2 { + // invariant: elements with index >= i have been locked in place. + i -= 1; + // lock element i in place. + values.swap(i, self.gen_range(0, i + 1)); + } + } +} + +impl<'a, R: ?Sized> Rng for &'a mut R where R: Rng { + fn next_u32(&mut self) -> u32 { + (**self).next_u32() + } + + fn next_u64(&mut self) -> u64 { + (**self).next_u64() + } + + fn next_f32(&mut self) -> f32 { + (**self).next_f32() + } + + fn next_f64(&mut self) -> f64 { + (**self).next_f64() + } + + fn fill_bytes(&mut self, dest: &mut [u8]) { + (**self).fill_bytes(dest) + } +} + +#[cfg(feature="std")] +impl<R: ?Sized> Rng for Box<R> where R: Rng { + fn next_u32(&mut self) -> u32 { + (**self).next_u32() + } + + fn next_u64(&mut self) -> u64 { + (**self).next_u64() + } + + fn next_f32(&mut self) -> f32 { + (**self).next_f32() + } + + fn next_f64(&mut self) -> f64 { + (**self).next_f64() + } + + fn fill_bytes(&mut self, dest: &mut [u8]) { + (**self).fill_bytes(dest) + } +} + +/// Iterator which will generate a stream of random items. +/// +/// This iterator is created via the [`gen_iter`] method on [`Rng`]. +/// +/// [`gen_iter`]: trait.Rng.html#method.gen_iter +/// [`Rng`]: trait.Rng.html +#[derive(Debug)] +pub struct Generator<'a, T, R:'a> { + rng: &'a mut R, + _marker: marker::PhantomData<fn() -> T>, +} + +impl<'a, T: Rand, R: Rng> Iterator for Generator<'a, T, R> { + type Item = T; + + fn next(&mut self) -> Option<T> { + Some(self.rng.gen()) + } +} + +/// Iterator which will continuously generate random ascii characters. +/// +/// This iterator is created via the [`gen_ascii_chars`] method on [`Rng`]. +/// +/// [`gen_ascii_chars`]: trait.Rng.html#method.gen_ascii_chars +/// [`Rng`]: trait.Rng.html +#[derive(Debug)] +pub struct AsciiGenerator<'a, R:'a> { + rng: &'a mut R, +} + +impl<'a, R: Rng> Iterator for AsciiGenerator<'a, R> { + type Item = char; + + fn next(&mut self) -> Option<char> { + const GEN_ASCII_STR_CHARSET: &'static [u8] = + b"ABCDEFGHIJKLMNOPQRSTUVWXYZ\ + abcdefghijklmnopqrstuvwxyz\ + 0123456789"; + Some(*self.rng.choose(GEN_ASCII_STR_CHARSET).unwrap() as char) + } +} + +/// A random number generator that can be explicitly seeded to produce +/// the same stream of randomness multiple times. +pub trait SeedableRng<Seed>: Rng { + /// Reseed an RNG with the given seed. + /// + /// # Example + /// + /// ```rust + /// use rand::{Rng, SeedableRng, StdRng}; + /// + /// let seed: &[_] = &[1, 2, 3, 4]; + /// let mut rng: StdRng = SeedableRng::from_seed(seed); + /// println!("{}", rng.gen::<f64>()); + /// rng.reseed(&[5, 6, 7, 8]); + /// println!("{}", rng.gen::<f64>()); + /// ``` + fn reseed(&mut self, Seed); + + /// Create a new RNG with the given seed. + /// + /// # Example + /// + /// ```rust + /// use rand::{Rng, SeedableRng, StdRng}; + /// + /// let seed: &[_] = &[1, 2, 3, 4]; + /// let mut rng: StdRng = SeedableRng::from_seed(seed); + /// println!("{}", rng.gen::<f64>()); + /// ``` + fn from_seed(seed: Seed) -> Self; +} + +/// A wrapper for generating floating point numbers uniformly in the +/// open interval `(0,1)` (not including either endpoint). +/// +/// Use `Closed01` for the closed interval `[0,1]`, and the default +/// `Rand` implementation for `f32` and `f64` for the half-open +/// `[0,1)`. +/// +/// # Example +/// ```rust +/// use rand::{random, Open01}; +/// +/// let Open01(val) = random::<Open01<f32>>(); +/// println!("f32 from (0,1): {}", val); +/// ``` +#[derive(Debug)] +pub struct Open01<F>(pub F); + +/// A wrapper for generating floating point numbers uniformly in the +/// closed interval `[0,1]` (including both endpoints). +/// +/// Use `Open01` for the closed interval `(0,1)`, and the default +/// `Rand` implementation of `f32` and `f64` for the half-open +/// `[0,1)`. +/// +/// # Example +/// +/// ```rust +/// use rand::{random, Closed01}; +/// +/// let Closed01(val) = random::<Closed01<f32>>(); +/// println!("f32 from [0,1]: {}", val); +/// ``` +#[derive(Debug)] +pub struct Closed01<F>(pub F); + +/// The standard RNG. This is designed to be efficient on the current +/// platform. +#[derive(Copy, Clone, Debug)] +pub struct StdRng { + rng: IsaacWordRng, +} + +impl StdRng { + /// Create a randomly seeded instance of `StdRng`. + /// + /// This is a very expensive operation as it has to read + /// randomness from the operating system and use this in an + /// expensive seeding operation. If one is only generating a small + /// number of random numbers, or doesn't need the utmost speed for + /// generating each number, `thread_rng` and/or `random` may be more + /// appropriate. + /// + /// Reading the randomness from the OS may fail, and any error is + /// propagated via the `io::Result` return value. + #[cfg(feature="std")] + pub fn new() -> io::Result<StdRng> { + match OsRng::new() { + Ok(mut r) => Ok(StdRng { rng: r.gen() }), + Err(e1) => { + match JitterRng::new() { + Ok(mut r) => Ok(StdRng { rng: r.gen() }), + Err(_) => { + Err(e1) + } + } + } + } + } +} + +impl Rng for StdRng { + #[inline] + fn next_u32(&mut self) -> u32 { + self.rng.next_u32() + } + + #[inline] + fn next_u64(&mut self) -> u64 { + self.rng.next_u64() + } +} + +impl<'a> SeedableRng<&'a [usize]> for StdRng { + fn reseed(&mut self, seed: &'a [usize]) { + // the internal RNG can just be seeded from the above + // randomness. + self.rng.reseed(unsafe {mem::transmute(seed)}) + } + + fn from_seed(seed: &'a [usize]) -> StdRng { + StdRng { rng: SeedableRng::from_seed(unsafe {mem::transmute(seed)}) } + } +} + +/// Create a weak random number generator with a default algorithm and seed. +/// +/// It returns the fastest `Rng` algorithm currently available in Rust without +/// consideration for cryptography or security. If you require a specifically +/// seeded `Rng` for consistency over time you should pick one algorithm and +/// create the `Rng` yourself. +/// +/// This will seed the generator with randomness from thread_rng. +#[cfg(feature="std")] +pub fn weak_rng() -> XorShiftRng { + thread_rng().gen() +} + +/// Controls how the thread-local RNG is reseeded. +#[cfg(feature="std")] +#[derive(Debug)] +struct ThreadRngReseeder; + +#[cfg(feature="std")] +impl reseeding::Reseeder<StdRng> for ThreadRngReseeder { + fn reseed(&mut self, rng: &mut StdRng) { + match StdRng::new() { + Ok(r) => *rng = r, + Err(e) => panic!("No entropy available: {}", e), + } + } +} +#[cfg(feature="std")] +const THREAD_RNG_RESEED_THRESHOLD: u64 = 32_768; +#[cfg(feature="std")] +type ThreadRngInner = reseeding::ReseedingRng<StdRng, ThreadRngReseeder>; + +/// The thread-local RNG. +#[cfg(feature="std")] +#[derive(Clone, Debug)] +pub struct ThreadRng { + rng: Rc<RefCell<ThreadRngInner>>, +} + +/// Retrieve the lazily-initialized thread-local random number +/// generator, seeded by the system. Intended to be used in method +/// chaining style, e.g. `thread_rng().gen::<i32>()`. +/// +/// After generating a certain amount of randomness, the RNG will reseed itself +/// from the operating system or, if the operating system RNG returns an error, +/// a seed based on the current system time. +/// +/// The internal RNG used is platform and architecture dependent, even +/// if the operating system random number generator is rigged to give +/// the same sequence always. If absolute consistency is required, +/// explicitly select an RNG, e.g. `IsaacRng` or `Isaac64Rng`. +#[cfg(feature="std")] +pub fn thread_rng() -> ThreadRng { + // used to make space in TLS for a random number generator + thread_local!(static THREAD_RNG_KEY: Rc<RefCell<ThreadRngInner>> = { + let r = match StdRng::new() { + Ok(r) => r, + Err(e) => panic!("No entropy available: {}", e), + }; + let rng = reseeding::ReseedingRng::new(r, + THREAD_RNG_RESEED_THRESHOLD, + ThreadRngReseeder); + Rc::new(RefCell::new(rng)) + }); + + ThreadRng { rng: THREAD_RNG_KEY.with(|t| t.clone()) } +} + +#[cfg(feature="std")] +impl Rng for ThreadRng { + fn next_u32(&mut self) -> u32 { + self.rng.borrow_mut().next_u32() + } + + fn next_u64(&mut self) -> u64 { + self.rng.borrow_mut().next_u64() + } + + #[inline] + fn fill_bytes(&mut self, bytes: &mut [u8]) { + self.rng.borrow_mut().fill_bytes(bytes) + } +} + +/// Generates a random value using the thread-local random number generator. +/// +/// `random()` can generate various types of random things, and so may require +/// type hinting to generate the specific type you want. +/// +/// This function uses the thread local random number generator. This means +/// that if you're calling `random()` in a loop, caching the generator can +/// increase performance. An example is shown below. +/// +/// # Examples +/// +/// ``` +/// let x = rand::random::<u8>(); +/// println!("{}", x); +/// +/// let y = rand::random::<f64>(); +/// println!("{}", y); +/// +/// if rand::random() { // generates a boolean +/// println!("Better lucky than good!"); +/// } +/// ``` +/// +/// Caching the thread local random number generator: +/// +/// ``` +/// use rand::Rng; +/// +/// let mut v = vec![1, 2, 3]; +/// +/// for x in v.iter_mut() { +/// *x = rand::random() +/// } +/// +/// // can be made faster by caching thread_rng +/// +/// let mut rng = rand::thread_rng(); +/// +/// for x in v.iter_mut() { +/// *x = rng.gen(); +/// } +/// ``` +#[cfg(feature="std")] +#[inline] +pub fn random<T: Rand>() -> T { + thread_rng().gen() +} + +/// DEPRECATED: use `seq::sample_iter` instead. +/// +/// Randomly sample up to `amount` elements from a finite iterator. +/// The order of elements in the sample is not random. +/// +/// # Example +/// +/// ```rust +/// use rand::{thread_rng, sample}; +/// +/// let mut rng = thread_rng(); +/// let sample = sample(&mut rng, 1..100, 5); +/// println!("{:?}", sample); +/// ``` +#[cfg(feature="std")] +#[inline(always)] +#[deprecated(since="0.4.0", note="renamed to seq::sample_iter")] +pub fn sample<T, I, R>(rng: &mut R, iterable: I, amount: usize) -> Vec<T> + where I: IntoIterator<Item=T>, + R: Rng, +{ + // the legacy sample didn't care whether amount was met + seq::sample_iter(rng, iterable, amount) + .unwrap_or_else(|e| e) +} + +#[cfg(test)] +mod test { + use super::{Rng, thread_rng, random, SeedableRng, StdRng, weak_rng}; + use std::iter::repeat; + + pub struct MyRng<R> { inner: R } + + impl<R: Rng> Rng for MyRng<R> { + fn next_u32(&mut self) -> u32 { + fn next<T: Rng>(t: &mut T) -> u32 { + t.next_u32() + } + next(&mut self.inner) + } + } + + pub fn rng() -> MyRng<::ThreadRng> { + MyRng { inner: ::thread_rng() } + } + + struct ConstRng { i: u64 } + impl Rng for ConstRng { + fn next_u32(&mut self) -> u32 { self.i as u32 } + fn next_u64(&mut self) -> u64 { self.i } + + // no fill_bytes on purpose + } + + pub fn iter_eq<I, J>(i: I, j: J) -> bool + where I: IntoIterator, + J: IntoIterator<Item=I::Item>, + I::Item: Eq + { + // make sure the iterators have equal length + let mut i = i.into_iter(); + let mut j = j.into_iter(); + loop { + match (i.next(), j.next()) { + (Some(ref ei), Some(ref ej)) if ei == ej => { } + (None, None) => return true, + _ => return false, + } + } + } + + #[test] + fn test_fill_bytes_default() { + let mut r = ConstRng { i: 0x11_22_33_44_55_66_77_88 }; + + // check every remainder mod 8, both in small and big vectors. + let lengths = [0, 1, 2, 3, 4, 5, 6, 7, + 80, 81, 82, 83, 84, 85, 86, 87]; + for &n in lengths.iter() { + let mut v = repeat(0u8).take(n).collect::<Vec<_>>(); + r.fill_bytes(&mut v); + + // use this to get nicer error messages. + for (i, &byte) in v.iter().enumerate() { + if byte == 0 { + panic!("byte {} of {} is zero", i, n) + } + } + } + } + + #[test] + fn test_gen_range() { + let mut r = thread_rng(); + for _ in 0..1000 { + let a = r.gen_range(-3, 42); + assert!(a >= -3 && a < 42); + assert_eq!(r.gen_range(0, 1), 0); + assert_eq!(r.gen_range(-12, -11), -12); + } + + for _ in 0..1000 { + let a = r.gen_range(10, 42); + assert!(a >= 10 && a < 42); + assert_eq!(r.gen_range(0, 1), 0); + assert_eq!(r.gen_range(3_000_000, 3_000_001), 3_000_000); + } + + } + + #[test] + #[should_panic] + fn test_gen_range_panic_int() { + let mut r = thread_rng(); + r.gen_range(5, -2); + } + + #[test] + #[should_panic] + fn test_gen_range_panic_usize() { + let mut r = thread_rng(); + r.gen_range(5, 2); + } + + #[test] + fn test_gen_weighted_bool() { + let mut r = thread_rng(); + assert_eq!(r.gen_weighted_bool(0), true); + assert_eq!(r.gen_weighted_bool(1), true); + } + + #[test] + fn test_gen_ascii_str() { + let mut r = thread_rng(); + assert_eq!(r.gen_ascii_chars().take(0).count(), 0); + assert_eq!(r.gen_ascii_chars().take(10).count(), 10); + assert_eq!(r.gen_ascii_chars().take(16).count(), 16); + } + + #[test] + fn test_gen_vec() { + let mut r = thread_rng(); + assert_eq!(r.gen_iter::<u8>().take(0).count(), 0); + assert_eq!(r.gen_iter::<u8>().take(10).count(), 10); + assert_eq!(r.gen_iter::<f64>().take(16).count(), 16); + } + + #[test] + fn test_choose() { + let mut r = thread_rng(); + assert_eq!(r.choose(&[1, 1, 1]).map(|&x|x), Some(1)); + + let v: &[isize] = &[]; + assert_eq!(r.choose(v), None); + } + + #[test] + fn test_shuffle() { + let mut r = thread_rng(); + let empty: &mut [isize] = &mut []; + r.shuffle(empty); + let mut one = [1]; + r.shuffle(&mut one); + let b: &[_] = &[1]; + assert_eq!(one, b); + + let mut two = [1, 2]; + r.shuffle(&mut two); + assert!(two == [1, 2] || two == [2, 1]); + + let mut x = [1, 1, 1]; + r.shuffle(&mut x); + let b: &[_] = &[1, 1, 1]; + assert_eq!(x, b); + } + + #[test] + fn test_thread_rng() { + let mut r = thread_rng(); + r.gen::<i32>(); + let mut v = [1, 1, 1]; + r.shuffle(&mut v); + let b: &[_] = &[1, 1, 1]; + assert_eq!(v, b); + assert_eq!(r.gen_range(0, 1), 0); + } + + #[test] + fn test_rng_trait_object() { + let mut rng = thread_rng(); + { + let mut r = &mut rng as &mut Rng; + r.next_u32(); + (&mut r).gen::<i32>(); + let mut v = [1, 1, 1]; + (&mut r).shuffle(&mut v); + let b: &[_] = &[1, 1, 1]; + assert_eq!(v, b); + assert_eq!((&mut r).gen_range(0, 1), 0); + } + { + let mut r = Box::new(rng) as Box<Rng>; + r.next_u32(); + r.gen::<i32>(); + let mut v = [1, 1, 1]; + r.shuffle(&mut v); + let b: &[_] = &[1, 1, 1]; + assert_eq!(v, b); + assert_eq!(r.gen_range(0, 1), 0); + } + } + + #[test] + fn test_random() { + // not sure how to test this aside from just getting some values + let _n : usize = random(); + let _f : f32 = random(); + let _o : Option<Option<i8>> = random(); + let _many : ((), + (usize, + isize, + Option<(u32, (bool,))>), + (u8, i8, u16, i16, u32, i32, u64, i64), + (f32, (f64, (f64,)))) = random(); + } + + #[test] + fn test_std_rng_seeded() { + let s = thread_rng().gen_iter::<usize>().take(256).collect::<Vec<usize>>(); + let mut ra: StdRng = SeedableRng::from_seed(&s[..]); + let mut rb: StdRng = SeedableRng::from_seed(&s[..]); + assert!(iter_eq(ra.gen_ascii_chars().take(100), + rb.gen_ascii_chars().take(100))); + } + + #[test] + fn test_std_rng_reseed() { + let s = thread_rng().gen_iter::<usize>().take(256).collect::<Vec<usize>>(); + let mut r: StdRng = SeedableRng::from_seed(&s[..]); + let string1 = r.gen_ascii_chars().take(100).collect::<String>(); + + r.reseed(&s); + + let string2 = r.gen_ascii_chars().take(100).collect::<String>(); + assert_eq!(string1, string2); + } + + #[test] + fn test_weak_rng() { + let s = weak_rng().gen_iter::<usize>().take(256).collect::<Vec<usize>>(); + let mut ra: StdRng = SeedableRng::from_seed(&s[..]); + let mut rb: StdRng = SeedableRng::from_seed(&s[..]); + assert!(iter_eq(ra.gen_ascii_chars().take(100), + rb.gen_ascii_chars().take(100))); + } +} diff --git a/rand/src/os.rs b/rand/src/os.rs new file mode 100644 index 0000000..10022fb --- /dev/null +++ b/rand/src/os.rs @@ -0,0 +1,617 @@ +// Copyright 2013-2015 The Rust Project Developers. See the COPYRIGHT +// file at the top-level directory of this distribution and at +// http://rust-lang.org/COPYRIGHT. +// +// Licensed under the Apache License, Version 2.0 <LICENSE-APACHE or +// http://www.apache.org/licenses/LICENSE-2.0> or the MIT license +// <LICENSE-MIT or http://opensource.org/licenses/MIT>, at your +// option. This file may not be copied, modified, or distributed +// except according to those terms. + +//! Interfaces to the operating system provided random number +//! generators. + +use std::{io, mem, fmt}; +use Rng; + +/// A random number generator that retrieves randomness straight from +/// the operating system. Platform sources: +/// +/// - Unix-like systems (Linux, Android, Mac OSX): read directly from +/// `/dev/urandom`, or from `getrandom(2)` system call if available. +/// - OpenBSD: calls `getentropy(2)` +/// - FreeBSD: uses the `kern.arandom` `sysctl(2)` mib +/// - Windows: calls `RtlGenRandom`, exported from `advapi32.dll` as +/// `SystemFunction036`. +/// - iOS: calls SecRandomCopyBytes as /dev/(u)random is sandboxed. +/// - PNaCl: calls into the `nacl-irt-random-0.1` IRT interface. +/// +/// This usually does not block. On some systems (e.g. FreeBSD, OpenBSD, +/// Max OS X, and modern Linux) this may block very early in the init +/// process, if the CSPRNG has not been seeded yet.[1] +/// +/// [1] See <https://www.python.org/dev/peps/pep-0524/> for a more +/// in-depth discussion. +pub struct OsRng(imp::OsRng); + +impl OsRng { + /// Create a new `OsRng`. + pub fn new() -> io::Result<OsRng> { + imp::OsRng::new().map(OsRng) + } +} + +impl Rng for OsRng { + fn next_u32(&mut self) -> u32 { self.0.next_u32() } + fn next_u64(&mut self) -> u64 { self.0.next_u64() } + fn fill_bytes(&mut self, v: &mut [u8]) { self.0.fill_bytes(v) } +} + +impl fmt::Debug for OsRng { + fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result { + write!(f, "OsRng {{}}") + } +} + +fn next_u32(fill_buf: &mut FnMut(&mut [u8])) -> u32 { + let mut buf: [u8; 4] = [0; 4]; + fill_buf(&mut buf); + unsafe { mem::transmute::<[u8; 4], u32>(buf) } +} + +fn next_u64(fill_buf: &mut FnMut(&mut [u8])) -> u64 { + let mut buf: [u8; 8] = [0; 8]; + fill_buf(&mut buf); + unsafe { mem::transmute::<[u8; 8], u64>(buf) } +} + +#[cfg(all(unix, not(target_os = "ios"), + not(target_os = "nacl"), + not(target_os = "freebsd"), + not(target_os = "fuchsia"), + not(target_os = "openbsd"), + not(target_os = "redox")))] +mod imp { + extern crate libc; + + use super::{next_u32, next_u64}; + use self::OsRngInner::*; + + use std::io; + use std::fs::File; + use Rng; + use read::ReadRng; + + #[cfg(all(target_os = "linux", + any(target_arch = "x86_64", + target_arch = "x86", + target_arch = "arm", + target_arch = "aarch64", + target_arch = "powerpc")))] + fn getrandom(buf: &mut [u8]) -> libc::c_long { + extern "C" { + fn syscall(number: libc::c_long, ...) -> libc::c_long; + } + + #[cfg(target_arch = "x86_64")] + const NR_GETRANDOM: libc::c_long = 318; + #[cfg(target_arch = "x86")] + const NR_GETRANDOM: libc::c_long = 355; + #[cfg(target_arch = "arm")] + const NR_GETRANDOM: libc::c_long = 384; + #[cfg(target_arch = "aarch64")] + const NR_GETRANDOM: libc::c_long = 278; + #[cfg(target_arch = "powerpc")] + const NR_GETRANDOM: libc::c_long = 359; + + unsafe { + syscall(NR_GETRANDOM, buf.as_mut_ptr(), buf.len(), 0) + } + } + + #[cfg(not(all(target_os = "linux", + any(target_arch = "x86_64", + target_arch = "x86", + target_arch = "arm", + target_arch = "aarch64", + target_arch = "powerpc"))))] + fn getrandom(_buf: &mut [u8]) -> libc::c_long { -1 } + + fn getrandom_fill_bytes(v: &mut [u8]) { + let mut read = 0; + let len = v.len(); + while read < len { + let result = getrandom(&mut v[read..]); + if result == -1 { + let err = io::Error::last_os_error(); + if err.kind() == io::ErrorKind::Interrupted { + continue + } else { + panic!("unexpected getrandom error: {}", err); + } + } else { + read += result as usize; + } + } + } + + #[cfg(all(target_os = "linux", + any(target_arch = "x86_64", + target_arch = "x86", + target_arch = "arm", + target_arch = "aarch64", + target_arch = "powerpc")))] + fn is_getrandom_available() -> bool { + use std::sync::atomic::{AtomicBool, ATOMIC_BOOL_INIT, Ordering}; + use std::sync::{Once, ONCE_INIT}; + + static CHECKER: Once = ONCE_INIT; + static AVAILABLE: AtomicBool = ATOMIC_BOOL_INIT; + + CHECKER.call_once(|| { + let mut buf: [u8; 0] = []; + let result = getrandom(&mut buf); + let available = if result == -1 { + let err = io::Error::last_os_error().raw_os_error(); + err != Some(libc::ENOSYS) + } else { + true + }; + AVAILABLE.store(available, Ordering::Relaxed); + }); + + AVAILABLE.load(Ordering::Relaxed) + } + + #[cfg(not(all(target_os = "linux", + any(target_arch = "x86_64", + target_arch = "x86", + target_arch = "arm", + target_arch = "aarch64", + target_arch = "powerpc"))))] + fn is_getrandom_available() -> bool { false } + + pub struct OsRng { + inner: OsRngInner, + } + + enum OsRngInner { + OsGetrandomRng, + OsReadRng(ReadRng<File>), + } + + impl OsRng { + pub fn new() -> io::Result<OsRng> { + if is_getrandom_available() { + return Ok(OsRng { inner: OsGetrandomRng }); + } + + let reader = try!(File::open("/dev/urandom")); + let reader_rng = ReadRng::new(reader); + + Ok(OsRng { inner: OsReadRng(reader_rng) }) + } + } + + impl Rng for OsRng { + fn next_u32(&mut self) -> u32 { + match self.inner { + OsGetrandomRng => next_u32(&mut getrandom_fill_bytes), + OsReadRng(ref mut rng) => rng.next_u32(), + } + } + fn next_u64(&mut self) -> u64 { + match self.inner { + OsGetrandomRng => next_u64(&mut getrandom_fill_bytes), + OsReadRng(ref mut rng) => rng.next_u64(), + } + } + fn fill_bytes(&mut self, v: &mut [u8]) { + match self.inner { + OsGetrandomRng => getrandom_fill_bytes(v), + OsReadRng(ref mut rng) => rng.fill_bytes(v) + } + } + } +} + +#[cfg(target_os = "ios")] +mod imp { + extern crate libc; + + use super::{next_u32, next_u64}; + + use std::io; + use Rng; + use self::libc::{c_int, size_t}; + + #[derive(Debug)] + pub struct OsRng; + + enum SecRandom {} + + #[allow(non_upper_case_globals)] + const kSecRandomDefault: *const SecRandom = 0 as *const SecRandom; + + #[link(name = "Security", kind = "framework")] + extern { + fn SecRandomCopyBytes(rnd: *const SecRandom, + count: size_t, bytes: *mut u8) -> c_int; + } + + impl OsRng { + pub fn new() -> io::Result<OsRng> { + Ok(OsRng) + } + } + + impl Rng for OsRng { + fn next_u32(&mut self) -> u32 { + next_u32(&mut |v| self.fill_bytes(v)) + } + fn next_u64(&mut self) -> u64 { + next_u64(&mut |v| self.fill_bytes(v)) + } + fn fill_bytes(&mut self, v: &mut [u8]) { + let ret = unsafe { + SecRandomCopyBytes(kSecRandomDefault, v.len() as size_t, v.as_mut_ptr()) + }; + if ret == -1 { + panic!("couldn't generate random bytes: {}", io::Error::last_os_error()); + } + } + } +} + +#[cfg(target_os = "freebsd")] +mod imp { + extern crate libc; + + use std::{io, ptr}; + use Rng; + + use super::{next_u32, next_u64}; + + #[derive(Debug)] + pub struct OsRng; + + impl OsRng { + pub fn new() -> io::Result<OsRng> { + Ok(OsRng) + } + } + + impl Rng for OsRng { + fn next_u32(&mut self) -> u32 { + next_u32(&mut |v| self.fill_bytes(v)) + } + fn next_u64(&mut self) -> u64 { + next_u64(&mut |v| self.fill_bytes(v)) + } + fn fill_bytes(&mut self, v: &mut [u8]) { + let mib = [libc::CTL_KERN, libc::KERN_ARND]; + // kern.arandom permits a maximum buffer size of 256 bytes + for s in v.chunks_mut(256) { + let mut s_len = s.len(); + let ret = unsafe { + libc::sysctl(mib.as_ptr(), mib.len() as libc::c_uint, + s.as_mut_ptr() as *mut _, &mut s_len, + ptr::null(), 0) + }; + if ret == -1 || s_len != s.len() { + panic!("kern.arandom sysctl failed! (returned {}, s.len() {}, oldlenp {})", + ret, s.len(), s_len); + } + } + } + } +} + +#[cfg(target_os = "openbsd")] +mod imp { + extern crate libc; + + use std::io; + use Rng; + + use super::{next_u32, next_u64}; + + #[derive(Debug)] + pub struct OsRng; + + impl OsRng { + pub fn new() -> io::Result<OsRng> { + Ok(OsRng) + } + } + + impl Rng for OsRng { + fn next_u32(&mut self) -> u32 { + next_u32(&mut |v| self.fill_bytes(v)) + } + fn next_u64(&mut self) -> u64 { + next_u64(&mut |v| self.fill_bytes(v)) + } + fn fill_bytes(&mut self, v: &mut [u8]) { + // getentropy(2) permits a maximum buffer size of 256 bytes + for s in v.chunks_mut(256) { + let ret = unsafe { + libc::getentropy(s.as_mut_ptr() as *mut libc::c_void, s.len()) + }; + if ret == -1 { + let err = io::Error::last_os_error(); + panic!("getentropy failed: {}", err); + } + } + } + } +} + +#[cfg(target_os = "redox")] +mod imp { + use std::io; + use std::fs::File; + use Rng; + use read::ReadRng; + + #[derive(Debug)] + pub struct OsRng { + inner: ReadRng<File>, + } + + impl OsRng { + pub fn new() -> io::Result<OsRng> { + let reader = try!(File::open("rand:")); + let reader_rng = ReadRng::new(reader); + + Ok(OsRng { inner: reader_rng }) + } + } + + impl Rng for OsRng { + fn next_u32(&mut self) -> u32 { + self.inner.next_u32() + } + fn next_u64(&mut self) -> u64 { + self.inner.next_u64() + } + fn fill_bytes(&mut self, v: &mut [u8]) { + self.inner.fill_bytes(v) + } + } +} + +#[cfg(target_os = "fuchsia")] +mod imp { + extern crate fuchsia_zircon; + + use std::io; + use Rng; + + use super::{next_u32, next_u64}; + + #[derive(Debug)] + pub struct OsRng; + + impl OsRng { + pub fn new() -> io::Result<OsRng> { + Ok(OsRng) + } + } + + impl Rng for OsRng { + fn next_u32(&mut self) -> u32 { + next_u32(&mut |v| self.fill_bytes(v)) + } + fn next_u64(&mut self) -> u64 { + next_u64(&mut |v| self.fill_bytes(v)) + } + fn fill_bytes(&mut self, v: &mut [u8]) { + for s in v.chunks_mut(fuchsia_zircon::sys::ZX_CPRNG_DRAW_MAX_LEN) { + let mut filled = 0; + while filled < s.len() { + match fuchsia_zircon::cprng_draw(&mut s[filled..]) { + Ok(actual) => filled += actual, + Err(e) => panic!("cprng_draw failed: {:?}", e), + }; + } + } + } + } +} + +#[cfg(windows)] +mod imp { + extern crate winapi; + + use std::io; + use Rng; + + use super::{next_u32, next_u64}; + + use self::winapi::shared::minwindef::ULONG; + use self::winapi::um::ntsecapi::RtlGenRandom; + use self::winapi::um::winnt::PVOID; + + #[derive(Debug)] + pub struct OsRng; + + impl OsRng { + pub fn new() -> io::Result<OsRng> { + Ok(OsRng) + } + } + + impl Rng for OsRng { + fn next_u32(&mut self) -> u32 { + next_u32(&mut |v| self.fill_bytes(v)) + } + fn next_u64(&mut self) -> u64 { + next_u64(&mut |v| self.fill_bytes(v)) + } + fn fill_bytes(&mut self, v: &mut [u8]) { + // RtlGenRandom takes an ULONG (u32) for the length so we need to + // split up the buffer. + for slice in v.chunks_mut(<ULONG>::max_value() as usize) { + let ret = unsafe { + RtlGenRandom(slice.as_mut_ptr() as PVOID, slice.len() as ULONG) + }; + if ret == 0 { + panic!("couldn't generate random bytes: {}", + io::Error::last_os_error()); + } + } + } + } +} + +#[cfg(target_os = "nacl")] +mod imp { + extern crate libc; + + use std::io; + use std::mem; + use Rng; + + use super::{next_u32, next_u64}; + + #[derive(Debug)] + pub struct OsRng(extern fn(dest: *mut libc::c_void, + bytes: libc::size_t, + read: *mut libc::size_t) -> libc::c_int); + + extern { + fn nacl_interface_query(name: *const libc::c_char, + table: *mut libc::c_void, + table_size: libc::size_t) -> libc::size_t; + } + + const INTERFACE: &'static [u8] = b"nacl-irt-random-0.1\0"; + + #[repr(C)] + struct NaClIRTRandom { + get_random_bytes: Option<extern fn(dest: *mut libc::c_void, + bytes: libc::size_t, + read: *mut libc::size_t) -> libc::c_int>, + } + + impl OsRng { + pub fn new() -> io::Result<OsRng> { + let mut iface = NaClIRTRandom { + get_random_bytes: None, + }; + let result = unsafe { + nacl_interface_query(INTERFACE.as_ptr() as *const _, + mem::transmute(&mut iface), + mem::size_of::<NaClIRTRandom>() as libc::size_t) + }; + if result != 0 { + assert!(iface.get_random_bytes.is_some()); + let result = OsRng(iface.get_random_bytes.take().unwrap()); + Ok(result) + } else { + let error = io::ErrorKind::NotFound; + let error = io::Error::new(error, "IRT random interface missing"); + Err(error) + } + } + } + + impl Rng for OsRng { + fn next_u32(&mut self) -> u32 { + next_u32(&mut |v| self.fill_bytes(v)) + } + fn next_u64(&mut self) -> u64 { + next_u64(&mut |v| self.fill_bytes(v)) + } + fn fill_bytes(&mut self, v: &mut [u8]) { + let mut read = 0; + loop { + let mut r: libc::size_t = 0; + let len = v.len(); + let error = (self.0)(v[read..].as_mut_ptr() as *mut _, + (len - read) as libc::size_t, + &mut r as *mut _); + assert!(error == 0, "`get_random_bytes` failed!"); + read += r as usize; + + if read >= v.len() { break; } + } + } + } +} + +#[cfg(all(target_arch = "wasm32", not(target_os = "emscripten")))] +mod imp { + use std::io; + use Rng; + + #[derive(Debug)] + pub struct OsRng; + + impl OsRng { + pub fn new() -> io::Result<OsRng> { + Err(io::Error::new(io::ErrorKind::Other, "Not supported")) + } + } + + impl Rng for OsRng { + fn next_u32(&mut self) -> u32 { + panic!("Not supported") + } + } +} + +#[cfg(test)] +mod test { + use std::sync::mpsc::channel; + use Rng; + use OsRng; + use std::thread; + + #[test] + fn test_os_rng() { + let mut r = OsRng::new().unwrap(); + + r.next_u32(); + r.next_u64(); + + let mut v = [0u8; 1000]; + r.fill_bytes(&mut v); + } + + #[test] + fn test_os_rng_tasks() { + + let mut txs = vec!(); + for _ in 0..20 { + let (tx, rx) = channel(); + txs.push(tx); + + thread::spawn(move|| { + // wait until all the tasks are ready to go. + rx.recv().unwrap(); + + // deschedule to attempt to interleave things as much + // as possible (XXX: is this a good test?) + let mut r = OsRng::new().unwrap(); + thread::yield_now(); + let mut v = [0u8; 1000]; + + for _ in 0..100 { + r.next_u32(); + thread::yield_now(); + r.next_u64(); + thread::yield_now(); + r.fill_bytes(&mut v); + thread::yield_now(); + } + }); + } + + // start all the tasks + for tx in txs.iter() { + tx.send(()).unwrap(); + } + } +} diff --git a/rand/src/prng/chacha.rs b/rand/src/prng/chacha.rs new file mode 100644 index 0000000..a73e8e7 --- /dev/null +++ b/rand/src/prng/chacha.rs @@ -0,0 +1,321 @@ +// Copyright 2014 The Rust Project Developers. See the COPYRIGHT +// file at the top-level directory of this distribution and at +// http://rust-lang.org/COPYRIGHT. +// +// Licensed under the Apache License, Version 2.0 <LICENSE-APACHE or +// http://www.apache.org/licenses/LICENSE-2.0> or the MIT license +// <LICENSE-MIT or http://opensource.org/licenses/MIT>, at your +// option. This file may not be copied, modified, or distributed +// except according to those terms. + +//! The ChaCha random number generator. + +use core::num::Wrapping as w; +use {Rng, SeedableRng, Rand}; + +#[allow(bad_style)] +type w32 = w<u32>; + +const KEY_WORDS : usize = 8; // 8 words for the 256-bit key +const STATE_WORDS : usize = 16; +const CHACHA_ROUNDS: u32 = 20; // Cryptographically secure from 8 upwards as of this writing + +/// A random number generator that uses the ChaCha20 algorithm [1]. +/// +/// The ChaCha algorithm is widely accepted as suitable for +/// cryptographic purposes, but this implementation has not been +/// verified as such. Prefer a generator like `OsRng` that defers to +/// the operating system for cases that need high security. +/// +/// [1]: D. J. Bernstein, [*ChaCha, a variant of +/// Salsa20*](http://cr.yp.to/chacha.html) +#[derive(Copy, Clone, Debug)] +pub struct ChaChaRng { + buffer: [w32; STATE_WORDS], // Internal buffer of output + state: [w32; STATE_WORDS], // Initial state + index: usize, // Index into state +} + +static EMPTY: ChaChaRng = ChaChaRng { + buffer: [w(0); STATE_WORDS], + state: [w(0); STATE_WORDS], + index: STATE_WORDS +}; + + +macro_rules! quarter_round{ + ($a: expr, $b: expr, $c: expr, $d: expr) => {{ + $a = $a + $b; $d = $d ^ $a; $d = w($d.0.rotate_left(16)); + $c = $c + $d; $b = $b ^ $c; $b = w($b.0.rotate_left(12)); + $a = $a + $b; $d = $d ^ $a; $d = w($d.0.rotate_left( 8)); + $c = $c + $d; $b = $b ^ $c; $b = w($b.0.rotate_left( 7)); + }} +} + +macro_rules! double_round{ + ($x: expr) => {{ + // Column round + quarter_round!($x[ 0], $x[ 4], $x[ 8], $x[12]); + quarter_round!($x[ 1], $x[ 5], $x[ 9], $x[13]); + quarter_round!($x[ 2], $x[ 6], $x[10], $x[14]); + quarter_round!($x[ 3], $x[ 7], $x[11], $x[15]); + // Diagonal round + quarter_round!($x[ 0], $x[ 5], $x[10], $x[15]); + quarter_round!($x[ 1], $x[ 6], $x[11], $x[12]); + quarter_round!($x[ 2], $x[ 7], $x[ 8], $x[13]); + quarter_round!($x[ 3], $x[ 4], $x[ 9], $x[14]); + }} +} + +#[inline] +fn core(output: &mut [w32; STATE_WORDS], input: &[w32; STATE_WORDS]) { + *output = *input; + + for _ in 0..CHACHA_ROUNDS / 2 { + double_round!(output); + } + + for i in 0..STATE_WORDS { + output[i] = output[i] + input[i]; + } +} + +impl ChaChaRng { + + /// Create an ChaCha random number generator using the default + /// fixed key of 8 zero words. + /// + /// # Examples + /// + /// ```rust + /// use rand::{Rng, ChaChaRng}; + /// + /// let mut ra = ChaChaRng::new_unseeded(); + /// println!("{:?}", ra.next_u32()); + /// println!("{:?}", ra.next_u32()); + /// ``` + /// + /// Since this equivalent to a RNG with a fixed seed, repeated executions + /// of an unseeded RNG will produce the same result. This code sample will + /// consistently produce: + /// + /// - 2917185654 + /// - 2419978656 + pub fn new_unseeded() -> ChaChaRng { + let mut rng = EMPTY; + rng.init(&[0; KEY_WORDS]); + rng + } + + /// Sets the internal 128-bit ChaCha counter to + /// a user-provided value. This permits jumping + /// arbitrarily ahead (or backwards) in the pseudorandom stream. + /// + /// Since the nonce words are used to extend the counter to 128 bits, + /// users wishing to obtain the conventional ChaCha pseudorandom stream + /// associated with a particular nonce can call this function with + /// arguments `0, desired_nonce`. + /// + /// # Examples + /// + /// ```rust + /// use rand::{Rng, ChaChaRng}; + /// + /// let mut ra = ChaChaRng::new_unseeded(); + /// ra.set_counter(0u64, 1234567890u64); + /// println!("{:?}", ra.next_u32()); + /// println!("{:?}", ra.next_u32()); + /// ``` + pub fn set_counter(&mut self, counter_low: u64, counter_high: u64) { + self.state[12] = w((counter_low >> 0) as u32); + self.state[13] = w((counter_low >> 32) as u32); + self.state[14] = w((counter_high >> 0) as u32); + self.state[15] = w((counter_high >> 32) as u32); + self.index = STATE_WORDS; // force recomputation + } + + /// Initializes `self.state` with the appropriate key and constants + /// + /// We deviate slightly from the ChaCha specification regarding + /// the nonce, which is used to extend the counter to 128 bits. + /// This is provably as strong as the original cipher, though, + /// since any distinguishing attack on our variant also works + /// against ChaCha with a chosen-nonce. See the XSalsa20 [1] + /// security proof for a more involved example of this. + /// + /// The modified word layout is: + /// ```text + /// constant constant constant constant + /// key key key key + /// key key key key + /// counter counter counter counter + /// ``` + /// [1]: Daniel J. Bernstein. [*Extending the Salsa20 + /// nonce.*](http://cr.yp.to/papers.html#xsalsa) + fn init(&mut self, key: &[u32; KEY_WORDS]) { + self.state[0] = w(0x61707865); + self.state[1] = w(0x3320646E); + self.state[2] = w(0x79622D32); + self.state[3] = w(0x6B206574); + + for i in 0..KEY_WORDS { + self.state[4+i] = w(key[i]); + } + + self.state[12] = w(0); + self.state[13] = w(0); + self.state[14] = w(0); + self.state[15] = w(0); + + self.index = STATE_WORDS; + } + + /// Refill the internal output buffer (`self.buffer`) + fn update(&mut self) { + core(&mut self.buffer, &self.state); + self.index = 0; + // update 128-bit counter + self.state[12] = self.state[12] + w(1); + if self.state[12] != w(0) { return }; + self.state[13] = self.state[13] + w(1); + if self.state[13] != w(0) { return }; + self.state[14] = self.state[14] + w(1); + if self.state[14] != w(0) { return }; + self.state[15] = self.state[15] + w(1); + } +} + +impl Rng for ChaChaRng { + #[inline] + fn next_u32(&mut self) -> u32 { + if self.index == STATE_WORDS { + self.update(); + } + + let value = self.buffer[self.index % STATE_WORDS]; + self.index += 1; + value.0 + } +} + +impl<'a> SeedableRng<&'a [u32]> for ChaChaRng { + + fn reseed(&mut self, seed: &'a [u32]) { + // reset state + self.init(&[0u32; KEY_WORDS]); + // set key in place + let key = &mut self.state[4 .. 4+KEY_WORDS]; + for (k, s) in key.iter_mut().zip(seed.iter()) { + *k = w(*s); + } + } + + /// Create a ChaCha generator from a seed, + /// obtained from a variable-length u32 array. + /// Only up to 8 words are used; if less than 8 + /// words are used, the remaining are set to zero. + fn from_seed(seed: &'a [u32]) -> ChaChaRng { + let mut rng = EMPTY; + rng.reseed(seed); + rng + } +} + +impl Rand for ChaChaRng { + fn rand<R: Rng>(other: &mut R) -> ChaChaRng { + let mut key : [u32; KEY_WORDS] = [0; KEY_WORDS]; + for word in key.iter_mut() { + *word = other.gen(); + } + SeedableRng::from_seed(&key[..]) + } +} + + +#[cfg(test)] +mod test { + use {Rng, SeedableRng}; + use super::ChaChaRng; + + #[test] + fn test_rng_rand_seeded() { + let s = ::test::rng().gen_iter::<u32>().take(8).collect::<Vec<u32>>(); + let mut ra: ChaChaRng = SeedableRng::from_seed(&s[..]); + let mut rb: ChaChaRng = SeedableRng::from_seed(&s[..]); + assert!(::test::iter_eq(ra.gen_ascii_chars().take(100), + rb.gen_ascii_chars().take(100))); + } + + #[test] + fn test_rng_seeded() { + let seed : &[_] = &[0,1,2,3,4,5,6,7]; + let mut ra: ChaChaRng = SeedableRng::from_seed(seed); + let mut rb: ChaChaRng = SeedableRng::from_seed(seed); + assert!(::test::iter_eq(ra.gen_ascii_chars().take(100), + rb.gen_ascii_chars().take(100))); + } + + #[test] + fn test_rng_reseed() { + let s = ::test::rng().gen_iter::<u32>().take(8).collect::<Vec<u32>>(); + let mut r: ChaChaRng = SeedableRng::from_seed(&s[..]); + let string1: String = r.gen_ascii_chars().take(100).collect(); + + r.reseed(&s); + + let string2: String = r.gen_ascii_chars().take(100).collect(); + assert_eq!(string1, string2); + } + + #[test] + fn test_rng_true_values() { + // Test vectors 1 and 2 from + // http://tools.ietf.org/html/draft-nir-cfrg-chacha20-poly1305-04 + let seed : &[_] = &[0u32; 8]; + let mut ra: ChaChaRng = SeedableRng::from_seed(seed); + + let v = (0..16).map(|_| ra.next_u32()).collect::<Vec<_>>(); + assert_eq!(v, + vec!(0xade0b876, 0x903df1a0, 0xe56a5d40, 0x28bd8653, + 0xb819d2bd, 0x1aed8da0, 0xccef36a8, 0xc70d778b, + 0x7c5941da, 0x8d485751, 0x3fe02477, 0x374ad8b8, + 0xf4b8436a, 0x1ca11815, 0x69b687c3, 0x8665eeb2)); + + let v = (0..16).map(|_| ra.next_u32()).collect::<Vec<_>>(); + assert_eq!(v, + vec!(0xbee7079f, 0x7a385155, 0x7c97ba98, 0x0d082d73, + 0xa0290fcb, 0x6965e348, 0x3e53c612, 0xed7aee32, + 0x7621b729, 0x434ee69c, 0xb03371d5, 0xd539d874, + 0x281fed31, 0x45fb0a51, 0x1f0ae1ac, 0x6f4d794b)); + + + let seed : &[_] = &[0,1,2,3,4,5,6,7]; + let mut ra: ChaChaRng = SeedableRng::from_seed(seed); + + // Store the 17*i-th 32-bit word, + // i.e., the i-th word of the i-th 16-word block + let mut v : Vec<u32> = Vec::new(); + for _ in 0..16 { + v.push(ra.next_u32()); + for _ in 0..16 { + ra.next_u32(); + } + } + + assert_eq!(v, + vec!(0xf225c81a, 0x6ab1be57, 0x04d42951, 0x70858036, + 0x49884684, 0x64efec72, 0x4be2d186, 0x3615b384, + 0x11cfa18e, 0xd3c50049, 0x75c775f6, 0x434c6530, + 0x2c5bad8f, 0x898881dc, 0x5f1c86d9, 0xc1f8e7f4)); + } + + #[test] + fn test_rng_clone() { + let seed : &[_] = &[0u32; 8]; + let mut rng: ChaChaRng = SeedableRng::from_seed(seed); + let mut clone = rng.clone(); + for _ in 0..16 { + assert_eq!(rng.next_u64(), clone.next_u64()); + } + } +} diff --git a/rand/src/prng/isaac.rs b/rand/src/prng/isaac.rs new file mode 100644 index 0000000..cf5eb67 --- /dev/null +++ b/rand/src/prng/isaac.rs @@ -0,0 +1,328 @@ +// Copyright 2013 The Rust Project Developers. See the COPYRIGHT +// file at the top-level directory of this distribution and at +// http://rust-lang.org/COPYRIGHT. +// +// Licensed under the Apache License, Version 2.0 <LICENSE-APACHE or +// http://www.apache.org/licenses/LICENSE-2.0> or the MIT license +// <LICENSE-MIT or http://opensource.org/licenses/MIT>, at your +// option. This file may not be copied, modified, or distributed +// except according to those terms. + +//! The ISAAC random number generator. + +#![allow(non_camel_case_types)] + +use core::slice; +use core::iter::repeat; +use core::num::Wrapping as w; +use core::fmt; + +use {Rng, SeedableRng, Rand}; + +#[allow(bad_style)] +type w32 = w<u32>; + +const RAND_SIZE_LEN: usize = 8; +const RAND_SIZE: u32 = 1 << RAND_SIZE_LEN; +const RAND_SIZE_USIZE: usize = 1 << RAND_SIZE_LEN; + +/// A random number generator that uses the ISAAC algorithm[1]. +/// +/// The ISAAC algorithm is generally accepted as suitable for +/// cryptographic purposes, but this implementation has not be +/// verified as such. Prefer a generator like `OsRng` that defers to +/// the operating system for cases that need high security. +/// +/// [1]: Bob Jenkins, [*ISAAC: A fast cryptographic random number +/// generator*](http://www.burtleburtle.net/bob/rand/isaacafa.html) +#[derive(Copy)] +pub struct IsaacRng { + cnt: u32, + rsl: [w32; RAND_SIZE_USIZE], + mem: [w32; RAND_SIZE_USIZE], + a: w32, + b: w32, + c: w32, +} + +static EMPTY: IsaacRng = IsaacRng { + cnt: 0, + rsl: [w(0); RAND_SIZE_USIZE], + mem: [w(0); RAND_SIZE_USIZE], + a: w(0), b: w(0), c: w(0), +}; + +impl IsaacRng { + + /// Create an ISAAC random number generator using the default + /// fixed seed. + pub fn new_unseeded() -> IsaacRng { + let mut rng = EMPTY; + rng.init(false); + rng + } + + /// Initialises `self`. If `use_rsl` is true, then use the current value + /// of `rsl` as a seed, otherwise construct one algorithmically (not + /// randomly). + fn init(&mut self, use_rsl: bool) { + let mut a = w(0x9e3779b9); + let mut b = a; + let mut c = a; + let mut d = a; + let mut e = a; + let mut f = a; + let mut g = a; + let mut h = a; + + macro_rules! mix { + () => {{ + a=a^(b<<11); d=d+a; b=b+c; + b=b^(c>>2); e=e+b; c=c+d; + c=c^(d<<8); f=f+c; d=d+e; + d=d^(e>>16); g=g+d; e=e+f; + e=e^(f<<10); h=h+e; f=f+g; + f=f^(g>>4); a=a+f; g=g+h; + g=g^(h<<8); b=b+g; h=h+a; + h=h^(a>>9); c=c+h; a=a+b; + }} + } + + for _ in 0..4 { + mix!(); + } + + if use_rsl { + macro_rules! memloop { + ($arr:expr) => {{ + for i in (0..RAND_SIZE_USIZE/8).map(|i| i * 8) { + a=a+$arr[i ]; b=b+$arr[i+1]; + c=c+$arr[i+2]; d=d+$arr[i+3]; + e=e+$arr[i+4]; f=f+$arr[i+5]; + g=g+$arr[i+6]; h=h+$arr[i+7]; + mix!(); + self.mem[i ]=a; self.mem[i+1]=b; + self.mem[i+2]=c; self.mem[i+3]=d; + self.mem[i+4]=e; self.mem[i+5]=f; + self.mem[i+6]=g; self.mem[i+7]=h; + } + }} + } + + memloop!(self.rsl); + memloop!(self.mem); + } else { + for i in (0..RAND_SIZE_USIZE/8).map(|i| i * 8) { + mix!(); + self.mem[i ]=a; self.mem[i+1]=b; + self.mem[i+2]=c; self.mem[i+3]=d; + self.mem[i+4]=e; self.mem[i+5]=f; + self.mem[i+6]=g; self.mem[i+7]=h; + } + } + + self.isaac(); + } + + /// Refills the output buffer (`self.rsl`) + #[inline] + fn isaac(&mut self) { + self.c = self.c + w(1); + // abbreviations + let mut a = self.a; + let mut b = self.b + self.c; + + const MIDPOINT: usize = RAND_SIZE_USIZE / 2; + + macro_rules! ind { + ($x:expr) => ( self.mem[($x >> 2usize).0 as usize & (RAND_SIZE_USIZE - 1)] ) + } + + let r = [(0, MIDPOINT), (MIDPOINT, 0)]; + for &(mr_offset, m2_offset) in r.iter() { + + macro_rules! rngstepp { + ($j:expr, $shift:expr) => {{ + let base = $j; + let mix = a << $shift; + + let x = self.mem[base + mr_offset]; + a = (a ^ mix) + self.mem[base + m2_offset]; + let y = ind!(x) + a + b; + self.mem[base + mr_offset] = y; + + b = ind!(y >> RAND_SIZE_LEN) + x; + self.rsl[base + mr_offset] = b; + }} + } + + macro_rules! rngstepn { + ($j:expr, $shift:expr) => {{ + let base = $j; + let mix = a >> $shift; + + let x = self.mem[base + mr_offset]; + a = (a ^ mix) + self.mem[base + m2_offset]; + let y = ind!(x) + a + b; + self.mem[base + mr_offset] = y; + + b = ind!(y >> RAND_SIZE_LEN) + x; + self.rsl[base + mr_offset] = b; + }} + } + + for i in (0..MIDPOINT/4).map(|i| i * 4) { + rngstepp!(i + 0, 13); + rngstepn!(i + 1, 6); + rngstepp!(i + 2, 2); + rngstepn!(i + 3, 16); + } + } + + self.a = a; + self.b = b; + self.cnt = RAND_SIZE; + } +} + +// Cannot be derived because [u32; 256] does not implement Clone +impl Clone for IsaacRng { + fn clone(&self) -> IsaacRng { + *self + } +} + +impl Rng for IsaacRng { + #[inline] + fn next_u32(&mut self) -> u32 { + if self.cnt == 0 { + // make some more numbers + self.isaac(); + } + self.cnt -= 1; + + // self.cnt is at most RAND_SIZE, but that is before the + // subtraction above. We want to index without bounds + // checking, but this could lead to incorrect code if someone + // misrefactors, so we check, sometimes. + // + // (Changes here should be reflected in Isaac64Rng.next_u64.) + debug_assert!(self.cnt < RAND_SIZE); + + // (the % is cheaply telling the optimiser that we're always + // in bounds, without unsafe. NB. this is a power of two, so + // it optimises to a bitwise mask). + self.rsl[(self.cnt % RAND_SIZE) as usize].0 + } +} + +impl<'a> SeedableRng<&'a [u32]> for IsaacRng { + fn reseed(&mut self, seed: &'a [u32]) { + // make the seed into [seed[0], seed[1], ..., seed[seed.len() + // - 1], 0, 0, ...], to fill rng.rsl. + let seed_iter = seed.iter().map(|&x| x).chain(repeat(0u32)); + + for (rsl_elem, seed_elem) in self.rsl.iter_mut().zip(seed_iter) { + *rsl_elem = w(seed_elem); + } + self.cnt = 0; + self.a = w(0); + self.b = w(0); + self.c = w(0); + + self.init(true); + } + + /// Create an ISAAC random number generator with a seed. This can + /// be any length, although the maximum number of elements used is + /// 256 and any more will be silently ignored. A generator + /// constructed with a given seed will generate the same sequence + /// of values as all other generators constructed with that seed. + fn from_seed(seed: &'a [u32]) -> IsaacRng { + let mut rng = EMPTY; + rng.reseed(seed); + rng + } +} + +impl Rand for IsaacRng { + fn rand<R: Rng>(other: &mut R) -> IsaacRng { + let mut ret = EMPTY; + unsafe { + let ptr = ret.rsl.as_mut_ptr() as *mut u8; + + let slice = slice::from_raw_parts_mut(ptr, RAND_SIZE_USIZE * 4); + other.fill_bytes(slice); + } + ret.cnt = 0; + ret.a = w(0); + ret.b = w(0); + ret.c = w(0); + + ret.init(true); + return ret; + } +} + +impl fmt::Debug for IsaacRng { + fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result { + write!(f, "IsaacRng {{}}") + } +} + +#[cfg(test)] +mod test { + use {Rng, SeedableRng}; + use super::IsaacRng; + + #[test] + fn test_rng_32_rand_seeded() { + let s = ::test::rng().gen_iter::<u32>().take(256).collect::<Vec<u32>>(); + let mut ra: IsaacRng = SeedableRng::from_seed(&s[..]); + let mut rb: IsaacRng = SeedableRng::from_seed(&s[..]); + assert!(::test::iter_eq(ra.gen_ascii_chars().take(100), + rb.gen_ascii_chars().take(100))); + } + + #[test] + fn test_rng_32_seeded() { + let seed: &[_] = &[1, 23, 456, 7890, 12345]; + let mut ra: IsaacRng = SeedableRng::from_seed(seed); + let mut rb: IsaacRng = SeedableRng::from_seed(seed); + assert!(::test::iter_eq(ra.gen_ascii_chars().take(100), + rb.gen_ascii_chars().take(100))); + } + + #[test] + fn test_rng_32_reseed() { + let s = ::test::rng().gen_iter::<u32>().take(256).collect::<Vec<u32>>(); + let mut r: IsaacRng = SeedableRng::from_seed(&s[..]); + let string1: String = r.gen_ascii_chars().take(100).collect(); + + r.reseed(&s[..]); + + let string2: String = r.gen_ascii_chars().take(100).collect(); + assert_eq!(string1, string2); + } + + #[test] + fn test_rng_32_true_values() { + let seed: &[_] = &[1, 23, 456, 7890, 12345]; + let mut ra: IsaacRng = SeedableRng::from_seed(seed); + // Regression test that isaac is actually using the above vector + let v = (0..10).map(|_| ra.next_u32()).collect::<Vec<_>>(); + assert_eq!(v, + vec!(2558573138, 873787463, 263499565, 2103644246, 3595684709, + 4203127393, 264982119, 2765226902, 2737944514, 3900253796)); + + let seed: &[_] = &[12345, 67890, 54321, 9876]; + let mut rb: IsaacRng = SeedableRng::from_seed(seed); + // skip forward to the 10000th number + for _ in 0..10000 { rb.next_u32(); } + + let v = (0..10).map(|_| rb.next_u32()).collect::<Vec<_>>(); + assert_eq!(v, + vec!(3676831399, 3183332890, 2834741178, 3854698763, 2717568474, + 1576568959, 3507990155, 179069555, 141456972, 2478885421)); + } +} diff --git a/rand/src/prng/isaac64.rs b/rand/src/prng/isaac64.rs new file mode 100644 index 0000000..b98e3fe --- /dev/null +++ b/rand/src/prng/isaac64.rs @@ -0,0 +1,340 @@ +// Copyright 2013 The Rust Project Developers. See the COPYRIGHT +// file at the top-level directory of this distribution and at +// http://rust-lang.org/COPYRIGHT. +// +// Licensed under the Apache License, Version 2.0 <LICENSE-APACHE or +// http://www.apache.org/licenses/LICENSE-2.0> or the MIT license +// <LICENSE-MIT or http://opensource.org/licenses/MIT>, at your +// option. This file may not be copied, modified, or distributed +// except according to those terms. + +//! The ISAAC-64 random number generator. + +use core::slice; +use core::iter::repeat; +use core::num::Wrapping as w; +use core::fmt; + +use {Rng, SeedableRng, Rand}; + +#[allow(bad_style)] +type w64 = w<u64>; + +const RAND_SIZE_64_LEN: usize = 8; +const RAND_SIZE_64: usize = 1 << RAND_SIZE_64_LEN; + +/// A random number generator that uses ISAAC-64[1], the 64-bit +/// variant of the ISAAC algorithm. +/// +/// The ISAAC algorithm is generally accepted as suitable for +/// cryptographic purposes, but this implementation has not be +/// verified as such. Prefer a generator like `OsRng` that defers to +/// the operating system for cases that need high security. +/// +/// [1]: Bob Jenkins, [*ISAAC: A fast cryptographic random number +/// generator*](http://www.burtleburtle.net/bob/rand/isaacafa.html) +#[derive(Copy)] +pub struct Isaac64Rng { + cnt: usize, + rsl: [w64; RAND_SIZE_64], + mem: [w64; RAND_SIZE_64], + a: w64, + b: w64, + c: w64, +} + +static EMPTY_64: Isaac64Rng = Isaac64Rng { + cnt: 0, + rsl: [w(0); RAND_SIZE_64], + mem: [w(0); RAND_SIZE_64], + a: w(0), b: w(0), c: w(0), +}; + +impl Isaac64Rng { + /// Create a 64-bit ISAAC random number generator using the + /// default fixed seed. + pub fn new_unseeded() -> Isaac64Rng { + let mut rng = EMPTY_64; + rng.init(false); + rng + } + + /// Initialises `self`. If `use_rsl` is true, then use the current value + /// of `rsl` as a seed, otherwise construct one algorithmically (not + /// randomly). + fn init(&mut self, use_rsl: bool) { + macro_rules! init { + ($var:ident) => ( + let mut $var = w(0x9e3779b97f4a7c13); + ) + } + init!(a); init!(b); init!(c); init!(d); + init!(e); init!(f); init!(g); init!(h); + + macro_rules! mix { + () => {{ + a=a-e; f=f^(h>>9); h=h+a; + b=b-f; g=g^(a<<9); a=a+b; + c=c-g; h=h^(b>>23); b=b+c; + d=d-h; a=a^(c<<15); c=c+d; + e=e-a; b=b^(d>>14); d=d+e; + f=f-b; c=c^(e<<20); e=e+f; + g=g-c; d=d^(f>>17); f=f+g; + h=h-d; e=e^(g<<14); g=g+h; + }} + } + + for _ in 0..4 { + mix!(); + } + + if use_rsl { + macro_rules! memloop { + ($arr:expr) => {{ + for i in (0..RAND_SIZE_64 / 8).map(|i| i * 8) { + a=a+$arr[i ]; b=b+$arr[i+1]; + c=c+$arr[i+2]; d=d+$arr[i+3]; + e=e+$arr[i+4]; f=f+$arr[i+5]; + g=g+$arr[i+6]; h=h+$arr[i+7]; + mix!(); + self.mem[i ]=a; self.mem[i+1]=b; + self.mem[i+2]=c; self.mem[i+3]=d; + self.mem[i+4]=e; self.mem[i+5]=f; + self.mem[i+6]=g; self.mem[i+7]=h; + } + }} + } + + memloop!(self.rsl); + memloop!(self.mem); + } else { + for i in (0..RAND_SIZE_64 / 8).map(|i| i * 8) { + mix!(); + self.mem[i ]=a; self.mem[i+1]=b; + self.mem[i+2]=c; self.mem[i+3]=d; + self.mem[i+4]=e; self.mem[i+5]=f; + self.mem[i+6]=g; self.mem[i+7]=h; + } + } + + self.isaac64(); + } + + /// Refills the output buffer (`self.rsl`) + fn isaac64(&mut self) { + self.c = self.c + w(1); + // abbreviations + let mut a = self.a; + let mut b = self.b + self.c; + const MIDPOINT: usize = RAND_SIZE_64 / 2; + const MP_VEC: [(usize, usize); 2] = [(0,MIDPOINT), (MIDPOINT, 0)]; + macro_rules! ind { + ($x:expr) => { + *self.mem.get_unchecked((($x >> 3usize).0 as usize) & (RAND_SIZE_64 - 1)) + } + } + + for &(mr_offset, m2_offset) in MP_VEC.iter() { + for base in (0..MIDPOINT / 4).map(|i| i * 4) { + + macro_rules! rngstepp { + ($j:expr, $shift:expr) => {{ + let base = base + $j; + let mix = a ^ (a << $shift); + let mix = if $j == 0 {!mix} else {mix}; + + unsafe { + let x = *self.mem.get_unchecked(base + mr_offset); + a = mix + *self.mem.get_unchecked(base + m2_offset); + let y = ind!(x) + a + b; + *self.mem.get_unchecked_mut(base + mr_offset) = y; + + b = ind!(y >> RAND_SIZE_64_LEN) + x; + *self.rsl.get_unchecked_mut(base + mr_offset) = b; + } + }} + } + + macro_rules! rngstepn { + ($j:expr, $shift:expr) => {{ + let base = base + $j; + let mix = a ^ (a >> $shift); + let mix = if $j == 0 {!mix} else {mix}; + + unsafe { + let x = *self.mem.get_unchecked(base + mr_offset); + a = mix + *self.mem.get_unchecked(base + m2_offset); + let y = ind!(x) + a + b; + *self.mem.get_unchecked_mut(base + mr_offset) = y; + + b = ind!(y >> RAND_SIZE_64_LEN) + x; + *self.rsl.get_unchecked_mut(base + mr_offset) = b; + } + }} + } + + rngstepp!(0, 21); + rngstepn!(1, 5); + rngstepp!(2, 12); + rngstepn!(3, 33); + } + } + + self.a = a; + self.b = b; + self.cnt = RAND_SIZE_64; + } +} + +// Cannot be derived because [u32; 256] does not implement Clone +impl Clone for Isaac64Rng { + fn clone(&self) -> Isaac64Rng { + *self + } +} + +impl Rng for Isaac64Rng { + #[inline] + fn next_u32(&mut self) -> u32 { + self.next_u64() as u32 + } + + #[inline] + fn next_u64(&mut self) -> u64 { + if self.cnt == 0 { + // make some more numbers + self.isaac64(); + } + self.cnt -= 1; + + // See corresponding location in IsaacRng.next_u32 for + // explanation. + debug_assert!(self.cnt < RAND_SIZE_64); + self.rsl[(self.cnt % RAND_SIZE_64) as usize].0 + } +} + +impl<'a> SeedableRng<&'a [u64]> for Isaac64Rng { + fn reseed(&mut self, seed: &'a [u64]) { + // make the seed into [seed[0], seed[1], ..., seed[seed.len() + // - 1], 0, 0, ...], to fill rng.rsl. + let seed_iter = seed.iter().map(|&x| x).chain(repeat(0u64)); + + for (rsl_elem, seed_elem) in self.rsl.iter_mut().zip(seed_iter) { + *rsl_elem = w(seed_elem); + } + self.cnt = 0; + self.a = w(0); + self.b = w(0); + self.c = w(0); + + self.init(true); + } + + /// Create an ISAAC random number generator with a seed. This can + /// be any length, although the maximum number of elements used is + /// 256 and any more will be silently ignored. A generator + /// constructed with a given seed will generate the same sequence + /// of values as all other generators constructed with that seed. + fn from_seed(seed: &'a [u64]) -> Isaac64Rng { + let mut rng = EMPTY_64; + rng.reseed(seed); + rng + } +} + +impl Rand for Isaac64Rng { + fn rand<R: Rng>(other: &mut R) -> Isaac64Rng { + let mut ret = EMPTY_64; + unsafe { + let ptr = ret.rsl.as_mut_ptr() as *mut u8; + + let slice = slice::from_raw_parts_mut(ptr, RAND_SIZE_64 * 8); + other.fill_bytes(slice); + } + ret.cnt = 0; + ret.a = w(0); + ret.b = w(0); + ret.c = w(0); + + ret.init(true); + return ret; + } +} + +impl fmt::Debug for Isaac64Rng { + fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result { + write!(f, "Isaac64Rng {{}}") + } +} + +#[cfg(test)] +mod test { + use {Rng, SeedableRng}; + use super::Isaac64Rng; + + #[test] + fn test_rng_64_rand_seeded() { + let s = ::test::rng().gen_iter::<u64>().take(256).collect::<Vec<u64>>(); + let mut ra: Isaac64Rng = SeedableRng::from_seed(&s[..]); + let mut rb: Isaac64Rng = SeedableRng::from_seed(&s[..]); + assert!(::test::iter_eq(ra.gen_ascii_chars().take(100), + rb.gen_ascii_chars().take(100))); + } + + #[test] + fn test_rng_64_seeded() { + let seed: &[_] = &[1, 23, 456, 7890, 12345]; + let mut ra: Isaac64Rng = SeedableRng::from_seed(seed); + let mut rb: Isaac64Rng = SeedableRng::from_seed(seed); + assert!(::test::iter_eq(ra.gen_ascii_chars().take(100), + rb.gen_ascii_chars().take(100))); + } + + #[test] + fn test_rng_64_reseed() { + let s = ::test::rng().gen_iter::<u64>().take(256).collect::<Vec<u64>>(); + let mut r: Isaac64Rng = SeedableRng::from_seed(&s[..]); + let string1: String = r.gen_ascii_chars().take(100).collect(); + + r.reseed(&s[..]); + + let string2: String = r.gen_ascii_chars().take(100).collect(); + assert_eq!(string1, string2); + } + + #[test] + fn test_rng_64_true_values() { + let seed: &[_] = &[1, 23, 456, 7890, 12345]; + let mut ra: Isaac64Rng = SeedableRng::from_seed(seed); + // Regression test that isaac is actually using the above vector + let v = (0..10).map(|_| ra.next_u64()).collect::<Vec<_>>(); + assert_eq!(v, + vec!(547121783600835980, 14377643087320773276, 17351601304698403469, + 1238879483818134882, 11952566807690396487, 13970131091560099343, + 4469761996653280935, 15552757044682284409, 6860251611068737823, + 13722198873481261842)); + + let seed: &[_] = &[12345, 67890, 54321, 9876]; + let mut rb: Isaac64Rng = SeedableRng::from_seed(seed); + // skip forward to the 10000th number + for _ in 0..10000 { rb.next_u64(); } + + let v = (0..10).map(|_| rb.next_u64()).collect::<Vec<_>>(); + assert_eq!(v, + vec!(18143823860592706164, 8491801882678285927, 2699425367717515619, + 17196852593171130876, 2606123525235546165, 15790932315217671084, + 596345674630742204, 9947027391921273664, 11788097613744130851, + 10391409374914919106)); + } + + #[test] + fn test_rng_clone() { + let seed: &[_] = &[1, 23, 456, 7890, 12345]; + let mut rng: Isaac64Rng = SeedableRng::from_seed(seed); + let mut clone = rng.clone(); + for _ in 0..16 { + assert_eq!(rng.next_u64(), clone.next_u64()); + } + } +} diff --git a/rand/src/prng/mod.rs b/rand/src/prng/mod.rs new file mode 100644 index 0000000..ed3e018 --- /dev/null +++ b/rand/src/prng/mod.rs @@ -0,0 +1,51 @@ +// Copyright 2017 The Rust Project Developers. See the COPYRIGHT +// file at the top-level directory of this distribution and at +// http://rust-lang.org/COPYRIGHT. +// +// Licensed under the Apache License, Version 2.0 <LICENSE-APACHE or +// http://www.apache.org/licenses/LICENSE-2.0> or the MIT license +// <LICENSE-MIT or http://opensource.org/licenses/MIT>, at your +// option. This file may not be copied, modified, or distributed +// except according to those terms. + +//! Pseudo random number generators are algorithms to produce *apparently +//! random* numbers deterministically, and usually fairly quickly. +//! +//! So long as the algorithm is computationally secure, is initialised with +//! sufficient entropy (i.e. unknown by an attacker), and its internal state is +//! also protected (unknown to an attacker), the output will also be +//! *computationally secure*. Computationally Secure Pseudo Random Number +//! Generators (CSPRNGs) are thus suitable sources of random numbers for +//! cryptography. There are a couple of gotchas here, however. First, the seed +//! used for initialisation must be unknown. Usually this should be provided by +//! the operating system and should usually be secure, however this may not +//! always be the case (especially soon after startup). Second, user-space +//! memory may be vulnerable, for example when written to swap space, and after +//! forking a child process should reinitialise any user-space PRNGs. For this +//! reason it may be preferable to source random numbers directly from the OS +//! for cryptographic applications. +//! +//! PRNGs are also widely used for non-cryptographic uses: randomised +//! algorithms, simulations, games. In these applications it is usually not +//! important for numbers to be cryptographically *unguessable*, but even +//! distribution and independence from other samples (from the point of view +//! of someone unaware of the algorithm used, at least) may still be important. +//! Good PRNGs should satisfy these properties, but do not take them for +//! granted; Wikipedia's article on +//! [Pseudorandom number generators](https://en.wikipedia.org/wiki/Pseudorandom_number_generator) +//! provides some background on this topic. +//! +//! Care should be taken when seeding (initialising) PRNGs. Some PRNGs have +//! short periods for some seeds. If one PRNG is seeded from another using the +//! same algorithm, it is possible that both will yield the same sequence of +//! values (with some lag). + +mod chacha; +mod isaac; +mod isaac64; +mod xorshift; + +pub use self::chacha::ChaChaRng; +pub use self::isaac::IsaacRng; +pub use self::isaac64::Isaac64Rng; +pub use self::xorshift::XorShiftRng; diff --git a/rand/src/prng/xorshift.rs b/rand/src/prng/xorshift.rs new file mode 100644 index 0000000..dd367e9 --- /dev/null +++ b/rand/src/prng/xorshift.rs @@ -0,0 +1,101 @@ +// Copyright 2017 The Rust Project Developers. See the COPYRIGHT +// file at the top-level directory of this distribution and at +// http://rust-lang.org/COPYRIGHT. +// +// Licensed under the Apache License, Version 2.0 <LICENSE-APACHE or +// http://www.apache.org/licenses/LICENSE-2.0> or the MIT license +// <LICENSE-MIT or http://opensource.org/licenses/MIT>, at your +// option. This file may not be copied, modified, or distributed +// except according to those terms. + +//! Xorshift generators + +use core::num::Wrapping as w; +use {Rng, SeedableRng, Rand}; + +/// An Xorshift[1] random number +/// generator. +/// +/// The Xorshift algorithm is not suitable for cryptographic purposes +/// but is very fast. If you do not know for sure that it fits your +/// requirements, use a more secure one such as `IsaacRng` or `OsRng`. +/// +/// [1]: Marsaglia, George (July 2003). ["Xorshift +/// RNGs"](http://www.jstatsoft.org/v08/i14/paper). *Journal of +/// Statistical Software*. Vol. 8 (Issue 14). +#[allow(missing_copy_implementations)] +#[derive(Clone, Debug)] +pub struct XorShiftRng { + x: w<u32>, + y: w<u32>, + z: w<u32>, + w: w<u32>, +} + +impl XorShiftRng { + /// Creates a new XorShiftRng instance which is not seeded. + /// + /// The initial values of this RNG are constants, so all generators created + /// by this function will yield the same stream of random numbers. It is + /// highly recommended that this is created through `SeedableRng` instead of + /// this function + pub fn new_unseeded() -> XorShiftRng { + XorShiftRng { + x: w(0x193a6754), + y: w(0xa8a7d469), + z: w(0x97830e05), + w: w(0x113ba7bb), + } + } +} + +impl Rng for XorShiftRng { + #[inline] + fn next_u32(&mut self) -> u32 { + let x = self.x; + let t = x ^ (x << 11); + self.x = self.y; + self.y = self.z; + self.z = self.w; + let w_ = self.w; + self.w = w_ ^ (w_ >> 19) ^ (t ^ (t >> 8)); + self.w.0 + } +} + +impl SeedableRng<[u32; 4]> for XorShiftRng { + /// Reseed an XorShiftRng. This will panic if `seed` is entirely 0. + fn reseed(&mut self, seed: [u32; 4]) { + assert!(!seed.iter().all(|&x| x == 0), + "XorShiftRng.reseed called with an all zero seed."); + + self.x = w(seed[0]); + self.y = w(seed[1]); + self.z = w(seed[2]); + self.w = w(seed[3]); + } + + /// Create a new XorShiftRng. This will panic if `seed` is entirely 0. + fn from_seed(seed: [u32; 4]) -> XorShiftRng { + assert!(!seed.iter().all(|&x| x == 0), + "XorShiftRng::from_seed called with an all zero seed."); + + XorShiftRng { + x: w(seed[0]), + y: w(seed[1]), + z: w(seed[2]), + w: w(seed[3]), + } + } +} + +impl Rand for XorShiftRng { + fn rand<R: Rng>(rng: &mut R) -> XorShiftRng { + let mut tuple: (u32, u32, u32, u32) = rng.gen(); + while tuple == (0, 0, 0, 0) { + tuple = rng.gen(); + } + let (x, y, z, w_) = tuple; + XorShiftRng { x: w(x), y: w(y), z: w(z), w: w(w_) } + } +} diff --git a/rand/src/rand_impls.rs b/rand/src/rand_impls.rs new file mode 100644 index 0000000..a865bb6 --- /dev/null +++ b/rand/src/rand_impls.rs @@ -0,0 +1,299 @@ +// Copyright 2013-2014 The Rust Project Developers. See the COPYRIGHT +// file at the top-level directory of this distribution and at +// http://rust-lang.org/COPYRIGHT. +// +// Licensed under the Apache License, Version 2.0 <LICENSE-APACHE or +// http://www.apache.org/licenses/LICENSE-2.0> or the MIT license +// <LICENSE-MIT or http://opensource.org/licenses/MIT>, at your +// option. This file may not be copied, modified, or distributed +// except according to those terms. + +//! The implementations of `Rand` for the built-in types. + +use core::{char, mem}; + +use {Rand,Rng}; + +impl Rand for isize { + #[inline] + fn rand<R: Rng>(rng: &mut R) -> isize { + if mem::size_of::<isize>() == 4 { + rng.gen::<i32>() as isize + } else { + rng.gen::<i64>() as isize + } + } +} + +impl Rand for i8 { + #[inline] + fn rand<R: Rng>(rng: &mut R) -> i8 { + rng.next_u32() as i8 + } +} + +impl Rand for i16 { + #[inline] + fn rand<R: Rng>(rng: &mut R) -> i16 { + rng.next_u32() as i16 + } +} + +impl Rand for i32 { + #[inline] + fn rand<R: Rng>(rng: &mut R) -> i32 { + rng.next_u32() as i32 + } +} + +impl Rand for i64 { + #[inline] + fn rand<R: Rng>(rng: &mut R) -> i64 { + rng.next_u64() as i64 + } +} + +#[cfg(feature = "i128_support")] +impl Rand for i128 { + #[inline] + fn rand<R: Rng>(rng: &mut R) -> i128 { + rng.gen::<u128>() as i128 + } +} + +impl Rand for usize { + #[inline] + fn rand<R: Rng>(rng: &mut R) -> usize { + if mem::size_of::<usize>() == 4 { + rng.gen::<u32>() as usize + } else { + rng.gen::<u64>() as usize + } + } +} + +impl Rand for u8 { + #[inline] + fn rand<R: Rng>(rng: &mut R) -> u8 { + rng.next_u32() as u8 + } +} + +impl Rand for u16 { + #[inline] + fn rand<R: Rng>(rng: &mut R) -> u16 { + rng.next_u32() as u16 + } +} + +impl Rand for u32 { + #[inline] + fn rand<R: Rng>(rng: &mut R) -> u32 { + rng.next_u32() + } +} + +impl Rand for u64 { + #[inline] + fn rand<R: Rng>(rng: &mut R) -> u64 { + rng.next_u64() + } +} + +#[cfg(feature = "i128_support")] +impl Rand for u128 { + #[inline] + fn rand<R: Rng>(rng: &mut R) -> u128 { + ((rng.next_u64() as u128) << 64) | (rng.next_u64() as u128) + } +} + + +macro_rules! float_impls { + ($mod_name:ident, $ty:ty, $mantissa_bits:expr, $method_name:ident) => { + mod $mod_name { + use {Rand, Rng, Open01, Closed01}; + + const SCALE: $ty = (1u64 << $mantissa_bits) as $ty; + + impl Rand for $ty { + /// Generate a floating point number in the half-open + /// interval `[0,1)`. + /// + /// See `Closed01` for the closed interval `[0,1]`, + /// and `Open01` for the open interval `(0,1)`. + #[inline] + fn rand<R: Rng>(rng: &mut R) -> $ty { + rng.$method_name() + } + } + impl Rand for Open01<$ty> { + #[inline] + fn rand<R: Rng>(rng: &mut R) -> Open01<$ty> { + // add a small amount (specifically 2 bits below + // the precision of f64/f32 at 1.0), so that small + // numbers are larger than 0, but large numbers + // aren't pushed to/above 1. + Open01(rng.$method_name() + 0.25 / SCALE) + } + } + impl Rand for Closed01<$ty> { + #[inline] + fn rand<R: Rng>(rng: &mut R) -> Closed01<$ty> { + // rescale so that 1.0 - epsilon becomes 1.0 + // precisely. + Closed01(rng.$method_name() * SCALE / (SCALE - 1.0)) + } + } + } + } +} +float_impls! { f64_rand_impls, f64, 53, next_f64 } +float_impls! { f32_rand_impls, f32, 24, next_f32 } + +impl Rand for char { + #[inline] + fn rand<R: Rng>(rng: &mut R) -> char { + // a char is 21 bits + const CHAR_MASK: u32 = 0x001f_ffff; + loop { + // Rejection sampling. About 0.2% of numbers with at most + // 21-bits are invalid codepoints (surrogates), so this + // will succeed first go almost every time. + match char::from_u32(rng.next_u32() & CHAR_MASK) { + Some(c) => return c, + None => {} + } + } + } +} + +impl Rand for bool { + #[inline] + fn rand<R: Rng>(rng: &mut R) -> bool { + rng.gen::<u8>() & 1 == 1 + } +} + +macro_rules! tuple_impl { + // use variables to indicate the arity of the tuple + ($($tyvar:ident),* ) => { + // the trailing commas are for the 1 tuple + impl< + $( $tyvar : Rand ),* + > Rand for ( $( $tyvar ),* , ) { + + #[inline] + fn rand<R: Rng>(_rng: &mut R) -> ( $( $tyvar ),* , ) { + ( + // use the $tyvar's to get the appropriate number of + // repeats (they're not actually needed) + $( + _rng.gen::<$tyvar>() + ),* + , + ) + } + } + } +} + +impl Rand for () { + #[inline] + fn rand<R: Rng>(_: &mut R) -> () { () } +} +tuple_impl!{A} +tuple_impl!{A, B} +tuple_impl!{A, B, C} +tuple_impl!{A, B, C, D} +tuple_impl!{A, B, C, D, E} +tuple_impl!{A, B, C, D, E, F} +tuple_impl!{A, B, C, D, E, F, G} +tuple_impl!{A, B, C, D, E, F, G, H} +tuple_impl!{A, B, C, D, E, F, G, H, I} +tuple_impl!{A, B, C, D, E, F, G, H, I, J} +tuple_impl!{A, B, C, D, E, F, G, H, I, J, K} +tuple_impl!{A, B, C, D, E, F, G, H, I, J, K, L} + +macro_rules! array_impl { + {$n:expr, $t:ident, $($ts:ident,)*} => { + array_impl!{($n - 1), $($ts,)*} + + impl<T> Rand for [T; $n] where T: Rand { + #[inline] + fn rand<R: Rng>(_rng: &mut R) -> [T; $n] { + [_rng.gen::<$t>(), $(_rng.gen::<$ts>()),*] + } + } + }; + {$n:expr,} => { + impl<T> Rand for [T; $n] { + fn rand<R: Rng>(_rng: &mut R) -> [T; $n] { [] } + } + }; +} + +array_impl!{32, T, T, T, T, T, T, T, T, T, T, T, T, T, T, T, T, T, T, T, T, T, T, T, T, T, T, T, T, T, T, T, T,} + +impl<T:Rand> Rand for Option<T> { + #[inline] + fn rand<R: Rng>(rng: &mut R) -> Option<T> { + if rng.gen() { + Some(rng.gen()) + } else { + None + } + } +} + +#[cfg(test)] +mod tests { + use {Rng, thread_rng, Open01, Closed01}; + + struct ConstantRng(u64); + impl Rng for ConstantRng { + fn next_u32(&mut self) -> u32 { + let ConstantRng(v) = *self; + v as u32 + } + fn next_u64(&mut self) -> u64 { + let ConstantRng(v) = *self; + v + } + } + + #[test] + fn floating_point_edge_cases() { + // the test for exact equality is correct here. + assert!(ConstantRng(0xffff_ffff).gen::<f32>() != 1.0); + assert!(ConstantRng(0xffff_ffff_ffff_ffff).gen::<f64>() != 1.0); + } + + #[test] + fn rand_open() { + // this is unlikely to catch an incorrect implementation that + // generates exactly 0 or 1, but it keeps it sane. + let mut rng = thread_rng(); + for _ in 0..1_000 { + // strict inequalities + let Open01(f) = rng.gen::<Open01<f64>>(); + assert!(0.0 < f && f < 1.0); + + let Open01(f) = rng.gen::<Open01<f32>>(); + assert!(0.0 < f && f < 1.0); + } + } + + #[test] + fn rand_closed() { + let mut rng = thread_rng(); + for _ in 0..1_000 { + // strict inequalities + let Closed01(f) = rng.gen::<Closed01<f64>>(); + assert!(0.0 <= f && f <= 1.0); + + let Closed01(f) = rng.gen::<Closed01<f32>>(); + assert!(0.0 <= f && f <= 1.0); + } + } +} diff --git a/rand/src/read.rs b/rand/src/read.rs new file mode 100644 index 0000000..c7351b7 --- /dev/null +++ b/rand/src/read.rs @@ -0,0 +1,123 @@ +// Copyright 2013 The Rust Project Developers. See the COPYRIGHT +// file at the top-level directory of this distribution and at +// http://rust-lang.org/COPYRIGHT. +// +// Licensed under the Apache License, Version 2.0 <LICENSE-APACHE or +// http://www.apache.org/licenses/LICENSE-2.0> or the MIT license +// <LICENSE-MIT or http://opensource.org/licenses/MIT>, at your +// option. This file may not be copied, modified, or distributed +// except according to those terms. + +//! A wrapper around any Read to treat it as an RNG. + +use std::io::{self, Read}; +use std::mem; +use Rng; + +/// An RNG that reads random bytes straight from a `Read`. This will +/// work best with an infinite reader, but this is not required. +/// +/// # Panics +/// +/// It will panic if it there is insufficient data to fulfill a request. +/// +/// # Example +/// +/// ```rust +/// use rand::{read, Rng}; +/// +/// let data = vec![1, 2, 3, 4, 5, 6, 7, 8]; +/// let mut rng = read::ReadRng::new(&data[..]); +/// println!("{:x}", rng.gen::<u32>()); +/// ``` +#[derive(Debug)] +pub struct ReadRng<R> { + reader: R +} + +impl<R: Read> ReadRng<R> { + /// Create a new `ReadRng` from a `Read`. + pub fn new(r: R) -> ReadRng<R> { + ReadRng { + reader: r + } + } +} + +impl<R: Read> Rng for ReadRng<R> { + fn next_u32(&mut self) -> u32 { + // This is designed for speed: reading a LE integer on a LE + // platform just involves blitting the bytes into the memory + // of the u32, similarly for BE on BE; avoiding byteswapping. + let mut buf = [0; 4]; + fill(&mut self.reader, &mut buf).unwrap(); + unsafe { *(buf.as_ptr() as *const u32) } + } + fn next_u64(&mut self) -> u64 { + // see above for explanation. + let mut buf = [0; 8]; + fill(&mut self.reader, &mut buf).unwrap(); + unsafe { *(buf.as_ptr() as *const u64) } + } + fn fill_bytes(&mut self, v: &mut [u8]) { + if v.len() == 0 { return } + fill(&mut self.reader, v).unwrap(); + } +} + +fn fill(r: &mut Read, mut buf: &mut [u8]) -> io::Result<()> { + while buf.len() > 0 { + match try!(r.read(buf)) { + 0 => return Err(io::Error::new(io::ErrorKind::Other, + "end of file reached")), + n => buf = &mut mem::replace(&mut buf, &mut [])[n..], + } + } + Ok(()) +} + +#[cfg(test)] +mod test { + use super::ReadRng; + use Rng; + + #[test] + fn test_reader_rng_u64() { + // transmute from the target to avoid endianness concerns. + let v = vec![0u8, 0, 0, 0, 0, 0, 0, 1, + 0 , 0, 0, 0, 0, 0, 0, 2, + 0, 0, 0, 0, 0, 0, 0, 3]; + let mut rng = ReadRng::new(&v[..]); + + assert_eq!(rng.next_u64(), 1_u64.to_be()); + assert_eq!(rng.next_u64(), 2_u64.to_be()); + assert_eq!(rng.next_u64(), 3_u64.to_be()); + } + #[test] + fn test_reader_rng_u32() { + let v = vec![0u8, 0, 0, 1, 0, 0, 0, 2, 0, 0, 0, 3]; + let mut rng = ReadRng::new(&v[..]); + + assert_eq!(rng.next_u32(), 1_u32.to_be()); + assert_eq!(rng.next_u32(), 2_u32.to_be()); + assert_eq!(rng.next_u32(), 3_u32.to_be()); + } + #[test] + fn test_reader_rng_fill_bytes() { + let v = [1u8, 2, 3, 4, 5, 6, 7, 8]; + let mut w = [0u8; 8]; + + let mut rng = ReadRng::new(&v[..]); + rng.fill_bytes(&mut w); + + assert!(v == w); + } + + #[test] + #[should_panic] + fn test_reader_rng_insufficient_bytes() { + let mut rng = ReadRng::new(&[][..]); + let mut v = [0u8; 3]; + rng.fill_bytes(&mut v); + } +} diff --git a/rand/src/reseeding.rs b/rand/src/reseeding.rs new file mode 100644 index 0000000..1f24e20 --- /dev/null +++ b/rand/src/reseeding.rs @@ -0,0 +1,229 @@ +// Copyright 2013 The Rust Project Developers. See the COPYRIGHT +// file at the top-level directory of this distribution and at +// http://rust-lang.org/COPYRIGHT. +// +// Licensed under the Apache License, Version 2.0 <LICENSE-APACHE or +// http://www.apache.org/licenses/LICENSE-2.0> or the MIT license +// <LICENSE-MIT or http://opensource.org/licenses/MIT>, at your +// option. This file may not be copied, modified, or distributed +// except according to those terms. + +//! A wrapper around another RNG that reseeds it after it +//! generates a certain number of random bytes. + +use core::default::Default; + +use {Rng, SeedableRng}; + +/// How many bytes of entropy the underling RNG is allowed to generate +/// before it is reseeded +const DEFAULT_GENERATION_THRESHOLD: u64 = 32 * 1024; + +/// A wrapper around any RNG which reseeds the underlying RNG after it +/// has generated a certain number of random bytes. +#[derive(Debug)] +pub struct ReseedingRng<R, Rsdr> { + rng: R, + generation_threshold: u64, + bytes_generated: u64, + /// Controls the behaviour when reseeding the RNG. + pub reseeder: Rsdr, +} + +impl<R: Rng, Rsdr: Reseeder<R>> ReseedingRng<R, Rsdr> { + /// Create a new `ReseedingRng` with the given parameters. + /// + /// # Arguments + /// + /// * `rng`: the random number generator to use. + /// * `generation_threshold`: the number of bytes of entropy at which to reseed the RNG. + /// * `reseeder`: the reseeding object to use. + pub fn new(rng: R, generation_threshold: u64, reseeder: Rsdr) -> ReseedingRng<R,Rsdr> { + ReseedingRng { + rng: rng, + generation_threshold: generation_threshold, + bytes_generated: 0, + reseeder: reseeder + } + } + + /// Reseed the internal RNG if the number of bytes that have been + /// generated exceed the threshold. + pub fn reseed_if_necessary(&mut self) { + if self.bytes_generated >= self.generation_threshold { + self.reseeder.reseed(&mut self.rng); + self.bytes_generated = 0; + } + } +} + + +impl<R: Rng, Rsdr: Reseeder<R>> Rng for ReseedingRng<R, Rsdr> { + fn next_u32(&mut self) -> u32 { + self.reseed_if_necessary(); + self.bytes_generated += 4; + self.rng.next_u32() + } + + fn next_u64(&mut self) -> u64 { + self.reseed_if_necessary(); + self.bytes_generated += 8; + self.rng.next_u64() + } + + fn fill_bytes(&mut self, dest: &mut [u8]) { + self.reseed_if_necessary(); + self.bytes_generated += dest.len() as u64; + self.rng.fill_bytes(dest) + } +} + +impl<S, R: SeedableRng<S>, Rsdr: Reseeder<R> + Default> + SeedableRng<(Rsdr, S)> for ReseedingRng<R, Rsdr> { + fn reseed(&mut self, (rsdr, seed): (Rsdr, S)) { + self.rng.reseed(seed); + self.reseeder = rsdr; + self.bytes_generated = 0; + } + + /// Create a new `ReseedingRng` from the given reseeder and + /// seed. This uses a default value for `generation_threshold`. + fn from_seed((rsdr, seed): (Rsdr, S)) -> ReseedingRng<R, Rsdr> { + ReseedingRng { + rng: SeedableRng::from_seed(seed), + generation_threshold: DEFAULT_GENERATION_THRESHOLD, + bytes_generated: 0, + reseeder: rsdr + } + } +} + +/// Something that can be used to reseed an RNG via `ReseedingRng`. +/// +/// # Example +/// +/// ```rust +/// use rand::{Rng, SeedableRng, StdRng}; +/// use rand::reseeding::{Reseeder, ReseedingRng}; +/// +/// struct TickTockReseeder { tick: bool } +/// impl Reseeder<StdRng> for TickTockReseeder { +/// fn reseed(&mut self, rng: &mut StdRng) { +/// let val = if self.tick {0} else {1}; +/// rng.reseed(&[val]); +/// self.tick = !self.tick; +/// } +/// } +/// fn main() { +/// let rsdr = TickTockReseeder { tick: true }; +/// +/// let inner = StdRng::new().unwrap(); +/// let mut rng = ReseedingRng::new(inner, 10, rsdr); +/// +/// // this will repeat, because it gets reseeded very regularly. +/// let s: String = rng.gen_ascii_chars().take(100).collect(); +/// println!("{}", s); +/// } +/// +/// ``` +pub trait Reseeder<R> { + /// Reseed the given RNG. + fn reseed(&mut self, rng: &mut R); +} + +/// Reseed an RNG using a `Default` instance. This reseeds by +/// replacing the RNG with the result of a `Default::default` call. +#[derive(Clone, Copy, Debug)] +pub struct ReseedWithDefault; + +impl<R: Rng + Default> Reseeder<R> for ReseedWithDefault { + fn reseed(&mut self, rng: &mut R) { + *rng = Default::default(); + } +} +impl Default for ReseedWithDefault { + fn default() -> ReseedWithDefault { ReseedWithDefault } +} + +#[cfg(test)] +mod test { + use std::default::Default; + use std::iter::repeat; + use super::{ReseedingRng, ReseedWithDefault}; + use {SeedableRng, Rng}; + + struct Counter { + i: u32 + } + + impl Rng for Counter { + fn next_u32(&mut self) -> u32 { + self.i += 1; + // very random + self.i - 1 + } + } + impl Default for Counter { + fn default() -> Counter { + Counter { i: 0 } + } + } + impl SeedableRng<u32> for Counter { + fn reseed(&mut self, seed: u32) { + self.i = seed; + } + fn from_seed(seed: u32) -> Counter { + Counter { i: seed } + } + } + type MyRng = ReseedingRng<Counter, ReseedWithDefault>; + + #[test] + fn test_reseeding() { + let mut rs = ReseedingRng::new(Counter {i:0}, 400, ReseedWithDefault); + + let mut i = 0; + for _ in 0..1000 { + assert_eq!(rs.next_u32(), i % 100); + i += 1; + } + } + + #[test] + fn test_rng_seeded() { + let mut ra: MyRng = SeedableRng::from_seed((ReseedWithDefault, 2)); + let mut rb: MyRng = SeedableRng::from_seed((ReseedWithDefault, 2)); + assert!(::test::iter_eq(ra.gen_ascii_chars().take(100), + rb.gen_ascii_chars().take(100))); + } + + #[test] + fn test_rng_reseed() { + let mut r: MyRng = SeedableRng::from_seed((ReseedWithDefault, 3)); + let string1: String = r.gen_ascii_chars().take(100).collect(); + + r.reseed((ReseedWithDefault, 3)); + + let string2: String = r.gen_ascii_chars().take(100).collect(); + assert_eq!(string1, string2); + } + + const FILL_BYTES_V_LEN: usize = 13579; + #[test] + fn test_rng_fill_bytes() { + let mut v = repeat(0u8).take(FILL_BYTES_V_LEN).collect::<Vec<_>>(); + ::test::rng().fill_bytes(&mut v); + + // Sanity test: if we've gotten here, `fill_bytes` has not infinitely + // recursed. + assert_eq!(v.len(), FILL_BYTES_V_LEN); + + // To test that `fill_bytes` actually did something, check that the + // average of `v` is not 0. + let mut sum = 0.0; + for &x in v.iter() { + sum += x as f64; + } + assert!(sum / v.len() as f64 != 0.0); + } +} diff --git a/rand/src/seq.rs b/rand/src/seq.rs new file mode 100644 index 0000000..a7889fe --- /dev/null +++ b/rand/src/seq.rs @@ -0,0 +1,337 @@ +// Copyright 2017 The Rust Project Developers. See the COPYRIGHT +// file at the top-level directory of this distribution and at +// http://rust-lang.org/COPYRIGHT. +// +// Licensed under the Apache License, Version 2.0 <LICENSE-APACHE or +// http://www.apache.org/licenses/LICENSE-2.0> or the MIT license +// <LICENSE-MIT or http://opensource.org/licenses/MIT>, at your +// option. This file may not be copied, modified, or distributed +// except according to those terms. + +//! Functions for randomly accessing and sampling sequences. + +use super::Rng; + +// This crate is only enabled when either std or alloc is available. +// BTreeMap is not as fast in tests, but better than nothing. +#[cfg(feature="std")] use std::collections::HashMap; +#[cfg(not(feature="std"))] use alloc::btree_map::BTreeMap; + +#[cfg(not(feature="std"))] use alloc::Vec; + +/// Randomly sample `amount` elements from a finite iterator. +/// +/// The following can be returned: +/// - `Ok`: `Vec` of `amount` non-repeating randomly sampled elements. The order is not random. +/// - `Err`: `Vec` of all the elements from `iterable` in sequential order. This happens when the +/// length of `iterable` was less than `amount`. This is considered an error since exactly +/// `amount` elements is typically expected. +/// +/// This implementation uses `O(len(iterable))` time and `O(amount)` memory. +/// +/// # Example +/// +/// ```rust +/// use rand::{thread_rng, seq}; +/// +/// let mut rng = thread_rng(); +/// let sample = seq::sample_iter(&mut rng, 1..100, 5).unwrap(); +/// println!("{:?}", sample); +/// ``` +pub fn sample_iter<T, I, R>(rng: &mut R, iterable: I, amount: usize) -> Result<Vec<T>, Vec<T>> + where I: IntoIterator<Item=T>, + R: Rng, +{ + let mut iter = iterable.into_iter(); + let mut reservoir = Vec::with_capacity(amount); + reservoir.extend(iter.by_ref().take(amount)); + + // Continue unless the iterator was exhausted + // + // note: this prevents iterators that "restart" from causing problems. + // If the iterator stops once, then so do we. + if reservoir.len() == amount { + for (i, elem) in iter.enumerate() { + let k = rng.gen_range(0, i + 1 + amount); + if let Some(spot) = reservoir.get_mut(k) { + *spot = elem; + } + } + Ok(reservoir) + } else { + // Don't hang onto extra memory. There is a corner case where + // `amount` was much less than `len(iterable)`. + reservoir.shrink_to_fit(); + Err(reservoir) + } +} + +/// Randomly sample exactly `amount` values from `slice`. +/// +/// The values are non-repeating and in random order. +/// +/// This implementation uses `O(amount)` time and memory. +/// +/// Panics if `amount > slice.len()` +/// +/// # Example +/// +/// ```rust +/// use rand::{thread_rng, seq}; +/// +/// let mut rng = thread_rng(); +/// let values = vec![5, 6, 1, 3, 4, 6, 7]; +/// println!("{:?}", seq::sample_slice(&mut rng, &values, 3)); +/// ``` +pub fn sample_slice<R, T>(rng: &mut R, slice: &[T], amount: usize) -> Vec<T> + where R: Rng, + T: Clone +{ + let indices = sample_indices(rng, slice.len(), amount); + + let mut out = Vec::with_capacity(amount); + out.extend(indices.iter().map(|i| slice[*i].clone())); + out +} + +/// Randomly sample exactly `amount` references from `slice`. +/// +/// The references are non-repeating and in random order. +/// +/// This implementation uses `O(amount)` time and memory. +/// +/// Panics if `amount > slice.len()` +/// +/// # Example +/// +/// ```rust +/// use rand::{thread_rng, seq}; +/// +/// let mut rng = thread_rng(); +/// let values = vec![5, 6, 1, 3, 4, 6, 7]; +/// println!("{:?}", seq::sample_slice_ref(&mut rng, &values, 3)); +/// ``` +pub fn sample_slice_ref<'a, R, T>(rng: &mut R, slice: &'a [T], amount: usize) -> Vec<&'a T> + where R: Rng +{ + let indices = sample_indices(rng, slice.len(), amount); + + let mut out = Vec::with_capacity(amount); + out.extend(indices.iter().map(|i| &slice[*i])); + out +} + +/// Randomly sample exactly `amount` indices from `0..length`. +/// +/// The values are non-repeating and in random order. +/// +/// This implementation uses `O(amount)` time and memory. +/// +/// This method is used internally by the slice sampling methods, but it can sometimes be useful to +/// have the indices themselves so this is provided as an alternative. +/// +/// Panics if `amount > length` +pub fn sample_indices<R>(rng: &mut R, length: usize, amount: usize) -> Vec<usize> + where R: Rng, +{ + if amount > length { + panic!("`amount` must be less than or equal to `slice.len()`"); + } + + // We are going to have to allocate at least `amount` for the output no matter what. However, + // if we use the `cached` version we will have to allocate `amount` as a HashMap as well since + // it inserts an element for every loop. + // + // Therefore, if `amount >= length / 2` then inplace will be both faster and use less memory. + // In fact, benchmarks show the inplace version is faster for length up to about 20 times + // faster than amount. + // + // TODO: there is probably even more fine-tuning that can be done here since + // `HashMap::with_capacity(amount)` probably allocates more than `amount` in practice, + // and a trade off could probably be made between memory/cpu, since hashmap operations + // are slower than array index swapping. + if amount >= length / 20 { + sample_indices_inplace(rng, length, amount) + } else { + sample_indices_cache(rng, length, amount) + } +} + +/// Sample an amount of indices using an inplace partial fisher yates method. +/// +/// This allocates the entire `length` of indices and randomizes only the first `amount`. +/// It then truncates to `amount` and returns. +/// +/// This is better than using a HashMap "cache" when `amount >= length / 2` since it does not +/// require allocating an extra cache and is much faster. +fn sample_indices_inplace<R>(rng: &mut R, length: usize, amount: usize) -> Vec<usize> + where R: Rng, +{ + debug_assert!(amount <= length); + let mut indices: Vec<usize> = Vec::with_capacity(length); + indices.extend(0..length); + for i in 0..amount { + let j: usize = rng.gen_range(i, length); + let tmp = indices[i]; + indices[i] = indices[j]; + indices[j] = tmp; + } + indices.truncate(amount); + debug_assert_eq!(indices.len(), amount); + indices +} + + +/// This method performs a partial fisher-yates on a range of indices using a HashMap +/// as a cache to record potential collisions. +/// +/// The cache avoids allocating the entire `length` of values. This is especially useful when +/// `amount <<< length`, i.e. select 3 non-repeating from 1_000_000 +fn sample_indices_cache<R>( + rng: &mut R, + length: usize, + amount: usize, +) -> Vec<usize> + where R: Rng, +{ + debug_assert!(amount <= length); + #[cfg(feature="std")] let mut cache = HashMap::with_capacity(amount); + #[cfg(not(feature="std"))] let mut cache = BTreeMap::new(); + let mut out = Vec::with_capacity(amount); + for i in 0..amount { + let j: usize = rng.gen_range(i, length); + + // equiv: let tmp = slice[i]; + let tmp = match cache.get(&i) { + Some(e) => *e, + None => i, + }; + + // equiv: slice[i] = slice[j]; + let x = match cache.get(&j) { + Some(x) => *x, + None => j, + }; + + // equiv: slice[j] = tmp; + cache.insert(j, tmp); + + // note that in the inplace version, slice[i] is automatically "returned" value + out.push(x); + } + debug_assert_eq!(out.len(), amount); + out +} + +#[cfg(test)] +mod test { + use super::*; + use {thread_rng, XorShiftRng, SeedableRng}; + + #[test] + fn test_sample_iter() { + let min_val = 1; + let max_val = 100; + + let mut r = thread_rng(); + let vals = (min_val..max_val).collect::<Vec<i32>>(); + let small_sample = sample_iter(&mut r, vals.iter(), 5).unwrap(); + let large_sample = sample_iter(&mut r, vals.iter(), vals.len() + 5).unwrap_err(); + + assert_eq!(small_sample.len(), 5); + assert_eq!(large_sample.len(), vals.len()); + // no randomization happens when amount >= len + assert_eq!(large_sample, vals.iter().collect::<Vec<_>>()); + + assert!(small_sample.iter().all(|e| { + **e >= min_val && **e <= max_val + })); + } + #[test] + fn test_sample_slice_boundaries() { + let empty: &[u8] = &[]; + + let mut r = thread_rng(); + + // sample 0 items + assert_eq!(sample_slice(&mut r, empty, 0), vec![]); + assert_eq!(sample_slice(&mut r, &[42, 2, 42], 0), vec![]); + + // sample 1 item + assert_eq!(sample_slice(&mut r, &[42], 1), vec![42]); + let v = sample_slice(&mut r, &[1, 42], 1)[0]; + assert!(v == 1 || v == 42); + + // sample "all" the items + let v = sample_slice(&mut r, &[42, 133], 2); + assert!(v == vec![42, 133] || v == vec![133, 42]); + + assert_eq!(sample_indices_inplace(&mut r, 0, 0), vec![]); + assert_eq!(sample_indices_inplace(&mut r, 1, 0), vec![]); + assert_eq!(sample_indices_inplace(&mut r, 1, 1), vec![0]); + + assert_eq!(sample_indices_cache(&mut r, 0, 0), vec![]); + assert_eq!(sample_indices_cache(&mut r, 1, 0), vec![]); + assert_eq!(sample_indices_cache(&mut r, 1, 1), vec![0]); + + // Make sure lucky 777's aren't lucky + let slice = &[42, 777]; + let mut num_42 = 0; + let total = 1000; + for _ in 0..total { + let v = sample_slice(&mut r, slice, 1); + assert_eq!(v.len(), 1); + let v = v[0]; + assert!(v == 42 || v == 777); + if v == 42 { + num_42 += 1; + } + } + let ratio_42 = num_42 as f64 / 1000 as f64; + assert!(0.4 <= ratio_42 || ratio_42 <= 0.6, "{}", ratio_42); + } + + #[test] + fn test_sample_slice() { + let xor_rng = XorShiftRng::from_seed; + + let max_range = 100; + let mut r = thread_rng(); + + for length in 1usize..max_range { + let amount = r.gen_range(0, length); + let seed: [u32; 4] = [ + r.next_u32(), r.next_u32(), r.next_u32(), r.next_u32() + ]; + + println!("Selecting indices: len={}, amount={}, seed={:?}", length, amount, seed); + + // assert that the two index methods give exactly the same result + let inplace = sample_indices_inplace( + &mut xor_rng(seed), length, amount); + let cache = sample_indices_cache( + &mut xor_rng(seed), length, amount); + assert_eq!(inplace, cache); + + // assert the basics work + let regular = sample_indices( + &mut xor_rng(seed), length, amount); + assert_eq!(regular.len(), amount); + assert!(regular.iter().all(|e| *e < length)); + assert_eq!(regular, inplace); + + // also test that sampling the slice works + let vec: Vec<usize> = (0..length).collect(); + { + let result = sample_slice(&mut xor_rng(seed), &vec, amount); + assert_eq!(result, regular); + } + + { + let result = sample_slice_ref(&mut xor_rng(seed), &vec, amount); + let expected = regular.iter().map(|v| v).collect::<Vec<_>>(); + assert_eq!(result, expected); + } + } + } +} |