1//===----------------------------------------------------------------------===//
2//
3// Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions.
4// See https://llvm.org/LICENSE.txt for license information.
5// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
6//
7//===----------------------------------------------------------------------===//
8//
9// REQUIRES: long_tests
10
11// <random>
12
13// template<class RealType = double>
14// class weibull_distribution
15
16// template<class _URNG> result_type operator()(_URNG& g, const param_type& parm);
17
18#include <random>
19#include <cassert>
20#include <vector>
21#include <numeric>
22#include <cstddef>
23
24#include "test_macros.h"
25
26template <class T>
27inline
28T
29sqr(T x)
30{
31 return x * x;
32}
33
34int main(int, char**)
35{
36 {
37 typedef std::weibull_distribution<> D;
38 typedef D::param_type P;
39 typedef std::mt19937 G;
40 G g;
41 D d(0.5, 2);
42 P p(1, .5);
43 const int N = 1000000;
44 std::vector<D::result_type> u;
45 for (int i = 0; i < N; ++i)
46 {
47 D::result_type v = d(g, p);
48 assert(d.min() <= v);
49 u.push_back(x: v);
50 }
51 double mean = std::accumulate(u.begin(), u.end(), 0.0) / u.size();
52 double var = 0;
53 double skew = 0;
54 double kurtosis = 0;
55 for (std::size_t i = 0; i < u.size(); ++i)
56 {
57 double dbl = (u[i] - mean);
58 double d2 = sqr(dbl);
59 var += d2;
60 skew += dbl * d2;
61 kurtosis += d2 * d2;
62 }
63 var /= u.size();
64 double dev = std::sqrt(x: var);
65 skew /= u.size() * dev * var;
66 kurtosis /= u.size() * var * var;
67 kurtosis -= 3;
68 double x_mean = p.b() * std::tgamma(1 + 1/p.a());
69 double x_var = sqr(p.b()) * std::tgamma(1 + 2/p.a()) - sqr(x_mean);
70 double x_skew = (sqr(p.b())*p.b() * std::tgamma(1 + 3/p.a()) -
71 3*x_mean*x_var - sqr(x_mean)*x_mean) /
72 (std::sqrt(x: x_var)*x_var);
73 double x_kurtosis = (sqr(sqr(p.b())) * std::tgamma(1 + 4/p.a()) -
74 4*x_skew*x_var*sqrt(x: x_var)*x_mean -
75 6*sqr(x_mean)*x_var - sqr(sqr(x_mean))) / sqr(x_var) - 3;
76 assert(std::abs((mean - x_mean) / x_mean) < 0.01);
77 assert(std::abs((var - x_var) / x_var) < 0.01);
78 assert(std::abs((skew - x_skew) / x_skew) < 0.01);
79 assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.01);
80 }
81 {
82 typedef std::weibull_distribution<> D;
83 typedef D::param_type P;
84 typedef std::mt19937 G;
85 G g;
86 D d(1, .5);
87 P p(2, 3);
88 const int N = 1000000;
89 std::vector<D::result_type> u;
90 for (int i = 0; i < N; ++i)
91 {
92 D::result_type v = d(g, p);
93 assert(d.min() <= v);
94 u.push_back(v);
95 }
96 double mean = std::accumulate(u.begin(), u.end(), 0.0) / u.size();
97 double var = 0;
98 double skew = 0;
99 double kurtosis = 0;
100 for (std::size_t i = 0; i < u.size(); ++i)
101 {
102 double dbl = (u[i] - mean);
103 double d2 = sqr(dbl);
104 var += d2;
105 skew += dbl * d2;
106 kurtosis += d2 * d2;
107 }
108 var /= u.size();
109 double dev = std::sqrt(x: var);
110 skew /= u.size() * dev * var;
111 kurtosis /= u.size() * var * var;
112 kurtosis -= 3;
113 double x_mean = p.b() * std::tgamma(1 + 1/p.a());
114 double x_var = sqr(p.b()) * std::tgamma(1 + 2/p.a()) - sqr(x_mean);
115 double x_skew = (sqr(p.b())*p.b() * std::tgamma(1 + 3/p.a()) -
116 3*x_mean*x_var - sqr(x_mean)*x_mean) /
117 (std::sqrt(x: x_var)*x_var);
118 double x_kurtosis = (sqr(sqr(p.b())) * std::tgamma(1 + 4/p.a()) -
119 4*x_skew*x_var*sqrt(x: x_var)*x_mean -
120 6*sqr(x_mean)*x_var - sqr(sqr(x_mean))) / sqr(x_var) - 3;
121 assert(std::abs((mean - x_mean) / x_mean) < 0.01);
122 assert(std::abs((var - x_var) / x_var) < 0.01);
123 assert(std::abs((skew - x_skew) / x_skew) < 0.01);
124 assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.03);
125 }
126 {
127 typedef std::weibull_distribution<> D;
128 typedef D::param_type P;
129 typedef std::mt19937 G;
130 G g;
131 D d(2, 3);
132 P p(.5, 2);
133 const int N = 1000000;
134 std::vector<D::result_type> u;
135 for (int i = 0; i < N; ++i)
136 {
137 D::result_type v = d(g, p);
138 assert(d.min() <= v);
139 u.push_back(v);
140 }
141 double mean = std::accumulate(u.begin(), u.end(), 0.0) / u.size();
142 double var = 0;
143 double skew = 0;
144 double kurtosis = 0;
145 for (std::size_t i = 0; i < u.size(); ++i)
146 {
147 double dbl = (u[i] - mean);
148 double d2 = sqr(dbl);
149 var += d2;
150 skew += dbl * d2;
151 kurtosis += d2 * d2;
152 }
153 var /= u.size();
154 double dev = std::sqrt(x: var);
155 skew /= u.size() * dev * var;
156 kurtosis /= u.size() * var * var;
157 kurtosis -= 3;
158 double x_mean = p.b() * std::tgamma(1 + 1/p.a());
159 double x_var = sqr(p.b()) * std::tgamma(1 + 2/p.a()) - sqr(x_mean);
160 double x_skew = (sqr(p.b())*p.b() * std::tgamma(1 + 3/p.a()) -
161 3*x_mean*x_var - sqr(x_mean)*x_mean) /
162 (std::sqrt(x: x_var)*x_var);
163 double x_kurtosis = (sqr(sqr(p.b())) * std::tgamma(1 + 4/p.a()) -
164 4*x_skew*x_var*sqrt(x: x_var)*x_mean -
165 6*sqr(x_mean)*x_var - sqr(sqr(x_mean))) / sqr(x_var) - 3;
166 assert(std::abs((mean - x_mean) / x_mean) < 0.01);
167 assert(std::abs((var - x_var) / x_var) < 0.01);
168 assert(std::abs((skew - x_skew) / x_skew) < 0.01);
169 assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.03);
170 }
171
172 return 0;
173}
174

source code of libcxx/test/std/numerics/rand/rand.dist/rand.dist.pois/rand.dist.pois.weibull/eval_param.pass.cpp