1bc8d3f97eb5c958007f2713238472e0c1c8fe02Howard Hinnant//===----------------------------------------------------------------------===// 2bc8d3f97eb5c958007f2713238472e0c1c8fe02Howard Hinnant// 3f5256e16dfc425c1d466f6308d4026d529ce9e0bHoward Hinnant// The LLVM Compiler Infrastructure 4bc8d3f97eb5c958007f2713238472e0c1c8fe02Howard Hinnant// 5b64f8b07c104c6cc986570ac8ee0ed16a9f23976Howard Hinnant// This file is dual licensed under the MIT and the University of Illinois Open 6b64f8b07c104c6cc986570ac8ee0ed16a9f23976Howard Hinnant// Source Licenses. See LICENSE.TXT for details. 7bc8d3f97eb5c958007f2713238472e0c1c8fe02Howard Hinnant// 8bc8d3f97eb5c958007f2713238472e0c1c8fe02Howard Hinnant//===----------------------------------------------------------------------===// 9bc8d3f97eb5c958007f2713238472e0c1c8fe02Howard Hinnant 10bc8d3f97eb5c958007f2713238472e0c1c8fe02Howard Hinnant// <random> 11bc8d3f97eb5c958007f2713238472e0c1c8fe02Howard Hinnant 12bc8d3f97eb5c958007f2713238472e0c1c8fe02Howard Hinnant// class bernoulli_distribution 13bc8d3f97eb5c958007f2713238472e0c1c8fe02Howard Hinnant 14bc8d3f97eb5c958007f2713238472e0c1c8fe02Howard Hinnant// template<class _URNG> result_type operator()(_URNG& g, const param_type& parm); 15bc8d3f97eb5c958007f2713238472e0c1c8fe02Howard Hinnant 16bc8d3f97eb5c958007f2713238472e0c1c8fe02Howard Hinnant#include <random> 176add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant#include <numeric> 186add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant#include <vector> 19bc8d3f97eb5c958007f2713238472e0c1c8fe02Howard Hinnant#include <cassert> 20bc8d3f97eb5c958007f2713238472e0c1c8fe02Howard Hinnant 216add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnanttemplate <class T> 226add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnantinline 236add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard HinnantT 246add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnantsqr(T x) 256add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant{ 266add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant return x * x; 276add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant} 286add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant 29bc8d3f97eb5c958007f2713238472e0c1c8fe02Howard Hinnantint main() 30bc8d3f97eb5c958007f2713238472e0c1c8fe02Howard Hinnant{ 31bc8d3f97eb5c958007f2713238472e0c1c8fe02Howard Hinnant { 32bc8d3f97eb5c958007f2713238472e0c1c8fe02Howard Hinnant typedef std::bernoulli_distribution D; 33bc8d3f97eb5c958007f2713238472e0c1c8fe02Howard Hinnant typedef D::param_type P; 346add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant typedef std::minstd_rand G; 35bc8d3f97eb5c958007f2713238472e0c1c8fe02Howard Hinnant G g; 36bc8d3f97eb5c958007f2713238472e0c1c8fe02Howard Hinnant D d(.75); 37bc8d3f97eb5c958007f2713238472e0c1c8fe02Howard Hinnant P p(.25); 386add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant const int N = 100000; 396add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant std::vector<D::result_type> u; 406add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant for (int i = 0; i < N; ++i) 416add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant u.push_back(d(g, p)); 426add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant double mean = std::accumulate(u.begin(), u.end(), 436add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant double(0)) / u.size(); 446add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant double var = 0; 45df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant double skew = 0; 46df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant double kurtosis = 0; 476add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant for (int i = 0; i < u.size(); ++i) 48df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant { 49df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant double d = (u[i] - mean); 50df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant double d2 = sqr(d); 51df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant var += d2; 52df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant skew += d * d2; 53df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant kurtosis += d2 * d2; 54df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant } 556add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant var /= u.size(); 56df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant double dev = std::sqrt(var); 57df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant skew /= u.size() * dev * var; 58df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant kurtosis /= u.size() * var * var; 59df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant kurtosis -= 3; 606add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant double x_mean = p.p(); 616add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant double x_var = p.p()*(1-p.p()); 62df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant double x_skew = (1 - 2 * p.p())/std::sqrt(x_var); 63df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant double x_kurtosis = (6 * sqr(p.p()) - 6 * p.p() + 1)/x_var; 64d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((mean - x_mean) / x_mean) < 0.01); 65d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((var - x_var) / x_var) < 0.01); 66d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((skew - x_skew) / x_skew) < 0.01); 67d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.02); 686add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant } 696add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant { 706add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant typedef std::bernoulli_distribution D; 716add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant typedef D::param_type P; 726add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant typedef std::minstd_rand G; 736add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant G g; 746add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant D d(.25); 756add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant P p(.75); 766add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant const int N = 100000; 776add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant std::vector<D::result_type> u; 786add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant for (int i = 0; i < N; ++i) 796add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant u.push_back(d(g, p)); 806add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant double mean = std::accumulate(u.begin(), u.end(), 816add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant double(0)) / u.size(); 826add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant double var = 0; 83df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant double skew = 0; 84df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant double kurtosis = 0; 856add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant for (int i = 0; i < u.size(); ++i) 86df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant { 87df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant double d = (u[i] - mean); 88df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant double d2 = sqr(d); 89df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant var += d2; 90df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant skew += d * d2; 91df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant kurtosis += d2 * d2; 92df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant } 936add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant var /= u.size(); 94df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant double dev = std::sqrt(var); 95df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant skew /= u.size() * dev * var; 96df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant kurtosis /= u.size() * var * var; 97df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant kurtosis -= 3; 986add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant double x_mean = p.p(); 996add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant double x_var = p.p()*(1-p.p()); 100df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant double x_skew = (1 - 2 * p.p())/std::sqrt(x_var); 101df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant double x_kurtosis = (6 * sqr(p.p()) - 6 * p.p() + 1)/x_var; 102d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((mean - x_mean) / x_mean) < 0.01); 103d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((var - x_var) / x_var) < 0.01); 104d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((skew - x_skew) / x_skew) < 0.01); 105d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.02); 106bc8d3f97eb5c958007f2713238472e0c1c8fe02Howard Hinnant } 107bc8d3f97eb5c958007f2713238472e0c1c8fe02Howard Hinnant} 108