diff options
Diffstat (limited to 'rand/src/distributions/exponential.rs')
-rw-r--r-- | rand/src/distributions/exponential.rs | 124 |
1 files changed, 124 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); + } +} |