rand/distributions/
float.rs

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
9//! Basic floating-point number distributions
10
11use crate::distributions::utils::FloatSIMDUtils;
12use crate::distributions::{Distribution, Standard};
13use crate::Rng;
14use core::mem;
15#[cfg(feature = "simd_support")] use packed_simd::*;
16
17#[cfg(feature = "serde1")]
18use serde::{Serialize, Deserialize};
19
20/// A distribution to sample floating point numbers uniformly in the half-open
21/// interval `(0, 1]`, i.e. including 1 but not 0.
22///
23/// All values that can be generated are of the form `n * ε/2`. For `f32`
24/// the 24 most significant random bits of a `u32` are used and for `f64` the
25/// 53 most significant bits of a `u64` are used. The conversion uses the
26/// multiplicative method.
27///
28/// See also: [`Standard`] which samples from `[0, 1)`, [`Open01`]
29/// which samples from `(0, 1)` and [`Uniform`] which samples from arbitrary
30/// ranges.
31///
32/// # Example
33/// ```
34/// use rand::{thread_rng, Rng};
35/// use rand::distributions::OpenClosed01;
36///
37/// let val: f32 = thread_rng().sample(OpenClosed01);
38/// println!("f32 from (0, 1): {}", val);
39/// ```
40///
41/// [`Standard`]: crate::distributions::Standard
42/// [`Open01`]: crate::distributions::Open01
43/// [`Uniform`]: crate::distributions::uniform::Uniform
44#[derive(Clone, Copy, Debug)]
45#[cfg_attr(feature = "serde1", derive(Serialize, Deserialize))]
46pub struct OpenClosed01;
47
48/// A distribution to sample floating point numbers uniformly in the open
49/// interval `(0, 1)`, i.e. not including either endpoint.
50///
51/// All values that can be generated are of the form `n * ε + ε/2`. For `f32`
52/// the 23 most significant random bits of an `u32` are used, for `f64` 52 from
53/// an `u64`. The conversion uses a transmute-based method.
54///
55/// See also: [`Standard`] which samples from `[0, 1)`, [`OpenClosed01`]
56/// which samples from `(0, 1]` and [`Uniform`] which samples from arbitrary
57/// ranges.
58///
59/// # Example
60/// ```
61/// use rand::{thread_rng, Rng};
62/// use rand::distributions::Open01;
63///
64/// let val: f32 = thread_rng().sample(Open01);
65/// println!("f32 from (0, 1): {}", val);
66/// ```
67///
68/// [`Standard`]: crate::distributions::Standard
69/// [`OpenClosed01`]: crate::distributions::OpenClosed01
70/// [`Uniform`]: crate::distributions::uniform::Uniform
71#[derive(Clone, Copy, Debug)]
72#[cfg_attr(feature = "serde1", derive(Serialize, Deserialize))]
73pub struct Open01;
74
75
76// This trait is needed by both this lib and rand_distr hence is a hidden export
77#[doc(hidden)]
78pub trait IntoFloat {
79    type F;
80
81    /// Helper method to combine the fraction and a constant exponent into a
82    /// float.
83    ///
84    /// Only the least significant bits of `self` may be set, 23 for `f32` and
85    /// 52 for `f64`.
86    /// The resulting value will fall in a range that depends on the exponent.
87    /// As an example the range with exponent 0 will be
88    /// [2<sup>0</sup>..2<sup>1</sup>), which is [1..2).
89    fn into_float_with_exponent(self, exponent: i32) -> Self::F;
90}
91
92macro_rules! float_impls {
93    ($ty:ident, $uty:ident, $f_scalar:ident, $u_scalar:ty,
94     $fraction_bits:expr, $exponent_bias:expr) => {
95        impl IntoFloat for $uty {
96            type F = $ty;
97            #[inline(always)]
98            fn into_float_with_exponent(self, exponent: i32) -> $ty {
99                // The exponent is encoded using an offset-binary representation
100                let exponent_bits: $u_scalar =
101                    (($exponent_bias + exponent) as $u_scalar) << $fraction_bits;
102                $ty::from_bits(self | exponent_bits)
103            }
104        }
105
106        impl Distribution<$ty> for Standard {
107            fn sample<R: Rng + ?Sized>(&self, rng: &mut R) -> $ty {
108                // Multiply-based method; 24/53 random bits; [0, 1) interval.
109                // We use the most significant bits because for simple RNGs
110                // those are usually more random.
111                let float_size = mem::size_of::<$f_scalar>() as u32 * 8;
112                let precision = $fraction_bits + 1;
113                let scale = 1.0 / ((1 as $u_scalar << precision) as $f_scalar);
114
115                let value: $uty = rng.gen();
116                let value = value >> (float_size - precision);
117                scale * $ty::cast_from_int(value)
118            }
119        }
120
121        impl Distribution<$ty> for OpenClosed01 {
122            fn sample<R: Rng + ?Sized>(&self, rng: &mut R) -> $ty {
123                // Multiply-based method; 24/53 random bits; (0, 1] interval.
124                // We use the most significant bits because for simple RNGs
125                // those are usually more random.
126                let float_size = mem::size_of::<$f_scalar>() as u32 * 8;
127                let precision = $fraction_bits + 1;
128                let scale = 1.0 / ((1 as $u_scalar << precision) as $f_scalar);
129
130                let value: $uty = rng.gen();
131                let value = value >> (float_size - precision);
132                // Add 1 to shift up; will not overflow because of right-shift:
133                scale * $ty::cast_from_int(value + 1)
134            }
135        }
136
137        impl Distribution<$ty> for Open01 {
138            fn sample<R: Rng + ?Sized>(&self, rng: &mut R) -> $ty {
139                // Transmute-based method; 23/52 random bits; (0, 1) interval.
140                // We use the most significant bits because for simple RNGs
141                // those are usually more random.
142                use core::$f_scalar::EPSILON;
143                let float_size = mem::size_of::<$f_scalar>() as u32 * 8;
144
145                let value: $uty = rng.gen();
146                let fraction = value >> (float_size - $fraction_bits);
147                fraction.into_float_with_exponent(0) - (1.0 - EPSILON / 2.0)
148            }
149        }
150    }
151}
152
153float_impls! { f32, u32, f32, u32, 23, 127 }
154float_impls! { f64, u64, f64, u64, 52, 1023 }
155
156#[cfg(feature = "simd_support")]
157float_impls! { f32x2, u32x2, f32, u32, 23, 127 }
158#[cfg(feature = "simd_support")]
159float_impls! { f32x4, u32x4, f32, u32, 23, 127 }
160#[cfg(feature = "simd_support")]
161float_impls! { f32x8, u32x8, f32, u32, 23, 127 }
162#[cfg(feature = "simd_support")]
163float_impls! { f32x16, u32x16, f32, u32, 23, 127 }
164
165#[cfg(feature = "simd_support")]
166float_impls! { f64x2, u64x2, f64, u64, 52, 1023 }
167#[cfg(feature = "simd_support")]
168float_impls! { f64x4, u64x4, f64, u64, 52, 1023 }
169#[cfg(feature = "simd_support")]
170float_impls! { f64x8, u64x8, f64, u64, 52, 1023 }
171
172
173#[cfg(test)]
174mod tests {
175    use super::*;
176    use crate::rngs::mock::StepRng;
177
178    const EPSILON32: f32 = ::core::f32::EPSILON;
179    const EPSILON64: f64 = ::core::f64::EPSILON;
180
181    macro_rules! test_f32 {
182        ($fnn:ident, $ty:ident, $ZERO:expr, $EPSILON:expr) => {
183            #[test]
184            fn $fnn() {
185                // Standard
186                let mut zeros = StepRng::new(0, 0);
187                assert_eq!(zeros.gen::<$ty>(), $ZERO);
188                let mut one = StepRng::new(1 << 8 | 1 << (8 + 32), 0);
189                assert_eq!(one.gen::<$ty>(), $EPSILON / 2.0);
190                let mut max = StepRng::new(!0, 0);
191                assert_eq!(max.gen::<$ty>(), 1.0 - $EPSILON / 2.0);
192
193                // OpenClosed01
194                let mut zeros = StepRng::new(0, 0);
195                assert_eq!(zeros.sample::<$ty, _>(OpenClosed01), 0.0 + $EPSILON / 2.0);
196                let mut one = StepRng::new(1 << 8 | 1 << (8 + 32), 0);
197                assert_eq!(one.sample::<$ty, _>(OpenClosed01), $EPSILON);
198                let mut max = StepRng::new(!0, 0);
199                assert_eq!(max.sample::<$ty, _>(OpenClosed01), $ZERO + 1.0);
200
201                // Open01
202                let mut zeros = StepRng::new(0, 0);
203                assert_eq!(zeros.sample::<$ty, _>(Open01), 0.0 + $EPSILON / 2.0);
204                let mut one = StepRng::new(1 << 9 | 1 << (9 + 32), 0);
205                assert_eq!(one.sample::<$ty, _>(Open01), $EPSILON / 2.0 * 3.0);
206                let mut max = StepRng::new(!0, 0);
207                assert_eq!(max.sample::<$ty, _>(Open01), 1.0 - $EPSILON / 2.0);
208            }
209        };
210    }
211    test_f32! { f32_edge_cases, f32, 0.0, EPSILON32 }
212    #[cfg(feature = "simd_support")]
213    test_f32! { f32x2_edge_cases, f32x2, f32x2::splat(0.0), f32x2::splat(EPSILON32) }
214    #[cfg(feature = "simd_support")]
215    test_f32! { f32x4_edge_cases, f32x4, f32x4::splat(0.0), f32x4::splat(EPSILON32) }
216    #[cfg(feature = "simd_support")]
217    test_f32! { f32x8_edge_cases, f32x8, f32x8::splat(0.0), f32x8::splat(EPSILON32) }
218    #[cfg(feature = "simd_support")]
219    test_f32! { f32x16_edge_cases, f32x16, f32x16::splat(0.0), f32x16::splat(EPSILON32) }
220
221    macro_rules! test_f64 {
222        ($fnn:ident, $ty:ident, $ZERO:expr, $EPSILON:expr) => {
223            #[test]
224            fn $fnn() {
225                // Standard
226                let mut zeros = StepRng::new(0, 0);
227                assert_eq!(zeros.gen::<$ty>(), $ZERO);
228                let mut one = StepRng::new(1 << 11, 0);
229                assert_eq!(one.gen::<$ty>(), $EPSILON / 2.0);
230                let mut max = StepRng::new(!0, 0);
231                assert_eq!(max.gen::<$ty>(), 1.0 - $EPSILON / 2.0);
232
233                // OpenClosed01
234                let mut zeros = StepRng::new(0, 0);
235                assert_eq!(zeros.sample::<$ty, _>(OpenClosed01), 0.0 + $EPSILON / 2.0);
236                let mut one = StepRng::new(1 << 11, 0);
237                assert_eq!(one.sample::<$ty, _>(OpenClosed01), $EPSILON);
238                let mut max = StepRng::new(!0, 0);
239                assert_eq!(max.sample::<$ty, _>(OpenClosed01), $ZERO + 1.0);
240
241                // Open01
242                let mut zeros = StepRng::new(0, 0);
243                assert_eq!(zeros.sample::<$ty, _>(Open01), 0.0 + $EPSILON / 2.0);
244                let mut one = StepRng::new(1 << 12, 0);
245                assert_eq!(one.sample::<$ty, _>(Open01), $EPSILON / 2.0 * 3.0);
246                let mut max = StepRng::new(!0, 0);
247                assert_eq!(max.sample::<$ty, _>(Open01), 1.0 - $EPSILON / 2.0);
248            }
249        };
250    }
251    test_f64! { f64_edge_cases, f64, 0.0, EPSILON64 }
252    #[cfg(feature = "simd_support")]
253    test_f64! { f64x2_edge_cases, f64x2, f64x2::splat(0.0), f64x2::splat(EPSILON64) }
254    #[cfg(feature = "simd_support")]
255    test_f64! { f64x4_edge_cases, f64x4, f64x4::splat(0.0), f64x4::splat(EPSILON64) }
256    #[cfg(feature = "simd_support")]
257    test_f64! { f64x8_edge_cases, f64x8, f64x8::splat(0.0), f64x8::splat(EPSILON64) }
258
259    #[test]
260    fn value_stability() {
261        fn test_samples<T: Copy + core::fmt::Debug + PartialEq, D: Distribution<T>>(
262            distr: &D, zero: T, expected: &[T],
263        ) {
264            let mut rng = crate::test::rng(0x6f44f5646c2a7334);
265            let mut buf = [zero; 3];
266            for x in &mut buf {
267                *x = rng.sample(&distr);
268            }
269            assert_eq!(&buf, expected);
270        }
271
272        test_samples(&Standard, 0f32, &[0.0035963655, 0.7346052, 0.09778172]);
273        test_samples(&Standard, 0f64, &[
274            0.7346051961657583,
275            0.20298547462974248,
276            0.8166436635290655,
277        ]);
278
279        test_samples(&OpenClosed01, 0f32, &[0.003596425, 0.73460525, 0.09778178]);
280        test_samples(&OpenClosed01, 0f64, &[
281            0.7346051961657584,
282            0.2029854746297426,
283            0.8166436635290656,
284        ]);
285
286        test_samples(&Open01, 0f32, &[0.0035963655, 0.73460525, 0.09778172]);
287        test_samples(&Open01, 0f64, &[
288            0.7346051961657584,
289            0.20298547462974248,
290            0.8166436635290656,
291        ]);
292
293        #[cfg(feature = "simd_support")]
294        {
295            // We only test a sub-set of types here. Values are identical to
296            // non-SIMD types; we assume this pattern continues across all
297            // SIMD types.
298
299            test_samples(&Standard, f32x2::new(0.0, 0.0), &[
300                f32x2::new(0.0035963655, 0.7346052),
301                f32x2::new(0.09778172, 0.20298547),
302                f32x2::new(0.34296435, 0.81664366),
303            ]);
304
305            test_samples(&Standard, f64x2::new(0.0, 0.0), &[
306                f64x2::new(0.7346051961657583, 0.20298547462974248),
307                f64x2::new(0.8166436635290655, 0.7423708925400552),
308                f64x2::new(0.16387782224016323, 0.9087068770169618),
309            ]);
310        }
311    }
312}
rand/distributions/float.rs

rand/distributions/
float.rs
