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, const param_type& parm); 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 { 362bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef std::lognormal_distribution<> D; 372bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef D::param_type P; 382bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef std::mt19937 G; 392bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant G g; 402bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant D d; 412bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant P p(-1./8192, 0.015625); 422bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant const int N = 1000000; 432bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant std::vector<D::result_type> u; 442bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant for (int i = 0; i < N; ++i) 452bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 462bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant D::result_type v = d(g, p); 472bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant assert(v > 0); 482bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant u.push_back(v); 492bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 502bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double mean = std::accumulate(u.begin(), u.end(), 0.0) / u.size(); 512bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double var = 0; 522bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double skew = 0; 532bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double kurtosis = 0; 542bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant for (int i = 0; i < u.size(); ++i) 552bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 562bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double d = (u[i] - mean); 572bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double d2 = sqr(d); 582bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant var += d2; 592bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant skew += d * d2; 602bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis += d2 * d2; 612bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 622bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant var /= u.size(); 632bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double dev = std::sqrt(var); 642bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant skew /= u.size() * dev * var; 652bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis /= u.size() * var * var; 662bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis -= 3; 672bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_mean = std::exp(p.m() + sqr(p.s())/2); 682bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_var = (std::exp(sqr(p.s())) - 1) * std::exp(2*p.m() + sqr(p.s())); 690e20cae1a5be18fba591cd884aa2a389b66a3f49Howard Hinnant double x_skew = (std::exp(sqr(p.s())) + 2) * 702bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant std::sqrt((std::exp(sqr(p.s())) - 1)); 712bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_kurtosis = std::exp(4*sqr(p.s())) + 2*std::exp(3*sqr(p.s())) + 722bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant 3*std::exp(2*sqr(p.s())) - 6; 73d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((mean - x_mean) / x_mean) < 0.01); 74d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((var - x_var) / x_var) < 0.01); 75d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((skew - x_skew) / x_skew) < 0.05); 76d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.25); 772bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 782bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 792bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef std::lognormal_distribution<> D; 802bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef D::param_type P; 812bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef std::mt19937 G; 822bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant G g; 832bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant D d; 842bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant P p(-1./32, 0.25); 852bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant const int N = 1000000; 862bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant std::vector<D::result_type> u; 872bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant for (int i = 0; i < N; ++i) 882bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 892bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant D::result_type v = d(g, p); 902bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant assert(v > 0); 912bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant u.push_back(v); 922bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 932bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double mean = std::accumulate(u.begin(), u.end(), 0.0) / u.size(); 942bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double var = 0; 952bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double skew = 0; 962bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double kurtosis = 0; 972bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant for (int i = 0; i < u.size(); ++i) 982bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 992bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double d = (u[i] - mean); 1002bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double d2 = sqr(d); 1012bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant var += d2; 1022bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant skew += d * d2; 1032bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis += d2 * d2; 1042bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 1052bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant var /= u.size(); 1062bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double dev = std::sqrt(var); 1072bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant skew /= u.size() * dev * var; 1082bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis /= u.size() * var * var; 1092bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis -= 3; 1102bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_mean = std::exp(p.m() + sqr(p.s())/2); 1112bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_var = (std::exp(sqr(p.s())) - 1) * std::exp(2*p.m() + sqr(p.s())); 1120e20cae1a5be18fba591cd884aa2a389b66a3f49Howard Hinnant double x_skew = (std::exp(sqr(p.s())) + 2) * 1132bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant std::sqrt((std::exp(sqr(p.s())) - 1)); 1142bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_kurtosis = std::exp(4*sqr(p.s())) + 2*std::exp(3*sqr(p.s())) + 1152bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant 3*std::exp(2*sqr(p.s())) - 6; 116d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((mean - x_mean) / x_mean) < 0.01); 117d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((var - x_var) / x_var) < 0.01); 118d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((skew - x_skew) / x_skew) < 0.01); 119d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.03); 1202bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 1212bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 1222bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef std::lognormal_distribution<> D; 1232bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef D::param_type P; 1242bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef std::mt19937 G; 1252bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant G g; 1262bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant D d; 1272bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant P p(-1./8, 0.5); 1282bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant const int N = 1000000; 1292bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant std::vector<D::result_type> u; 1302bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant for (int i = 0; i < N; ++i) 1312bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 1322bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant D::result_type v = d(g, p); 1332bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant assert(v > 0); 1342bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant u.push_back(v); 1352bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 1362bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double mean = std::accumulate(u.begin(), u.end(), 0.0) / u.size(); 1372bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double var = 0; 1382bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double skew = 0; 1392bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double kurtosis = 0; 1402bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant for (int i = 0; i < u.size(); ++i) 1412bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 1422bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double d = (u[i] - mean); 1432bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double d2 = sqr(d); 1442bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant var += d2; 1452bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant skew += d * d2; 1462bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis += d2 * d2; 1472bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 1482bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant var /= u.size(); 1492bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double dev = std::sqrt(var); 1502bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant skew /= u.size() * dev * var; 1512bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis /= u.size() * var * var; 1522bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis -= 3; 1532bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_mean = std::exp(p.m() + sqr(p.s())/2); 1542bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_var = (std::exp(sqr(p.s())) - 1) * std::exp(2*p.m() + sqr(p.s())); 1550e20cae1a5be18fba591cd884aa2a389b66a3f49Howard Hinnant double x_skew = (std::exp(sqr(p.s())) + 2) * 1562bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant std::sqrt((std::exp(sqr(p.s())) - 1)); 1572bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_kurtosis = std::exp(4*sqr(p.s())) + 2*std::exp(3*sqr(p.s())) + 1582bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant 3*std::exp(2*sqr(p.s())) - 6; 159d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((mean - x_mean) / x_mean) < 0.01); 160d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((var - x_var) / x_var) < 0.01); 161d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((skew - x_skew) / x_skew) < 0.02); 162d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.05); 1632bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 1642bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 1652bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef std::lognormal_distribution<> D; 1662bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef D::param_type P; 1672bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef std::mt19937 G; 1682bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant G g; 1692bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant D d(3, 4); 1702bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant P p; 1712bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant const int N = 1000000; 1722bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant std::vector<D::result_type> u; 1732bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant for (int i = 0; i < N; ++i) 1742bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 1752bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant D::result_type v = d(g, p); 1762bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant assert(v > 0); 1772bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant u.push_back(v); 1782bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 1792bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double mean = std::accumulate(u.begin(), u.end(), 0.0) / u.size(); 1802bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double var = 0; 1812bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double skew = 0; 1822bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double kurtosis = 0; 1832bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant for (int i = 0; i < u.size(); ++i) 1842bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 1852bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double d = (u[i] - mean); 1862bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double d2 = sqr(d); 1872bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant var += d2; 1882bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant skew += d * d2; 1892bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis += d2 * d2; 1902bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 1912bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant var /= u.size(); 1922bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double dev = std::sqrt(var); 1932bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant skew /= u.size() * dev * var; 1942bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis /= u.size() * var * var; 1952bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis -= 3; 1962bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_mean = std::exp(p.m() + sqr(p.s())/2); 1972bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_var = (std::exp(sqr(p.s())) - 1) * std::exp(2*p.m() + sqr(p.s())); 1980e20cae1a5be18fba591cd884aa2a389b66a3f49Howard Hinnant double x_skew = (std::exp(sqr(p.s())) + 2) * 1992bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant std::sqrt((std::exp(sqr(p.s())) - 1)); 2002bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_kurtosis = std::exp(4*sqr(p.s())) + 2*std::exp(3*sqr(p.s())) + 2012bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant 3*std::exp(2*sqr(p.s())) - 6; 202d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((mean - x_mean) / x_mean) < 0.01); 203d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((var - x_var) / x_var) < 0.02); 204d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((skew - x_skew) / x_skew) < 0.08); 205d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.4); 2062bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 2072bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 2082bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef std::lognormal_distribution<> D; 2092bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef D::param_type P; 2102bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef std::mt19937 G; 2112bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant G g; 2122bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant D d; 2132bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant P p(-0.78125, 1.25); 2142bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant const int N = 1000000; 2152bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant std::vector<D::result_type> u; 2162bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant for (int i = 0; i < N; ++i) 2172bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 2182bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant D::result_type v = d(g, p); 2192bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant assert(v > 0); 2202bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant u.push_back(v); 2212bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 2222bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double mean = std::accumulate(u.begin(), u.end(), 0.0) / u.size(); 2232bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double var = 0; 2242bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double skew = 0; 2252bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double kurtosis = 0; 2262bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant for (int i = 0; i < u.size(); ++i) 2272bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 2282bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double d = (u[i] - mean); 2292bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double d2 = sqr(d); 2302bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant var += d2; 2312bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant skew += d * d2; 2322bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis += d2 * d2; 2332bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 2342bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant var /= u.size(); 2352bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double dev = std::sqrt(var); 2362bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant skew /= u.size() * dev * var; 2372bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis /= u.size() * var * var; 2382bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis -= 3; 2392bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_mean = std::exp(p.m() + sqr(p.s())/2); 2402bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_var = (std::exp(sqr(p.s())) - 1) * std::exp(2*p.m() + sqr(p.s())); 2410e20cae1a5be18fba591cd884aa2a389b66a3f49Howard Hinnant double x_skew = (std::exp(sqr(p.s())) + 2) * 2422bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant std::sqrt((std::exp(sqr(p.s())) - 1)); 2432bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_kurtosis = std::exp(4*sqr(p.s())) + 2*std::exp(3*sqr(p.s())) + 2442bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant 3*std::exp(2*sqr(p.s())) - 6; 245d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((mean - x_mean) / x_mean) < 0.01); 246d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((var - x_var) / x_var) < 0.04); 247d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((skew - x_skew) / x_skew) < 0.2); 248d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.7); 2492bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 2502bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant} 251