Skip to content

Commit 4dfec8a

Browse files
committed
Implement sampling from the unit circle
1 parent 3841dc8 commit 4dfec8a

File tree

3 files changed

+99
-1
lines changed

3 files changed

+99
-1
lines changed

src/distributions/mod.rs

+4
Original file line numberDiff line numberDiff line change
@@ -101,6 +101,7 @@
101101
//! - Multivariate probability distributions
102102
//! - [`Dirichlet`] distribution
103103
//! - [`UnitSphereSurface`] distribution
104+
//! - [`UnitCircle`] distribution
104105
//!
105106
//! # Examples
106107
//!
@@ -171,6 +172,7 @@
171172
//! [`Uniform::new`]: struct.Uniform.html#method.new
172173
//! [`Uniform::new_inclusive`]: struct.Uniform.html#method.new_inclusive
173174
//! [`UnitSphereSurface`]: struct.UnitSphereSurface.html
175+
//! [`UnitCircle`]: struct.UnitCircle.html
174176
//! [`WeightedIndex`]: struct.WeightedIndex.html
175177
176178
use Rng;
@@ -181,6 +183,7 @@ pub use self::float::{OpenClosed01, Open01};
181183
pub use self::bernoulli::Bernoulli;
182184
#[cfg(feature="alloc")] pub use self::weighted::{WeightedIndex, WeightedError};
183185
#[cfg(feature="std")] pub use self::unit_sphere::UnitSphereSurface;
186+
#[cfg(feature="std")] pub use self::unit_circle::UnitCircle;
184187
#[cfg(feature="std")] pub use self::gamma::{Gamma, ChiSquared, FisherF, StudentT};
185188
#[cfg(feature="std")] pub use self::normal::{Normal, LogNormal, StandardNormal};
186189
#[cfg(feature="std")] pub use self::exponential::{Exp, Exp1};
@@ -194,6 +197,7 @@ pub mod uniform;
194197
mod bernoulli;
195198
#[cfg(feature="alloc")] mod weighted;
196199
#[cfg(feature="std")] mod unit_sphere;
200+
#[cfg(feature="std")] mod unit_circle;
197201
#[cfg(feature="std")] mod gamma;
198202
#[cfg(feature="std")] mod normal;
199203
#[cfg(feature="std")] mod exponential;

src/distributions/unit_circle.rs

+94
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,94 @@
1+
use Rng;
2+
use distributions::{Distribution, Uniform};
3+
4+
/// Samples uniformly from the edge of the unit circle in two dimensions.
5+
///
6+
/// Implemented via a method by von Neumann[^1].
7+
///
8+
///
9+
/// # Example
10+
///
11+
/// ```
12+
/// use rand::distributions::{UnitCircle, Distribution};
13+
///
14+
/// let circle = UnitCircle::new();
15+
/// let v = circle.sample(&mut rand::thread_rng());
16+
/// println!("{:?} is from the unit circle.", v)
17+
/// ```
18+
///
19+
/// [^1]: von Neumann, J. (1951) [*Various Techniques Used in Connection with
20+
/// Random Digits.*](https://mcnp.lanl.gov/pdf_files/nbs_vonneumann.pdf)
21+
/// NBS Appl. Math. Ser., No. 12. Washington, DC: U.S. Government Printing
22+
/// Office, pp. 36-38.
23+
#[derive(Clone, Copy, Debug)]
24+
pub struct UnitCircle {
25+
uniform: Uniform<f64>,
26+
}
27+
28+
impl UnitCircle {
29+
/// Construct a new `UnitCircle` distribution.
30+
#[inline]
31+
pub fn new() -> UnitCircle {
32+
UnitCircle { uniform: Uniform::new(-1., 1.) }
33+
}
34+
}
35+
36+
impl Distribution<[f64; 2]> for UnitCircle {
37+
#[inline]
38+
fn sample<R: Rng + ?Sized>(&self, rng: &mut R) -> [f64; 2] {
39+
let mut x1;
40+
let mut x2;
41+
let mut sum;
42+
loop {
43+
x1 = self.uniform.sample(rng);
44+
x2 = self.uniform.sample(rng);
45+
sum = x1*x1 + x2*x2;
46+
if sum < 1. {
47+
break;
48+
}
49+
}
50+
let diff = x1*x1 - x2*x2;
51+
[diff / sum, 2.*x1*x2 / sum]
52+
}
53+
}
54+
55+
#[cfg(test)]
56+
mod tests {
57+
use distributions::Distribution;
58+
use super::UnitCircle;
59+
60+
/// Assert that two numbers are almost equal to each other.
61+
///
62+
/// On panic, this macro will print the values of the expressions with their
63+
/// debug representations.
64+
macro_rules! assert_almost_eq {
65+
($a:expr, $b:expr, $prec:expr) => (
66+
let diff = ($a - $b).abs();
67+
if diff > $prec {
68+
panic!(format!(
69+
"assertion failed: `abs(left - right) = {:.1e} < {:e}`, \
70+
(left: `{}`, right: `{}`)",
71+
diff, $prec, $a, $b));
72+
}
73+
);
74+
}
75+
76+
#[test]
77+
fn norm() {
78+
let mut rng = ::test::rng(1);
79+
let dist = UnitCircle::new();
80+
for _ in 0..1000 {
81+
let x = dist.sample(&mut rng);
82+
assert_almost_eq!(x[0]*x[0] + x[1]*x[1], 1., 1e-15);
83+
}
84+
}
85+
86+
#[test]
87+
fn value_stability() {
88+
let mut rng = ::test::rng(2);
89+
let dist = UnitCircle::new();
90+
assert_eq!(dist.sample(&mut rng), [-0.8150602311723979, 0.5793762331690843]);
91+
assert_eq!(dist.sample(&mut rng), [-0.056204569973983196, 0.998419273809375]);
92+
assert_eq!(dist.sample(&mut rng), [0.7761923749562624, -0.630496151502733]);
93+
}
94+
}

src/distributions/unit_sphere.rs

+1-1
Original file line numberDiff line numberDiff line change
@@ -18,7 +18,7 @@ use distributions::{Distribution, Uniform};
1818
///
1919
/// [^1]: Marsaglia, George (1972). [*Choosing a Point from the Surface of a
2020
/// Sphere.*](https://doi.org/10.1214/aoms/1177692644)
21-
/// Ann. Math. Statist. 43 (1972), no. 2, 645--646.
21+
/// Ann. Math. Statist. 43, no. 2, 645--646.
2222
#[derive(Clone, Copy, Debug)]
2323
pub struct UnitSphereSurface {
2424
uniform: Uniform<f64>,

0 commit comments

Comments
 (0)