12bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant//===----------------------------------------------------------------------===// 22bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant// 32bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant// The LLVM Compiler Infrastructure 42bc36fcff3de1ace5c74bb7c1459def41a67e862Howard 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. 72bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant// 82bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant//===----------------------------------------------------------------------===// 906086258d3d8c48a916ec51c33e1ad8f46821b81Dan Albert// 1006086258d3d8c48a916ec51c33e1ad8f46821b81Dan Albert// REQUIRES: long_tests 112bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant 122bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant// <random> 132bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant 142bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant// template<class RealType = double> 152bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant// class lognormal_distribution 162bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant 172bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant// template<class _URNG> result_type operator()(_URNG& g); 182bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant 192bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant#include <random> 202bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant#include <cassert> 212bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant#include <vector> 222bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant#include <numeric> 232bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant 242bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnanttemplate <class T> 252bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnantinline 262bc36fcff3de1ace5c74bb7c1459def41a67e862Howard HinnantT 272bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnantsqr(T x) 282bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant{ 292bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant return x * x; 302bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant} 312bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant 322bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnantint main() 332bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant{ 342bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 352bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef std::lognormal_distribution<> D; 362bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef D::param_type P; 372bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef std::mt19937 G; 382bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant G g; 392bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant D d(-1./8192, 0.015625); 402bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant const int N = 1000000; 412bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant std::vector<D::result_type> u; 422bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant for (int i = 0; i < N; ++i) 432bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 442bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant D::result_type v = d(g); 452bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant assert(v > 0); 462bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant u.push_back(v); 472bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 482bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double mean = std::accumulate(u.begin(), u.end(), 0.0) / u.size(); 492bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double var = 0; 502bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double skew = 0; 512bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double kurtosis = 0; 522bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant for (int i = 0; i < u.size(); ++i) 532bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 542bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double d = (u[i] - mean); 552bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double d2 = sqr(d); 562bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant var += d2; 572bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant skew += d * d2; 582bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis += d2 * d2; 592bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 602bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant var /= u.size(); 612bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double dev = std::sqrt(var); 622bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant skew /= u.size() * dev * var; 632bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis /= u.size() * var * var; 642bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis -= 3; 652bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_mean = std::exp(d.m() + sqr(d.s())/2); 662bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_var = (std::exp(sqr(d.s())) - 1) * std::exp(2*d.m() + sqr(d.s())); 670e20cae1a5be18fba591cd884aa2a389b66a3f49Howard Hinnant double x_skew = (std::exp(sqr(d.s())) + 2) * 682bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant std::sqrt((std::exp(sqr(d.s())) - 1)); 692bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_kurtosis = std::exp(4*sqr(d.s())) + 2*std::exp(3*sqr(d.s())) + 702bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant 3*std::exp(2*sqr(d.s())) - 6; 71d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((mean - x_mean) / x_mean) < 0.01); 72d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((var - x_var) / x_var) < 0.01); 73d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((skew - x_skew) / x_skew) < 0.05); 74d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.25); 752bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 762bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 772bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef std::lognormal_distribution<> D; 782bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef D::param_type P; 792bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef std::mt19937 G; 802bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant G g; 812bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant D d(-1./32, 0.25); 822bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant const int N = 1000000; 832bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant std::vector<D::result_type> u; 842bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant for (int i = 0; i < N; ++i) 852bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 862bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant D::result_type v = d(g); 872bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant assert(v > 0); 882bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant u.push_back(v); 892bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 902bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double mean = std::accumulate(u.begin(), u.end(), 0.0) / u.size(); 912bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double var = 0; 922bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double skew = 0; 932bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double kurtosis = 0; 942bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant for (int i = 0; i < u.size(); ++i) 952bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 962bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double d = (u[i] - mean); 972bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double d2 = sqr(d); 982bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant var += d2; 992bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant skew += d * d2; 1002bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis += d2 * d2; 1012bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 1022bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant var /= u.size(); 1032bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double dev = std::sqrt(var); 1042bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant skew /= u.size() * dev * var; 1052bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis /= u.size() * var * var; 1062bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis -= 3; 1072bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_mean = std::exp(d.m() + sqr(d.s())/2); 1082bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_var = (std::exp(sqr(d.s())) - 1) * std::exp(2*d.m() + sqr(d.s())); 1090e20cae1a5be18fba591cd884aa2a389b66a3f49Howard Hinnant double x_skew = (std::exp(sqr(d.s())) + 2) * 1102bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant std::sqrt((std::exp(sqr(d.s())) - 1)); 1112bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_kurtosis = std::exp(4*sqr(d.s())) + 2*std::exp(3*sqr(d.s())) + 1122bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant 3*std::exp(2*sqr(d.s())) - 6; 113d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((mean - x_mean) / x_mean) < 0.01); 114d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((var - x_var) / x_var) < 0.01); 115d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((skew - x_skew) / x_skew) < 0.01); 116d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.03); 1172bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 1182bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 1192bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef std::lognormal_distribution<> D; 1202bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef D::param_type P; 1212bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef std::mt19937 G; 1222bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant G g; 1232bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant D d(-1./8, 0.5); 1242bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant const int N = 1000000; 1252bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant std::vector<D::result_type> u; 1262bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant for (int i = 0; i < N; ++i) 1272bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 1282bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant D::result_type v = d(g); 1292bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant assert(v > 0); 1302bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant u.push_back(v); 1312bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 1322bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double mean = std::accumulate(u.begin(), u.end(), 0.0) / u.size(); 1332bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double var = 0; 1342bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double skew = 0; 1352bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double kurtosis = 0; 1362bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant for (int i = 0; i < u.size(); ++i) 1372bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 1382bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double d = (u[i] - mean); 1392bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double d2 = sqr(d); 1402bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant var += d2; 1412bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant skew += d * d2; 1422bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis += d2 * d2; 1432bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 1442bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant var /= u.size(); 1452bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double dev = std::sqrt(var); 1462bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant skew /= u.size() * dev * var; 1472bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis /= u.size() * var * var; 1482bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis -= 3; 1492bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_mean = std::exp(d.m() + sqr(d.s())/2); 1502bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_var = (std::exp(sqr(d.s())) - 1) * std::exp(2*d.m() + sqr(d.s())); 1510e20cae1a5be18fba591cd884aa2a389b66a3f49Howard Hinnant double x_skew = (std::exp(sqr(d.s())) + 2) * 1522bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant std::sqrt((std::exp(sqr(d.s())) - 1)); 1532bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_kurtosis = std::exp(4*sqr(d.s())) + 2*std::exp(3*sqr(d.s())) + 1542bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant 3*std::exp(2*sqr(d.s())) - 6; 155d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((mean - x_mean) / x_mean) < 0.01); 156d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((var - x_var) / x_var) < 0.01); 157d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((skew - x_skew) / x_skew) < 0.02); 158d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.05); 1592bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 1602bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 1612bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef std::lognormal_distribution<> D; 1622bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef D::param_type P; 1632bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef std::mt19937 G; 1642bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant G g; 1652bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant D d; 1662bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant const int N = 1000000; 1672bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant std::vector<D::result_type> u; 1682bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant for (int i = 0; i < N; ++i) 1692bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 1702bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant D::result_type v = d(g); 1712bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant assert(v > 0); 1722bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant u.push_back(v); 1732bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 1742bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double mean = std::accumulate(u.begin(), u.end(), 0.0) / u.size(); 1752bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double var = 0; 1762bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double skew = 0; 1772bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double kurtosis = 0; 1782bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant for (int i = 0; i < u.size(); ++i) 1792bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 1802bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double d = (u[i] - mean); 1812bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double d2 = sqr(d); 1822bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant var += d2; 1832bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant skew += d * d2; 1842bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis += d2 * d2; 1852bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 1862bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant var /= u.size(); 1872bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double dev = std::sqrt(var); 1882bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant skew /= u.size() * dev * var; 1892bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis /= u.size() * var * var; 1902bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis -= 3; 1912bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_mean = std::exp(d.m() + sqr(d.s())/2); 1922bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_var = (std::exp(sqr(d.s())) - 1) * std::exp(2*d.m() + sqr(d.s())); 1930e20cae1a5be18fba591cd884aa2a389b66a3f49Howard Hinnant double x_skew = (std::exp(sqr(d.s())) + 2) * 1942bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant std::sqrt((std::exp(sqr(d.s())) - 1)); 1952bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_kurtosis = std::exp(4*sqr(d.s())) + 2*std::exp(3*sqr(d.s())) + 1962bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant 3*std::exp(2*sqr(d.s())) - 6; 197d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((mean - x_mean) / x_mean) < 0.01); 198d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((var - x_var) / x_var) < 0.02); 199d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((skew - x_skew) / x_skew) < 0.08); 200d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.4); 2012bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 2022bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 2032bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef std::lognormal_distribution<> D; 2042bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef D::param_type P; 2052bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef std::mt19937 G; 2062bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant G g; 2072bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant D d(-0.78125, 1.25); 2082bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant const int N = 1000000; 2092bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant std::vector<D::result_type> u; 2102bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant for (int i = 0; i < N; ++i) 2112bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 2122bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant D::result_type v = d(g); 2132bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant assert(v > 0); 2142bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant u.push_back(v); 2152bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 2162bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double mean = std::accumulate(u.begin(), u.end(), 0.0) / u.size(); 2172bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double var = 0; 2182bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double skew = 0; 2192bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double kurtosis = 0; 2202bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant for (int i = 0; i < u.size(); ++i) 2212bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 2222bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double d = (u[i] - mean); 2232bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double d2 = sqr(d); 2242bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant var += d2; 2252bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant skew += d * d2; 2262bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis += d2 * d2; 2272bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 2282bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant var /= u.size(); 2292bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double dev = std::sqrt(var); 2302bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant skew /= u.size() * dev * var; 2312bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis /= u.size() * var * var; 2322bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis -= 3; 2332bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_mean = std::exp(d.m() + sqr(d.s())/2); 2342bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_var = (std::exp(sqr(d.s())) - 1) * std::exp(2*d.m() + sqr(d.s())); 2350e20cae1a5be18fba591cd884aa2a389b66a3f49Howard Hinnant double x_skew = (std::exp(sqr(d.s())) + 2) * 2362bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant std::sqrt((std::exp(sqr(d.s())) - 1)); 2372bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_kurtosis = std::exp(4*sqr(d.s())) + 2*std::exp(3*sqr(d.s())) + 2382bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant 3*std::exp(2*sqr(d.s())) - 6; 239d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((mean - x_mean) / x_mean) < 0.01); 240d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((var - x_var) / x_var) < 0.04); 241d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((skew - x_skew) / x_skew) < 0.2); 242d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.7); 2432bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 2442bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant} 245