eval.pass.cpp revision b64f8b07c104c6cc986570ac8ee0ed16a9f23976
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//===----------------------------------------------------------------------===// 92bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant 102bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant// <random> 112bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant 122bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant// template<class RealType = double> 132bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant// class lognormal_distribution 142bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant 152bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant// template<class _URNG> result_type operator()(_URNG& g); 162bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant 172bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant#include <random> 182bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant#include <cassert> 192bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant#include <vector> 202bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant#include <numeric> 212bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant 222bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnanttemplate <class T> 232bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnantinline 242bc36fcff3de1ace5c74bb7c1459def41a67e862Howard HinnantT 252bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnantsqr(T x) 262bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant{ 272bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant return x * x; 282bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant} 292bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant 302bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnantint main() 312bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant{ 322bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 332bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef std::lognormal_distribution<> D; 342bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef D::param_type P; 352bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef std::mt19937 G; 362bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant G g; 372bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant D d(-1./8192, 0.015625); 382bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant const int N = 1000000; 392bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant std::vector<D::result_type> u; 402bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant for (int i = 0; i < N; ++i) 412bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 422bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant D::result_type v = d(g); 432bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant assert(v > 0); 442bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant u.push_back(v); 452bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 462bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double mean = std::accumulate(u.begin(), u.end(), 0.0) / u.size(); 472bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double var = 0; 482bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double skew = 0; 492bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double kurtosis = 0; 502bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant for (int i = 0; i < u.size(); ++i) 512bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 522bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double d = (u[i] - mean); 532bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double d2 = sqr(d); 542bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant var += d2; 552bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant skew += d * d2; 562bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis += d2 * d2; 572bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 582bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant var /= u.size(); 592bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double dev = std::sqrt(var); 602bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant skew /= u.size() * dev * var; 612bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis /= u.size() * var * var; 622bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis -= 3; 632bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_mean = std::exp(d.m() + sqr(d.s())/2); 642bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_var = (std::exp(sqr(d.s())) - 1) * std::exp(2*d.m() + sqr(d.s())); 650e20cae1a5be18fba591cd884aa2a389b66a3f49Howard Hinnant double x_skew = (std::exp(sqr(d.s())) + 2) * 662bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant std::sqrt((std::exp(sqr(d.s())) - 1)); 672bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_kurtosis = std::exp(4*sqr(d.s())) + 2*std::exp(3*sqr(d.s())) + 682bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant 3*std::exp(2*sqr(d.s())) - 6; 69d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((mean - x_mean) / x_mean) < 0.01); 70d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((var - x_var) / x_var) < 0.01); 71d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((skew - x_skew) / x_skew) < 0.05); 72d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.25); 732bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 742bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 752bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef std::lognormal_distribution<> D; 762bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef D::param_type P; 772bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef std::mt19937 G; 782bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant G g; 792bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant D d(-1./32, 0.25); 802bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant const int N = 1000000; 812bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant std::vector<D::result_type> u; 822bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant for (int i = 0; i < N; ++i) 832bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 842bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant D::result_type v = d(g); 852bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant assert(v > 0); 862bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant u.push_back(v); 872bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 882bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double mean = std::accumulate(u.begin(), u.end(), 0.0) / u.size(); 892bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double var = 0; 902bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double skew = 0; 912bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double kurtosis = 0; 922bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant for (int i = 0; i < u.size(); ++i) 932bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 942bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double d = (u[i] - mean); 952bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double d2 = sqr(d); 962bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant var += d2; 972bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant skew += d * d2; 982bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis += d2 * d2; 992bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 1002bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant var /= u.size(); 1012bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double dev = std::sqrt(var); 1022bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant skew /= u.size() * dev * var; 1032bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis /= u.size() * var * var; 1042bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis -= 3; 1052bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_mean = std::exp(d.m() + sqr(d.s())/2); 1062bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_var = (std::exp(sqr(d.s())) - 1) * std::exp(2*d.m() + sqr(d.s())); 1070e20cae1a5be18fba591cd884aa2a389b66a3f49Howard Hinnant double x_skew = (std::exp(sqr(d.s())) + 2) * 1082bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant std::sqrt((std::exp(sqr(d.s())) - 1)); 1092bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_kurtosis = std::exp(4*sqr(d.s())) + 2*std::exp(3*sqr(d.s())) + 1102bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant 3*std::exp(2*sqr(d.s())) - 6; 111d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((mean - x_mean) / x_mean) < 0.01); 112d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((var - x_var) / x_var) < 0.01); 113d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((skew - x_skew) / x_skew) < 0.01); 114d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.03); 1152bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 1162bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 1172bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef std::lognormal_distribution<> D; 1182bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef D::param_type P; 1192bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef std::mt19937 G; 1202bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant G g; 1212bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant D d(-1./8, 0.5); 1222bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant const int N = 1000000; 1232bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant std::vector<D::result_type> u; 1242bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant for (int i = 0; i < N; ++i) 1252bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 1262bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant D::result_type v = d(g); 1272bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant assert(v > 0); 1282bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant u.push_back(v); 1292bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 1302bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double mean = std::accumulate(u.begin(), u.end(), 0.0) / u.size(); 1312bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double var = 0; 1322bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double skew = 0; 1332bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double kurtosis = 0; 1342bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant for (int i = 0; i < u.size(); ++i) 1352bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 1362bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double d = (u[i] - mean); 1372bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double d2 = sqr(d); 1382bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant var += d2; 1392bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant skew += d * d2; 1402bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis += d2 * d2; 1412bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 1422bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant var /= u.size(); 1432bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double dev = std::sqrt(var); 1442bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant skew /= u.size() * dev * var; 1452bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis /= u.size() * var * var; 1462bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis -= 3; 1472bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_mean = std::exp(d.m() + sqr(d.s())/2); 1482bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_var = (std::exp(sqr(d.s())) - 1) * std::exp(2*d.m() + sqr(d.s())); 1490e20cae1a5be18fba591cd884aa2a389b66a3f49Howard Hinnant double x_skew = (std::exp(sqr(d.s())) + 2) * 1502bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant std::sqrt((std::exp(sqr(d.s())) - 1)); 1512bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_kurtosis = std::exp(4*sqr(d.s())) + 2*std::exp(3*sqr(d.s())) + 1522bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant 3*std::exp(2*sqr(d.s())) - 6; 153d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((mean - x_mean) / x_mean) < 0.01); 154d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((var - x_var) / x_var) < 0.01); 155d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((skew - x_skew) / x_skew) < 0.02); 156d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.05); 1572bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 1582bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 1592bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef std::lognormal_distribution<> D; 1602bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef D::param_type P; 1612bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef std::mt19937 G; 1622bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant G g; 1632bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant D d; 1642bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant const int N = 1000000; 1652bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant std::vector<D::result_type> u; 1662bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant for (int i = 0; i < N; ++i) 1672bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 1682bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant D::result_type v = d(g); 1692bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant assert(v > 0); 1702bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant u.push_back(v); 1712bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 1722bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double mean = std::accumulate(u.begin(), u.end(), 0.0) / u.size(); 1732bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double var = 0; 1742bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double skew = 0; 1752bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double kurtosis = 0; 1762bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant for (int i = 0; i < u.size(); ++i) 1772bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 1782bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double d = (u[i] - mean); 1792bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double d2 = sqr(d); 1802bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant var += d2; 1812bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant skew += d * d2; 1822bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis += d2 * d2; 1832bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 1842bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant var /= u.size(); 1852bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double dev = std::sqrt(var); 1862bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant skew /= u.size() * dev * var; 1872bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis /= u.size() * var * var; 1882bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis -= 3; 1892bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_mean = std::exp(d.m() + sqr(d.s())/2); 1902bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_var = (std::exp(sqr(d.s())) - 1) * std::exp(2*d.m() + sqr(d.s())); 1910e20cae1a5be18fba591cd884aa2a389b66a3f49Howard Hinnant double x_skew = (std::exp(sqr(d.s())) + 2) * 1922bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant std::sqrt((std::exp(sqr(d.s())) - 1)); 1932bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_kurtosis = std::exp(4*sqr(d.s())) + 2*std::exp(3*sqr(d.s())) + 1942bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant 3*std::exp(2*sqr(d.s())) - 6; 195d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((mean - x_mean) / x_mean) < 0.01); 196d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((var - x_var) / x_var) < 0.02); 197d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((skew - x_skew) / x_skew) < 0.08); 198d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.4); 1992bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 2002bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 2012bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef std::lognormal_distribution<> D; 2022bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef D::param_type P; 2032bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant typedef std::mt19937 G; 2042bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant G g; 2052bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant D d(-0.78125, 1.25); 2062bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant const int N = 1000000; 2072bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant std::vector<D::result_type> u; 2082bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant for (int i = 0; i < N; ++i) 2092bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 2102bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant D::result_type v = d(g); 2112bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant assert(v > 0); 2122bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant u.push_back(v); 2132bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 2142bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double mean = std::accumulate(u.begin(), u.end(), 0.0) / u.size(); 2152bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double var = 0; 2162bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double skew = 0; 2172bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double kurtosis = 0; 2182bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant for (int i = 0; i < u.size(); ++i) 2192bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant { 2202bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double d = (u[i] - mean); 2212bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double d2 = sqr(d); 2222bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant var += d2; 2232bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant skew += d * d2; 2242bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis += d2 * d2; 2252bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 2262bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant var /= u.size(); 2272bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double dev = std::sqrt(var); 2282bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant skew /= u.size() * dev * var; 2292bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis /= u.size() * var * var; 2302bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant kurtosis -= 3; 2312bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_mean = std::exp(d.m() + sqr(d.s())/2); 2322bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_var = (std::exp(sqr(d.s())) - 1) * std::exp(2*d.m() + sqr(d.s())); 2330e20cae1a5be18fba591cd884aa2a389b66a3f49Howard Hinnant double x_skew = (std::exp(sqr(d.s())) + 2) * 2342bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant std::sqrt((std::exp(sqr(d.s())) - 1)); 2352bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant double x_kurtosis = std::exp(4*sqr(d.s())) + 2*std::exp(3*sqr(d.s())) + 2362bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant 3*std::exp(2*sqr(d.s())) - 6; 237d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((mean - x_mean) / x_mean) < 0.01); 238d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((var - x_var) / x_var) < 0.04); 239d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((skew - x_skew) / x_skew) < 0.2); 240d6d1171f2c3f254582ae1d5b9e14cea0ea8e701bHoward Hinnant assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.7); 2412bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant } 2422bc36fcff3de1ace5c74bb7c1459def41a67e862Howard Hinnant} 243