]>
Commit | Line | Data |
---|---|---|
0731742a XL |
1 | // Copyright 2018 Developers of the Rand project. |
2 | // | |
3 | // Licensed under the Apache License, Version 2.0 <LICENSE-APACHE or | |
4 | // https://www.apache.org/licenses/LICENSE-2.0> or the MIT license | |
5 | // <LICENSE-MIT or https://opensource.org/licenses/MIT>, at your | |
6 | // option. This file may not be copied, modified, or distributed | |
7 | // except according to those terms. | |
8 | ||
416331ca | 9 | #![allow(deprecated)] |
dfeec247 | 10 | #![allow(clippy::all)] |
416331ca | 11 | |
416331ca | 12 | use crate::distributions::{Distribution, Uniform}; |
dfeec247 | 13 | use crate::Rng; |
0731742a XL |
14 | |
15 | /// Samples uniformly from the surface of the unit sphere in three dimensions. | |
16 | /// | |
17 | /// Implemented via a method by Marsaglia[^1]. | |
18 | /// | |
0731742a XL |
19 | /// [^1]: Marsaglia, George (1972). [*Choosing a Point from the Surface of a |
20 | /// Sphere.*](https://doi.org/10.1214/aoms/1177692644) | |
21 | /// Ann. Math. Statist. 43, no. 2, 645--646. | |
dfeec247 | 22 | #[deprecated(since = "0.7.0", note = "moved to rand_distr crate")] |
0731742a | 23 | #[derive(Clone, Copy, Debug)] |
416331ca | 24 | pub struct UnitSphereSurface; |
0731742a XL |
25 | |
26 | impl UnitSphereSurface { | |
27 | /// Construct a new `UnitSphereSurface` distribution. | |
28 | #[inline] | |
29 | pub fn new() -> UnitSphereSurface { | |
416331ca | 30 | UnitSphereSurface |
0731742a XL |
31 | } |
32 | } | |
33 | ||
34 | impl Distribution<[f64; 3]> for UnitSphereSurface { | |
35 | #[inline] | |
36 | fn sample<R: Rng + ?Sized>(&self, rng: &mut R) -> [f64; 3] { | |
416331ca | 37 | let uniform = Uniform::new(-1., 1.); |
0731742a | 38 | loop { |
416331ca | 39 | let (x1, x2) = (uniform.sample(rng), uniform.sample(rng)); |
dfeec247 | 40 | let sum = x1 * x1 + x2 * x2; |
0731742a XL |
41 | if sum >= 1. { |
42 | continue; | |
43 | } | |
44 | let factor = 2. * (1.0_f64 - sum).sqrt(); | |
dfeec247 | 45 | return [x1 * factor, x2 * factor, 1. - 2. * sum]; |
0731742a XL |
46 | } |
47 | } | |
48 | } | |
49 | ||
50 | #[cfg(test)] | |
51 | mod tests { | |
0731742a | 52 | use super::UnitSphereSurface; |
dfeec247 | 53 | use crate::distributions::Distribution; |
0731742a XL |
54 | |
55 | /// Assert that two numbers are almost equal to each other. | |
56 | /// | |
57 | /// On panic, this macro will print the values of the expressions with their | |
58 | /// debug representations. | |
59 | macro_rules! assert_almost_eq { | |
dfeec247 | 60 | ($a:expr, $b:expr, $prec:expr) => { |
0731742a XL |
61 | let diff = ($a - $b).abs(); |
62 | if diff > $prec { | |
63 | panic!(format!( | |
64 | "assertion failed: `abs(left - right) = {:.1e} < {:e}`, \ | |
65 | (left: `{}`, right: `{}`)", | |
dfeec247 XL |
66 | diff, $prec, $a, $b |
67 | )); | |
0731742a | 68 | } |
dfeec247 | 69 | }; |
0731742a XL |
70 | } |
71 | ||
72 | #[test] | |
73 | fn norm() { | |
416331ca | 74 | let mut rng = crate::test::rng(1); |
0731742a XL |
75 | let dist = UnitSphereSurface::new(); |
76 | for _ in 0..1000 { | |
77 | let x = dist.sample(&mut rng); | |
dfeec247 | 78 | assert_almost_eq!(x[0] * x[0] + x[1] * x[1] + x[2] * x[2], 1., 1e-15); |
0731742a XL |
79 | } |
80 | } | |
81 | ||
82 | #[test] | |
83 | fn value_stability() { | |
416331ca XL |
84 | let mut rng = crate::test::rng(2); |
85 | let expected = [ | |
dfeec247 XL |
86 | [0.03247542860231647, -0.7830477442152738, 0.6211131755296027], |
87 | [-0.09978440840914075, 0.9706650829833128, -0.21875184231323952], | |
88 | [0.2735582468624679, 0.9435374242279655, -0.1868234852870203], | |
89 | ]; | |
416331ca | 90 | let samples = [ |
dfeec247 XL |
91 | UnitSphereSurface.sample(&mut rng), |
92 | UnitSphereSurface.sample(&mut rng), | |
93 | UnitSphereSurface.sample(&mut rng), | |
94 | ]; | |
416331ca | 95 | assert_eq!(samples, expected); |
0731742a XL |
96 | } |
97 | } |