blob: a258bae1613c54d4b600145e3175d3ac963bf13c [file] [log] [blame]
<!DOCTYPE html><html lang="en"><head><meta charset="utf-8"><meta name="viewport" content="width=device-width, initial-scale=1.0"><meta name="generator" content="rustdoc"><meta name="description" content="Source of the Rust file `/root/.cargo/registry/src/github.com-1ecc6299db9ec823/rand_distr-0.4.3/src/geometric.rs`."><meta name="keywords" content="rust, rustlang, rust-lang"><title>geometric.rs - source</title><link rel="preload" as="font" type="font/woff2" crossorigin href="../../SourceSerif4-Regular.ttf.woff2"><link rel="preload" as="font" type="font/woff2" crossorigin href="../../FiraSans-Regular.woff2"><link rel="preload" as="font" type="font/woff2" crossorigin href="../../FiraSans-Medium.woff2"><link rel="preload" as="font" type="font/woff2" crossorigin href="../../SourceCodePro-Regular.ttf.woff2"><link rel="preload" as="font" type="font/woff2" crossorigin href="../../SourceSerif4-Bold.ttf.woff2"><link rel="preload" as="font" type="font/woff2" crossorigin href="../../SourceCodePro-Semibold.ttf.woff2"><link rel="stylesheet" href="../../normalize.css"><link rel="stylesheet" href="../../rustdoc.css" id="mainThemeStyle"><link rel="stylesheet" href="../../ayu.css" disabled><link rel="stylesheet" href="../../dark.css" disabled><link rel="stylesheet" href="../../light.css" id="themeStyle"><script id="default-settings" ></script><script src="../../storage.js"></script><script defer src="../../source-script.js"></script><script defer src="../../source-files.js"></script><script defer src="../../main.js"></script><noscript><link rel="stylesheet" href="../../noscript.css"></noscript><link rel="icon" href="https://www.rust-lang.org/favicon.ico"></head><body class="rustdoc source"><!--[if lte IE 11]><div class="warning">This old browser is unsupported and will most likely display funky things.</div><![endif]--><nav class="sidebar"><a class="sidebar-logo" href="../../rand_distr/index.html"><div class="logo-container">
<img src="https://www.rust-lang.org/logos/rust-logo-128x128-blk.png" alt="logo"></div></a></nav><main><div class="width-limiter"><nav class="sub"><a class="sub-logo-container" href="../../rand_distr/index.html">
<img src="https://www.rust-lang.org/logos/rust-logo-128x128-blk.png" alt="logo"></a><form class="search-form"><div class="search-container"><span></span><input class="search-input" name="search" autocomplete="off" spellcheck="false" placeholder="Click or press ‘S’ to search, ‘?’ for more options…" type="search"><div id="help-button" title="help" tabindex="-1"><a href="../../help.html">?</a></div><div id="settings-menu" tabindex="-1"><a href="../../settings.html" title="settings"><img width="22" height="22" alt="Change settings" src="../../wheel.svg"></a></div></div></form></nav><section id="main-content" class="content"><div class="example-wrap"><pre class="src-line-numbers"><span id="1">1</span>
<span id="2">2</span>
<span id="3">3</span>
<span id="4">4</span>
<span id="5">5</span>
<span id="6">6</span>
<span id="7">7</span>
<span id="8">8</span>
<span id="9">9</span>
<span id="10">10</span>
<span id="11">11</span>
<span id="12">12</span>
<span id="13">13</span>
<span id="14">14</span>
<span id="15">15</span>
<span id="16">16</span>
<span id="17">17</span>
<span id="18">18</span>
<span id="19">19</span>
<span id="20">20</span>
<span id="21">21</span>
<span id="22">22</span>
<span id="23">23</span>
<span id="24">24</span>
<span id="25">25</span>
<span id="26">26</span>
<span id="27">27</span>
<span id="28">28</span>
<span id="29">29</span>
<span id="30">30</span>
<span id="31">31</span>
<span id="32">32</span>
<span id="33">33</span>
<span id="34">34</span>
<span id="35">35</span>
<span id="36">36</span>
<span id="37">37</span>
<span id="38">38</span>
<span id="39">39</span>
<span id="40">40</span>
<span id="41">41</span>
<span id="42">42</span>
<span id="43">43</span>
<span id="44">44</span>
<span id="45">45</span>
<span id="46">46</span>
<span id="47">47</span>
<span id="48">48</span>
<span id="49">49</span>
<span id="50">50</span>
<span id="51">51</span>
<span id="52">52</span>
<span id="53">53</span>
<span id="54">54</span>
<span id="55">55</span>
<span id="56">56</span>
<span id="57">57</span>
<span id="58">58</span>
<span id="59">59</span>
<span id="60">60</span>
<span id="61">61</span>
<span id="62">62</span>
<span id="63">63</span>
<span id="64">64</span>
<span id="65">65</span>
<span id="66">66</span>
<span id="67">67</span>
<span id="68">68</span>
<span id="69">69</span>
<span id="70">70</span>
<span id="71">71</span>
<span id="72">72</span>
<span id="73">73</span>
<span id="74">74</span>
<span id="75">75</span>
<span id="76">76</span>
<span id="77">77</span>
<span id="78">78</span>
<span id="79">79</span>
<span id="80">80</span>
<span id="81">81</span>
<span id="82">82</span>
<span id="83">83</span>
<span id="84">84</span>
<span id="85">85</span>
<span id="86">86</span>
<span id="87">87</span>
<span id="88">88</span>
<span id="89">89</span>
<span id="90">90</span>
<span id="91">91</span>
<span id="92">92</span>
<span id="93">93</span>
<span id="94">94</span>
<span id="95">95</span>
<span id="96">96</span>
<span id="97">97</span>
<span id="98">98</span>
<span id="99">99</span>
<span id="100">100</span>
<span id="101">101</span>
<span id="102">102</span>
<span id="103">103</span>
<span id="104">104</span>
<span id="105">105</span>
<span id="106">106</span>
<span id="107">107</span>
<span id="108">108</span>
<span id="109">109</span>
<span id="110">110</span>
<span id="111">111</span>
<span id="112">112</span>
<span id="113">113</span>
<span id="114">114</span>
<span id="115">115</span>
<span id="116">116</span>
<span id="117">117</span>
<span id="118">118</span>
<span id="119">119</span>
<span id="120">120</span>
<span id="121">121</span>
<span id="122">122</span>
<span id="123">123</span>
<span id="124">124</span>
<span id="125">125</span>
<span id="126">126</span>
<span id="127">127</span>
<span id="128">128</span>
<span id="129">129</span>
<span id="130">130</span>
<span id="131">131</span>
<span id="132">132</span>
<span id="133">133</span>
<span id="134">134</span>
<span id="135">135</span>
<span id="136">136</span>
<span id="137">137</span>
<span id="138">138</span>
<span id="139">139</span>
<span id="140">140</span>
<span id="141">141</span>
<span id="142">142</span>
<span id="143">143</span>
<span id="144">144</span>
<span id="145">145</span>
<span id="146">146</span>
<span id="147">147</span>
<span id="148">148</span>
<span id="149">149</span>
<span id="150">150</span>
<span id="151">151</span>
<span id="152">152</span>
<span id="153">153</span>
<span id="154">154</span>
<span id="155">155</span>
<span id="156">156</span>
<span id="157">157</span>
<span id="158">158</span>
<span id="159">159</span>
<span id="160">160</span>
<span id="161">161</span>
<span id="162">162</span>
<span id="163">163</span>
<span id="164">164</span>
<span id="165">165</span>
<span id="166">166</span>
<span id="167">167</span>
<span id="168">168</span>
<span id="169">169</span>
<span id="170">170</span>
<span id="171">171</span>
<span id="172">172</span>
<span id="173">173</span>
<span id="174">174</span>
<span id="175">175</span>
<span id="176">176</span>
<span id="177">177</span>
<span id="178">178</span>
<span id="179">179</span>
<span id="180">180</span>
<span id="181">181</span>
<span id="182">182</span>
<span id="183">183</span>
<span id="184">184</span>
<span id="185">185</span>
<span id="186">186</span>
<span id="187">187</span>
<span id="188">188</span>
<span id="189">189</span>
<span id="190">190</span>
<span id="191">191</span>
<span id="192">192</span>
<span id="193">193</span>
<span id="194">194</span>
<span id="195">195</span>
<span id="196">196</span>
<span id="197">197</span>
<span id="198">198</span>
<span id="199">199</span>
<span id="200">200</span>
<span id="201">201</span>
<span id="202">202</span>
<span id="203">203</span>
<span id="204">204</span>
<span id="205">205</span>
<span id="206">206</span>
<span id="207">207</span>
<span id="208">208</span>
<span id="209">209</span>
<span id="210">210</span>
<span id="211">211</span>
<span id="212">212</span>
<span id="213">213</span>
<span id="214">214</span>
<span id="215">215</span>
<span id="216">216</span>
<span id="217">217</span>
<span id="218">218</span>
<span id="219">219</span>
<span id="220">220</span>
<span id="221">221</span>
<span id="222">222</span>
<span id="223">223</span>
<span id="224">224</span>
<span id="225">225</span>
<span id="226">226</span>
<span id="227">227</span>
<span id="228">228</span>
<span id="229">229</span>
<span id="230">230</span>
<span id="231">231</span>
<span id="232">232</span>
<span id="233">233</span>
<span id="234">234</span>
<span id="235">235</span>
<span id="236">236</span>
<span id="237">237</span>
<span id="238">238</span>
</pre><pre class="rust"><code><span class="doccomment">//! The geometric distribution.
</span><span class="kw">use </span><span class="kw">crate</span>::Distribution;
<span class="kw">use </span>rand::Rng;
<span class="kw">use </span>core::fmt;
<span class="attribute">#[allow(unused_imports)]
</span><span class="kw">use </span>num_traits::Float;
<span class="doccomment">/// The geometric distribution `Geometric(p)` bounded to `[0, u64::MAX]`.
///
/// This is the probability distribution of the number of failures before the
/// first success in a series of Bernoulli trials. It has the density function
/// `f(k) = (1 - p)^k p` for `k &gt;= 0`, where `p` is the probability of success
/// on each trial.
///
/// This is the discrete analogue of the [exponential distribution](crate::Exp).
///
/// Note that [`StandardGeometric`](crate::StandardGeometric) is an optimised
/// implementation for `p = 0.5`.
///
/// # Example
///
/// ```
/// use rand_distr::{Geometric, Distribution};
///
/// let geo = Geometric::new(0.25).unwrap();
/// let v = geo.sample(&amp;mut rand::thread_rng());
/// println!(&quot;{} is from a Geometric(0.25) distribution&quot;, v);
/// ```
</span><span class="attribute">#[derive(Copy, Clone, Debug)]
#[cfg_attr(feature = <span class="string">&quot;serde1&quot;</span>, derive(serde::Serialize, serde::Deserialize))]
</span><span class="kw">pub struct </span>Geometric
{
p: f64,
pi: f64,
k: u64
}
<span class="doccomment">/// Error type returned from `Geometric::new`.
</span><span class="attribute">#[derive(Clone, Copy, Debug, PartialEq, Eq)]
</span><span class="kw">pub enum </span>Error {
<span class="doccomment">/// `p &lt; 0 || p &gt; 1` or `nan`
</span>InvalidProbability,
}
<span class="kw">impl </span>fmt::Display <span class="kw">for </span>Error {
<span class="kw">fn </span>fmt(<span class="kw-2">&amp;</span><span class="self">self</span>, f: <span class="kw-2">&amp;mut </span>fmt::Formatter&lt;<span class="lifetime">&#39;_</span>&gt;) -&gt; fmt::Result {
f.write_str(<span class="kw">match </span><span class="self">self </span>{
Error::InvalidProbability =&gt; <span class="string">&quot;p is NaN or outside the interval [0, 1] in geometric distribution&quot;</span>,
})
}
}
<span class="attribute">#[cfg(feature = <span class="string">&quot;std&quot;</span>)]
#[cfg_attr(doc_cfg, doc(cfg(feature = <span class="string">&quot;std&quot;</span>)))]
</span><span class="kw">impl </span>std::error::Error <span class="kw">for </span>Error {}
<span class="kw">impl </span>Geometric {
<span class="doccomment">/// Construct a new `Geometric` with the given shape parameter `p`
/// (probability of success on each trial).
</span><span class="kw">pub fn </span>new(p: f64) -&gt; <span class="prelude-ty">Result</span>&lt;<span class="self">Self</span>, Error&gt; {
<span class="kw">if </span>!p.is_finite() || p &lt; <span class="number">0.0 </span>|| p &gt; <span class="number">1.0 </span>{
<span class="prelude-val">Err</span>(Error::InvalidProbability)
} <span class="kw">else if </span>p == <span class="number">0.0 </span>|| p &gt;= <span class="number">2.0 </span>/ <span class="number">3.0 </span>{
<span class="prelude-val">Ok</span>(Geometric { p, pi: p, k: <span class="number">0 </span>})
} <span class="kw">else </span>{
<span class="kw">let </span>(pi, k) = {
<span class="comment">// choose smallest k such that pi = (1 - p)^(2^k) &lt;= 0.5
</span><span class="kw">let </span><span class="kw-2">mut </span>k = <span class="number">1</span>;
<span class="kw">let </span><span class="kw-2">mut </span>pi = (<span class="number">1.0 </span>- p).powi(<span class="number">2</span>);
<span class="kw">while </span>pi &gt; <span class="number">0.5 </span>{
k += <span class="number">1</span>;
pi = pi * pi;
}
(pi, k)
};
<span class="prelude-val">Ok</span>(Geometric { p, pi, k })
}
}
}
<span class="kw">impl </span>Distribution&lt;u64&gt; <span class="kw">for </span>Geometric
{
<span class="kw">fn </span>sample&lt;R: Rng + <span class="question-mark">?</span>Sized&gt;(<span class="kw-2">&amp;</span><span class="self">self</span>, rng: <span class="kw-2">&amp;mut </span>R) -&gt; u64 {
<span class="kw">if </span><span class="self">self</span>.p &gt;= <span class="number">2.0 </span>/ <span class="number">3.0 </span>{
<span class="comment">// use the trivial algorithm:
</span><span class="kw">let </span><span class="kw-2">mut </span>failures = <span class="number">0</span>;
<span class="kw">loop </span>{
<span class="kw">let </span>u = rng.gen::&lt;f64&gt;();
<span class="kw">if </span>u &lt;= <span class="self">self</span>.p { <span class="kw">break</span>; }
failures += <span class="number">1</span>;
}
<span class="kw">return </span>failures;
}
<span class="kw">if </span><span class="self">self</span>.p == <span class="number">0.0 </span>{ <span class="kw">return </span>core::u64::MAX; }
<span class="kw">let </span>Geometric { p, pi, k } = <span class="kw-2">*</span><span class="self">self</span>;
<span class="comment">// Based on the algorithm presented in section 3 of
// Karl Bringmann and Tobias Friedrich (July 2013) - Exact and Efficient
// Generation of Geometric Random Variates and Random Graphs, published
// in International Colloquium on Automata, Languages and Programming
// (pp.267-278)
// https://people.mpi-inf.mpg.de/~kbringma/paper/2013ICALP-1.pdf
// Use the trivial algorithm to sample D from Geo(pi) = Geo(p) / 2^k:
</span><span class="kw">let </span>d = {
<span class="kw">let </span><span class="kw-2">mut </span>failures = <span class="number">0</span>;
<span class="kw">while </span>rng.gen::&lt;f64&gt;() &lt; pi {
failures += <span class="number">1</span>;
}
failures
};
<span class="comment">// Use rejection sampling for the remainder M from Geo(p) % 2^k:
// choose M uniformly from [0, 2^k), but reject with probability (1 - p)^M
// NOTE: The paper suggests using bitwise sampling here, which is
// currently unsupported, but should improve performance by requiring
// fewer iterations on average. ~ October 28, 2020
</span><span class="kw">let </span>m = <span class="kw">loop </span>{
<span class="kw">let </span>m = rng.gen::&lt;u64&gt;() &amp; ((<span class="number">1 </span>&lt;&lt; k) - <span class="number">1</span>);
<span class="kw">let </span>p_reject = <span class="kw">if </span>m &lt;= core::i32::MAX <span class="kw">as </span>u64 {
(<span class="number">1.0 </span>- p).powi(m <span class="kw">as </span>i32)
} <span class="kw">else </span>{
(<span class="number">1.0 </span>- p).powf(m <span class="kw">as </span>f64)
};
<span class="kw">let </span>u = rng.gen::&lt;f64&gt;();
<span class="kw">if </span>u &lt; p_reject {
<span class="kw">break </span>m;
}
};
(d &lt;&lt; k) + m
}
}
<span class="doccomment">/// Samples integers according to the geometric distribution with success
/// probability `p = 0.5`. This is equivalent to `Geometeric::new(0.5)`,
/// but faster.
///
/// See [`Geometric`](crate::Geometric) for the general geometric distribution.
///
/// Implemented via iterated [Rng::gen::&lt;u64&gt;().leading_zeros()].
///
/// # Example
/// ```
/// use rand::prelude::*;
/// use rand_distr::StandardGeometric;
///
/// let v = StandardGeometric.sample(&amp;mut thread_rng());
/// println!(&quot;{} is from a Geometric(0.5) distribution&quot;, v);
/// ```
</span><span class="attribute">#[derive(Copy, Clone, Debug)]
#[cfg_attr(feature = <span class="string">&quot;serde1&quot;</span>, derive(serde::Serialize, serde::Deserialize))]
</span><span class="kw">pub struct </span>StandardGeometric;
<span class="kw">impl </span>Distribution&lt;u64&gt; <span class="kw">for </span>StandardGeometric {
<span class="kw">fn </span>sample&lt;R: Rng + <span class="question-mark">?</span>Sized&gt;(<span class="kw-2">&amp;</span><span class="self">self</span>, rng: <span class="kw-2">&amp;mut </span>R) -&gt; u64 {
<span class="kw">let </span><span class="kw-2">mut </span>result = <span class="number">0</span>;
<span class="kw">loop </span>{
<span class="kw">let </span>x = rng.gen::&lt;u64&gt;().leading_zeros() <span class="kw">as </span>u64;
result += x;
<span class="kw">if </span>x &lt; <span class="number">64 </span>{ <span class="kw">break</span>; }
}
result
}
}
<span class="attribute">#[cfg(test)]
</span><span class="kw">mod </span>test {
<span class="kw">use super</span>::<span class="kw-2">*</span>;
<span class="attribute">#[test]
</span><span class="kw">fn </span>test_geo_invalid_p() {
<span class="macro">assert!</span>(Geometric::new(core::f64::NAN).is_err());
<span class="macro">assert!</span>(Geometric::new(core::f64::INFINITY).is_err());
<span class="macro">assert!</span>(Geometric::new(core::f64::NEG_INFINITY).is_err());
<span class="macro">assert!</span>(Geometric::new(-<span class="number">0.5</span>).is_err());
<span class="macro">assert!</span>(Geometric::new(<span class="number">0.0</span>).is_ok());
<span class="macro">assert!</span>(Geometric::new(<span class="number">1.0</span>).is_ok());
<span class="macro">assert!</span>(Geometric::new(<span class="number">2.0</span>).is_err());
}
<span class="kw">fn </span>test_geo_mean_and_variance&lt;R: Rng&gt;(p: f64, rng: <span class="kw-2">&amp;mut </span>R) {
<span class="kw">let </span>distr = Geometric::new(p).unwrap();
<span class="kw">let </span>expected_mean = (<span class="number">1.0 </span>- p) / p;
<span class="kw">let </span>expected_variance = (<span class="number">1.0 </span>- p) / (p * p);
<span class="kw">let </span><span class="kw-2">mut </span>results = [<span class="number">0.0</span>; <span class="number">10000</span>];
<span class="kw">for </span>i <span class="kw">in </span>results.iter_mut() {
<span class="kw-2">*</span>i = distr.sample(rng) <span class="kw">as </span>f64;
}
<span class="kw">let </span>mean = results.iter().sum::&lt;f64&gt;() / results.len() <span class="kw">as </span>f64;
<span class="macro">assert!</span>((mean <span class="kw">as </span>f64 - expected_mean).abs() &lt; expected_mean / <span class="number">40.0</span>);
<span class="kw">let </span>variance =
results.iter().map(|x| (x - mean) * (x - mean)).sum::&lt;f64&gt;() / results.len() <span class="kw">as </span>f64;
<span class="macro">assert!</span>((variance - expected_variance).abs() &lt; expected_variance / <span class="number">10.0</span>);
}
<span class="attribute">#[test]
</span><span class="kw">fn </span>test_geometric() {
<span class="kw">let </span><span class="kw-2">mut </span>rng = <span class="kw">crate</span>::test::rng(<span class="number">12345</span>);
test_geo_mean_and_variance(<span class="number">0.10</span>, <span class="kw-2">&amp;mut </span>rng);
test_geo_mean_and_variance(<span class="number">0.25</span>, <span class="kw-2">&amp;mut </span>rng);
test_geo_mean_and_variance(<span class="number">0.50</span>, <span class="kw-2">&amp;mut </span>rng);
test_geo_mean_and_variance(<span class="number">0.75</span>, <span class="kw-2">&amp;mut </span>rng);
test_geo_mean_and_variance(<span class="number">0.90</span>, <span class="kw-2">&amp;mut </span>rng);
}
<span class="attribute">#[test]
</span><span class="kw">fn </span>test_standard_geometric() {
<span class="kw">let </span><span class="kw-2">mut </span>rng = <span class="kw">crate</span>::test::rng(<span class="number">654321</span>);
<span class="kw">let </span>distr = StandardGeometric;
<span class="kw">let </span>expected_mean = <span class="number">1.0</span>;
<span class="kw">let </span>expected_variance = <span class="number">2.0</span>;
<span class="kw">let </span><span class="kw-2">mut </span>results = [<span class="number">0.0</span>; <span class="number">1000</span>];
<span class="kw">for </span>i <span class="kw">in </span>results.iter_mut() {
<span class="kw-2">*</span>i = distr.sample(<span class="kw-2">&amp;mut </span>rng) <span class="kw">as </span>f64;
}
<span class="kw">let </span>mean = results.iter().sum::&lt;f64&gt;() / results.len() <span class="kw">as </span>f64;
<span class="macro">assert!</span>((mean <span class="kw">as </span>f64 - expected_mean).abs() &lt; expected_mean / <span class="number">50.0</span>);
<span class="kw">let </span>variance =
results.iter().map(|x| (x - mean) * (x - mean)).sum::&lt;f64&gt;() / results.len() <span class="kw">as </span>f64;
<span class="macro">assert!</span>((variance - expected_variance).abs() &lt; expected_variance / <span class="number">10.0</span>);
}
}
</code></pre></div>
</section></div></main><div id="rustdoc-vars" data-root-path="../../" data-current-crate="rand_distr" data-themes="ayu,dark,light" data-resource-suffix="" data-rustdoc-version="1.66.0-nightly (5c8bff74b 2022-10-21)" ></div></body></html>