blob: fee9d8ce78fbb4c32be751ca2edf088a563b7ed2 [file] [log] [blame]
 // Copyright 2018 Developers of the Rand project. // // Licensed under the Apache License, Version 2.0 or the MIT license // , at your // option. This file may not be copied, modified, or distributed // except according to those terms. use num_traits::Float; use crate::{uniform::SampleUniform, Distribution, Uniform}; use rand::Rng; /// Samples uniformly from the edge of the unit circle in two dimensions. /// /// Implemented via a method by von Neumann[^1]. /// /// /// # Example /// /// ``` /// use rand_distr::{UnitCircle, Distribution}; /// /// let v: [f64; 2] = UnitCircle.sample(&mut rand::thread_rng()); /// println!("{:?} is from the unit circle.", v) /// ``` /// /// [^1]: von Neumann, J. (1951) [*Various Techniques Used in Connection with /// Random Digits.*](https://mcnp.lanl.gov/pdf_files/nbs_vonneumann.pdf) /// NBS Appl. Math. Ser., No. 12. Washington, DC: U.S. Government Printing /// Office, pp. 36-38. #[derive(Clone, Copy, Debug)] pub struct UnitCircle; impl Distribution<[F; 2]> for UnitCircle { #[inline] fn sample(&self, rng: &mut R) -> [F; 2] { let uniform = Uniform::new(F::from(-1.).unwrap(), F::from(1.).unwrap()); let mut x1; let mut x2; let mut sum; loop { x1 = uniform.sample(rng); x2 = uniform.sample(rng); sum = x1 * x1 + x2 * x2; if sum < F::from(1.).unwrap() { break; } } let diff = x1 * x1 - x2 * x2; [diff / sum, F::from(2.).unwrap() * x1 * x2 / sum] } } #[cfg(test)] mod tests { use super::UnitCircle; use crate::Distribution; /// Assert that two numbers are almost equal to each other. /// /// On panic, this macro will print the values of the expressions with their /// debug representations. macro_rules! assert_almost_eq { (\$a:expr, \$b:expr, \$prec:expr) => { let diff = (\$a - \$b).abs(); if diff > \$prec { panic!( "assertion failed: `abs(left - right) = {:.1e} < {:e}`, \ (left: `{}`, right: `{}`)", diff, \$prec, \$a, \$b ); } }; } #[test] fn norm() { let mut rng = crate::test::rng(1); for _ in 0..1000 { let x: [f64; 2] = UnitCircle.sample(&mut rng); assert_almost_eq!(x[0] * x[0] + x[1] * x[1], 1., 1e-15); } } }