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); 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 { 3203aad81afdf77c64e50d0df3f1a8767ea24aa827Howard Hinnant typedef std::bernoulli_distribution D; 336add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant typedef std::minstd_rand G; 34bc8d3f97eb5c958007f2713238472e0c1c8fe02Howard Hinnant G g; 35bc8d3f97eb5c958007f2713238472e0c1c8fe02Howard Hinnant D d(.75); 366add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant const int N = 100000; 376add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant std::vector<D::result_type> u; 386add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant for (int i = 0; i < N; ++i) 396add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant u.push_back(d(g)); 406add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant double mean = std::accumulate(u.begin(), u.end(), 416add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant double(0)) / u.size(); 426add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant double var = 0; 43df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant double skew = 0; 44df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant double kurtosis = 0; 456add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant for (int i = 0; i < u.size(); ++i) 46df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant { 47df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant double d = (u[i] - mean); 48df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant double d2 = sqr(d); 49df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant var += d2; 50df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant skew += d * d2; 51df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant kurtosis += d2 * d2; 52df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant } 536add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant var /= u.size(); 54df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant double dev = std::sqrt(var); 55df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant skew /= u.size() * dev * var; 56df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant kurtosis /= u.size() * var * var; 57df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant kurtosis -= 3; 586add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant double x_mean = d.p(); 596add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant double x_var = d.p()*(1-d.p()); 60df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant double x_skew = (1 - 2 * d.p())/std::sqrt(x_var); 61df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant double x_kurtosis = (6 * sqr(d.p()) - 6 * d.p() + 1)/x_var; 62d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((mean - x_mean) / x_mean) < 0.01); 63d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((var - x_var) / x_var) < 0.01); 64d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((skew - x_skew) / x_skew) < 0.01); 65d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.02); 666add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant } 676add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant { 686add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant typedef std::bernoulli_distribution D; 696add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant typedef std::minstd_rand G; 706add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant G g; 716add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant D d(.25); 726add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant const int N = 100000; 736add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant std::vector<D::result_type> u; 746add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant for (int i = 0; i < N; ++i) 756add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant u.push_back(d(g)); 766add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant double mean = std::accumulate(u.begin(), u.end(), 776add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant double(0)) / u.size(); 786add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant double var = 0; 79df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant double skew = 0; 80df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant double kurtosis = 0; 816add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant for (int i = 0; i < u.size(); ++i) 82df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant { 83df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant double d = (u[i] - mean); 84df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant double d2 = sqr(d); 85df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant var += d2; 86df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant skew += d * d2; 87df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant kurtosis += d2 * d2; 88df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant } 896add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant var /= u.size(); 90df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant double dev = std::sqrt(var); 91df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant skew /= u.size() * dev * var; 92df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant kurtosis /= u.size() * var * var; 93df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant kurtosis -= 3; 946add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant double x_mean = d.p(); 956add8ddfef4081e35ecfb7e67a457bfab97e1fe3Howard Hinnant double x_var = d.p()*(1-d.p()); 96df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant double x_skew = (1 - 2 * d.p())/std::sqrt(x_var); 97df40dc6c1a8ca0bf00fb6aec030f69042f61d974Howard Hinnant double x_kurtosis = (6 * sqr(d.p()) - 6 * d.p() + 1)/x_var; 98d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((mean - x_mean) / x_mean) < 0.01); 99d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((var - x_var) / x_var) < 0.01); 100d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((skew - x_skew) / x_skew) < 0.01); 101d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.02); 102bc8d3f97eb5c958007f2713238472e0c1c8fe02Howard Hinnant } 103bc8d3f97eb5c958007f2713238472e0c1c8fe02Howard Hinnant} 104