1 //===----------------------------------------------------------------------===//
2 //
3 // The LLVM Compiler Infrastructure
4 //
5 // This file is dual licensed under the MIT and the University of Illinois Open
6 // Source Licenses. See LICENSE.TXT for details.
7 //
8 //===----------------------------------------------------------------------===//
9
10 // <random>
11
12 // template<class RealType = double>
13 // class lognormal_distribution
14
15 // template<class _URNG> result_type operator()(_URNG& g, const param_type& parm);
16
17 #include <random>
18 #include <cassert>
19 #include <vector>
20 #include <numeric>
21
22 template <class T>
23 inline
24 T
sqr(T x)25 sqr(T x)
26 {
27 return x * x;
28 }
29
main()30 int main()
31 {
32
33 {
34 typedef std::lognormal_distribution<> D;
35 typedef D::param_type P;
36 typedef std::mt19937 G;
37 G g;
38 D d;
39 P p(-1./8192, 0.015625);
40 const int N = 1000000;
41 std::vector<D::result_type> u;
42 for (int i = 0; i < N; ++i)
43 {
44 D::result_type v = d(g, p);
45 assert(v > 0);
46 u.push_back(v);
47 }
48 double mean = std::accumulate(u.begin(), u.end(), 0.0) / u.size();
49 double var = 0;
50 double skew = 0;
51 double kurtosis = 0;
52 for (int i = 0; i < u.size(); ++i)
53 {
54 double d = (u[i] - mean);
55 double d2 = sqr(d);
56 var += d2;
57 skew += d * d2;
58 kurtosis += d2 * d2;
59 }
60 var /= u.size();
61 double dev = std::sqrt(var);
62 skew /= u.size() * dev * var;
63 kurtosis /= u.size() * var * var;
64 kurtosis -= 3;
65 double x_mean = std::exp(p.m() + sqr(p.s())/2);
66 double x_var = (std::exp(sqr(p.s())) - 1) * std::exp(2*p.m() + sqr(p.s()));
67 double x_skew = (std::exp(sqr(p.s())) + 2) *
68 std::sqrt((std::exp(sqr(p.s())) - 1));
69 double x_kurtosis = std::exp(4*sqr(p.s())) + 2*std::exp(3*sqr(p.s())) +
70 3*std::exp(2*sqr(p.s())) - 6;
71 assert(std::abs((mean - x_mean) / x_mean) < 0.01);
72 assert(std::abs((var - x_var) / x_var) < 0.01);
73 assert(std::abs((skew - x_skew) / x_skew) < 0.05);
74 assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.25);
75 }
76 {
77 typedef std::lognormal_distribution<> D;
78 typedef D::param_type P;
79 typedef std::mt19937 G;
80 G g;
81 D d;
82 P p(-1./32, 0.25);
83 const int N = 1000000;
84 std::vector<D::result_type> u;
85 for (int i = 0; i < N; ++i)
86 {
87 D::result_type v = d(g, p);
88 assert(v > 0);
89 u.push_back(v);
90 }
91 double mean = std::accumulate(u.begin(), u.end(), 0.0) / u.size();
92 double var = 0;
93 double skew = 0;
94 double kurtosis = 0;
95 for (int i = 0; i < u.size(); ++i)
96 {
97 double d = (u[i] - mean);
98 double d2 = sqr(d);
99 var += d2;
100 skew += d * d2;
101 kurtosis += d2 * d2;
102 }
103 var /= u.size();
104 double dev = std::sqrt(var);
105 skew /= u.size() * dev * var;
106 kurtosis /= u.size() * var * var;
107 kurtosis -= 3;
108 double x_mean = std::exp(p.m() + sqr(p.s())/2);
109 double x_var = (std::exp(sqr(p.s())) - 1) * std::exp(2*p.m() + sqr(p.s()));
110 double x_skew = (std::exp(sqr(p.s())) + 2) *
111 std::sqrt((std::exp(sqr(p.s())) - 1));
112 double x_kurtosis = std::exp(4*sqr(p.s())) + 2*std::exp(3*sqr(p.s())) +
113 3*std::exp(2*sqr(p.s())) - 6;
114 assert(std::abs((mean - x_mean) / x_mean) < 0.01);
115 assert(std::abs((var - x_var) / x_var) < 0.01);
116 assert(std::abs((skew - x_skew) / x_skew) < 0.01);
117 assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.03);
118 }
119 {
120 typedef std::lognormal_distribution<> D;
121 typedef D::param_type P;
122 typedef std::mt19937 G;
123 G g;
124 D d;
125 P p(-1./8, 0.5);
126 const int N = 1000000;
127 std::vector<D::result_type> u;
128 for (int i = 0; i < N; ++i)
129 {
130 D::result_type v = d(g, p);
131 assert(v > 0);
132 u.push_back(v);
133 }
134 double mean = std::accumulate(u.begin(), u.end(), 0.0) / u.size();
135 double var = 0;
136 double skew = 0;
137 double kurtosis = 0;
138 for (int i = 0; i < u.size(); ++i)
139 {
140 double d = (u[i] - mean);
141 double d2 = sqr(d);
142 var += d2;
143 skew += d * d2;
144 kurtosis += d2 * d2;
145 }
146 var /= u.size();
147 double dev = std::sqrt(var);
148 skew /= u.size() * dev * var;
149 kurtosis /= u.size() * var * var;
150 kurtosis -= 3;
151 double x_mean = std::exp(p.m() + sqr(p.s())/2);
152 double x_var = (std::exp(sqr(p.s())) - 1) * std::exp(2*p.m() + sqr(p.s()));
153 double x_skew = (std::exp(sqr(p.s())) + 2) *
154 std::sqrt((std::exp(sqr(p.s())) - 1));
155 double x_kurtosis = std::exp(4*sqr(p.s())) + 2*std::exp(3*sqr(p.s())) +
156 3*std::exp(2*sqr(p.s())) - 6;
157 assert(std::abs((mean - x_mean) / x_mean) < 0.01);
158 assert(std::abs((var - x_var) / x_var) < 0.01);
159 assert(std::abs((skew - x_skew) / x_skew) < 0.02);
160 assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.05);
161 }
162 {
163 typedef std::lognormal_distribution<> D;
164 typedef D::param_type P;
165 typedef std::mt19937 G;
166 G g;
167 D d(3, 4);
168 P p;
169 const int N = 1000000;
170 std::vector<D::result_type> u;
171 for (int i = 0; i < N; ++i)
172 {
173 D::result_type v = d(g, p);
174 assert(v > 0);
175 u.push_back(v);
176 }
177 double mean = std::accumulate(u.begin(), u.end(), 0.0) / u.size();
178 double var = 0;
179 double skew = 0;
180 double kurtosis = 0;
181 for (int i = 0; i < u.size(); ++i)
182 {
183 double d = (u[i] - mean);
184 double d2 = sqr(d);
185 var += d2;
186 skew += d * d2;
187 kurtosis += d2 * d2;
188 }
189 var /= u.size();
190 double dev = std::sqrt(var);
191 skew /= u.size() * dev * var;
192 kurtosis /= u.size() * var * var;
193 kurtosis -= 3;
194 double x_mean = std::exp(p.m() + sqr(p.s())/2);
195 double x_var = (std::exp(sqr(p.s())) - 1) * std::exp(2*p.m() + sqr(p.s()));
196 double x_skew = (std::exp(sqr(p.s())) + 2) *
197 std::sqrt((std::exp(sqr(p.s())) - 1));
198 double x_kurtosis = std::exp(4*sqr(p.s())) + 2*std::exp(3*sqr(p.s())) +
199 3*std::exp(2*sqr(p.s())) - 6;
200 assert(std::abs((mean - x_mean) / x_mean) < 0.01);
201 assert(std::abs((var - x_var) / x_var) < 0.02);
202 assert(std::abs((skew - x_skew) / x_skew) < 0.08);
203 assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.4);
204 }
205 {
206 typedef std::lognormal_distribution<> D;
207 typedef D::param_type P;
208 typedef std::mt19937 G;
209 G g;
210 D d;
211 P p(-0.78125, 1.25);
212 const int N = 1000000;
213 std::vector<D::result_type> u;
214 for (int i = 0; i < N; ++i)
215 {
216 D::result_type v = d(g, p);
217 assert(v > 0);
218 u.push_back(v);
219 }
220 double mean = std::accumulate(u.begin(), u.end(), 0.0) / u.size();
221 double var = 0;
222 double skew = 0;
223 double kurtosis = 0;
224 for (int i = 0; i < u.size(); ++i)
225 {
226 double d = (u[i] - mean);
227 double d2 = sqr(d);
228 var += d2;
229 skew += d * d2;
230 kurtosis += d2 * d2;
231 }
232 var /= u.size();
233 double dev = std::sqrt(var);
234 skew /= u.size() * dev * var;
235 kurtosis /= u.size() * var * var;
236 kurtosis -= 3;
237 double x_mean = std::exp(p.m() + sqr(p.s())/2);
238 double x_var = (std::exp(sqr(p.s())) - 1) * std::exp(2*p.m() + sqr(p.s()));
239 double x_skew = (std::exp(sqr(p.s())) + 2) *
240 std::sqrt((std::exp(sqr(p.s())) - 1));
241 double x_kurtosis = std::exp(4*sqr(p.s())) + 2*std::exp(3*sqr(p.s())) +
242 3*std::exp(2*sqr(p.s())) - 6;
243 assert(std::abs((mean - x_mean) / x_mean) < 0.01);
244 assert(std::abs((var - x_var) / x_var) < 0.04);
245 assert(std::abs((skew - x_skew) / x_skew) < 0.2);
246 assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.7);
247 }
248 }
249