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 | #ifndef _LIBCPP___RANDOM_PIECEWISE_CONSTANT_DISTRIBUTION_H |
10 | #define _LIBCPP___RANDOM_PIECEWISE_CONSTANT_DISTRIBUTION_H |
11 | |
12 | #include <__algorithm/upper_bound.h> |
13 | #include <__config> |
14 | #include <__random/is_valid.h> |
15 | #include <__random/uniform_real_distribution.h> |
16 | #include <iosfwd> |
17 | #include <numeric> |
18 | #include <vector> |
19 | |
20 | #if !defined(_LIBCPP_HAS_NO_PRAGMA_SYSTEM_HEADER) |
21 | # pragma GCC system_header |
22 | #endif |
23 | |
24 | _LIBCPP_PUSH_MACROS |
25 | #include <__undef_macros> |
26 | |
27 | _LIBCPP_BEGIN_NAMESPACE_STD |
28 | |
29 | template <class _RealType = double> |
30 | class _LIBCPP_TEMPLATE_VIS piecewise_constant_distribution { |
31 | static_assert(__libcpp_random_is_valid_realtype<_RealType>::value, |
32 | "RealType must be a supported floating-point type" ); |
33 | |
34 | public: |
35 | // types |
36 | typedef _RealType result_type; |
37 | |
38 | class _LIBCPP_TEMPLATE_VIS param_type { |
39 | vector<result_type> __b_; |
40 | vector<result_type> __densities_; |
41 | vector<result_type> __areas_; |
42 | |
43 | public: |
44 | typedef piecewise_constant_distribution distribution_type; |
45 | |
46 | _LIBCPP_HIDE_FROM_ABI param_type(); |
47 | template <class _InputIteratorB, class _InputIteratorW> |
48 | _LIBCPP_HIDE_FROM_ABI param_type(_InputIteratorB __f_b, _InputIteratorB __l_b, _InputIteratorW __f_w); |
49 | #ifndef _LIBCPP_CXX03_LANG |
50 | template <class _UnaryOperation> |
51 | _LIBCPP_HIDE_FROM_ABI param_type(initializer_list<result_type> __bl, _UnaryOperation __fw); |
52 | #endif // _LIBCPP_CXX03_LANG |
53 | template <class _UnaryOperation> |
54 | _LIBCPP_HIDE_FROM_ABI param_type(size_t __nw, result_type __xmin, result_type __xmax, _UnaryOperation __fw); |
55 | _LIBCPP_HIDE_FROM_ABI param_type(param_type const&) = default; |
56 | _LIBCPP_HIDE_FROM_ABI param_type& operator=(const param_type& __rhs); |
57 | |
58 | _LIBCPP_HIDE_FROM_ABI vector<result_type> intervals() const { return __b_; } |
59 | _LIBCPP_HIDE_FROM_ABI vector<result_type> densities() const { return __densities_; } |
60 | |
61 | friend _LIBCPP_HIDE_FROM_ABI bool operator==(const param_type& __x, const param_type& __y) { |
62 | return __x.__densities_ == __y.__densities_ && __x.__b_ == __y.__b_; |
63 | } |
64 | friend _LIBCPP_HIDE_FROM_ABI bool operator!=(const param_type& __x, const param_type& __y) { return !(__x == __y); } |
65 | |
66 | private: |
67 | _LIBCPP_HIDE_FROM_ABI void __init(); |
68 | |
69 | friend class piecewise_constant_distribution; |
70 | |
71 | template <class _CharT, class _Traits, class _RT> |
72 | friend basic_ostream<_CharT, _Traits>& |
73 | operator<<(basic_ostream<_CharT, _Traits>& __os, const piecewise_constant_distribution<_RT>& __x); |
74 | |
75 | template <class _CharT, class _Traits, class _RT> |
76 | friend basic_istream<_CharT, _Traits>& |
77 | operator>>(basic_istream<_CharT, _Traits>& __is, piecewise_constant_distribution<_RT>& __x); |
78 | }; |
79 | |
80 | private: |
81 | param_type __p_; |
82 | |
83 | public: |
84 | // constructor and reset functions |
85 | _LIBCPP_HIDE_FROM_ABI piecewise_constant_distribution() {} |
86 | template <class _InputIteratorB, class _InputIteratorW> |
87 | _LIBCPP_HIDE_FROM_ABI |
88 | piecewise_constant_distribution(_InputIteratorB __f_b, _InputIteratorB __l_b, _InputIteratorW __f_w) |
89 | : __p_(__f_b, __l_b, __f_w) {} |
90 | |
91 | #ifndef _LIBCPP_CXX03_LANG |
92 | template <class _UnaryOperation> |
93 | _LIBCPP_HIDE_FROM_ABI piecewise_constant_distribution(initializer_list<result_type> __bl, _UnaryOperation __fw) |
94 | : __p_(__bl, __fw) {} |
95 | #endif // _LIBCPP_CXX03_LANG |
96 | |
97 | template <class _UnaryOperation> |
98 | _LIBCPP_HIDE_FROM_ABI |
99 | piecewise_constant_distribution(size_t __nw, result_type __xmin, result_type __xmax, _UnaryOperation __fw) |
100 | : __p_(__nw, __xmin, __xmax, __fw) {} |
101 | |
102 | _LIBCPP_HIDE_FROM_ABI explicit piecewise_constant_distribution(const param_type& __p) : __p_(__p) {} |
103 | |
104 | _LIBCPP_HIDE_FROM_ABI void reset() {} |
105 | |
106 | // generating functions |
107 | template <class _URNG> |
108 | _LIBCPP_HIDE_FROM_ABI result_type operator()(_URNG& __g) { |
109 | return (*this)(__g, __p_); |
110 | } |
111 | template <class _URNG> |
112 | _LIBCPP_HIDE_FROM_ABI result_type operator()(_URNG& __g, const param_type& __p); |
113 | |
114 | // property functions |
115 | _LIBCPP_HIDE_FROM_ABI vector<result_type> intervals() const { return __p_.intervals(); } |
116 | _LIBCPP_HIDE_FROM_ABI vector<result_type> densities() const { return __p_.densities(); } |
117 | |
118 | _LIBCPP_HIDE_FROM_ABI param_type param() const { return __p_; } |
119 | _LIBCPP_HIDE_FROM_ABI void param(const param_type& __p) { __p_ = __p; } |
120 | |
121 | _LIBCPP_HIDE_FROM_ABI result_type min() const { return __p_.__b_.front(); } |
122 | _LIBCPP_HIDE_FROM_ABI result_type max() const { return __p_.__b_.back(); } |
123 | |
124 | friend _LIBCPP_HIDE_FROM_ABI bool |
125 | operator==(const piecewise_constant_distribution& __x, const piecewise_constant_distribution& __y) { |
126 | return __x.__p_ == __y.__p_; |
127 | } |
128 | friend _LIBCPP_HIDE_FROM_ABI bool |
129 | operator!=(const piecewise_constant_distribution& __x, const piecewise_constant_distribution& __y) { |
130 | return !(__x == __y); |
131 | } |
132 | |
133 | template <class _CharT, class _Traits, class _RT> |
134 | friend basic_ostream<_CharT, _Traits>& |
135 | operator<<(basic_ostream<_CharT, _Traits>& __os, const piecewise_constant_distribution<_RT>& __x); |
136 | |
137 | template <class _CharT, class _Traits, class _RT> |
138 | friend basic_istream<_CharT, _Traits>& |
139 | operator>>(basic_istream<_CharT, _Traits>& __is, piecewise_constant_distribution<_RT>& __x); |
140 | }; |
141 | |
142 | template <class _RealType> |
143 | typename piecewise_constant_distribution<_RealType>::param_type& |
144 | piecewise_constant_distribution<_RealType>::param_type::operator=(const param_type& __rhs) { |
145 | // These can throw |
146 | __b_.reserve(__rhs.__b_.size()); |
147 | __densities_.reserve(__rhs.__densities_.size()); |
148 | __areas_.reserve(__rhs.__areas_.size()); |
149 | |
150 | // These can not throw |
151 | __b_ = __rhs.__b_; |
152 | __densities_ = __rhs.__densities_; |
153 | __areas_ = __rhs.__areas_; |
154 | return *this; |
155 | } |
156 | |
157 | template <class _RealType> |
158 | void piecewise_constant_distribution<_RealType>::param_type::__init() { |
159 | // __densities_ contains non-normalized areas |
160 | result_type __total_area = std::accumulate(__densities_.begin(), __densities_.end(), result_type()); |
161 | for (size_t __i = 0; __i < __densities_.size(); ++__i) |
162 | __densities_[__i] /= __total_area; |
163 | // __densities_ contains normalized areas |
164 | __areas_.assign(__densities_.size(), result_type()); |
165 | std::partial_sum(__densities_.begin(), __densities_.end() - 1, __areas_.begin() + 1); |
166 | // __areas_ contains partial sums of normalized areas: [0, __densities_ - 1] |
167 | __densities_.back() = 1 - __areas_.back(); // correct round off error |
168 | for (size_t __i = 0; __i < __densities_.size(); ++__i) |
169 | __densities_[__i] /= (__b_[__i + 1] - __b_[__i]); |
170 | // __densities_ now contains __densities_ |
171 | } |
172 | |
173 | template <class _RealType> |
174 | piecewise_constant_distribution<_RealType>::param_type::param_type() : __b_(2), __densities_(1, 1.0), __areas_(1, 0.0) { |
175 | __b_[1] = 1; |
176 | } |
177 | |
178 | template <class _RealType> |
179 | template <class _InputIteratorB, class _InputIteratorW> |
180 | piecewise_constant_distribution<_RealType>::param_type::param_type( |
181 | _InputIteratorB __f_b, _InputIteratorB __l_b, _InputIteratorW __f_w) |
182 | : __b_(__f_b, __l_b) { |
183 | if (__b_.size() < 2) { |
184 | __b_.resize(2); |
185 | __b_[0] = 0; |
186 | __b_[1] = 1; |
187 | __densities_.assign(1, 1.0); |
188 | __areas_.assign(1, 0.0); |
189 | } else { |
190 | __densities_.reserve(__b_.size() - 1); |
191 | for (size_t __i = 0; __i < __b_.size() - 1; ++__i, ++__f_w) |
192 | __densities_.push_back(*__f_w); |
193 | __init(); |
194 | } |
195 | } |
196 | |
197 | #ifndef _LIBCPP_CXX03_LANG |
198 | |
199 | template <class _RealType> |
200 | template <class _UnaryOperation> |
201 | piecewise_constant_distribution<_RealType>::param_type::param_type( |
202 | initializer_list<result_type> __bl, _UnaryOperation __fw) |
203 | : __b_(__bl.begin(), __bl.end()) { |
204 | if (__b_.size() < 2) { |
205 | __b_.resize(2); |
206 | __b_[0] = 0; |
207 | __b_[1] = 1; |
208 | __densities_.assign(1, 1.0); |
209 | __areas_.assign(1, 0.0); |
210 | } else { |
211 | __densities_.reserve(__b_.size() - 1); |
212 | for (size_t __i = 0; __i < __b_.size() - 1; ++__i) |
213 | __densities_.push_back(__fw((__b_[__i + 1] + __b_[__i]) * .5)); |
214 | __init(); |
215 | } |
216 | } |
217 | |
218 | #endif // _LIBCPP_CXX03_LANG |
219 | |
220 | template <class _RealType> |
221 | template <class _UnaryOperation> |
222 | piecewise_constant_distribution<_RealType>::param_type::param_type( |
223 | size_t __nw, result_type __xmin, result_type __xmax, _UnaryOperation __fw) |
224 | : __b_(__nw == 0 ? 2 : __nw + 1) { |
225 | size_t __n = __b_.size() - 1; |
226 | result_type __d = (__xmax - __xmin) / __n; |
227 | __densities_.reserve(__n); |
228 | for (size_t __i = 0; __i < __n; ++__i) { |
229 | __b_[__i] = __xmin + __i * __d; |
230 | __densities_.push_back(__fw(__b_[__i] + __d * .5)); |
231 | } |
232 | __b_[__n] = __xmax; |
233 | __init(); |
234 | } |
235 | |
236 | template <class _RealType> |
237 | template <class _URNG> |
238 | _RealType piecewise_constant_distribution<_RealType>::operator()(_URNG& __g, const param_type& __p) { |
239 | static_assert(__libcpp_random_is_valid_urng<_URNG>::value, "" ); |
240 | typedef uniform_real_distribution<result_type> _Gen; |
241 | result_type __u = _Gen()(__g); |
242 | ptrdiff_t __k = std::upper_bound(__p.__areas_.begin(), __p.__areas_.end(), __u) - __p.__areas_.begin() - 1; |
243 | return (__u - __p.__areas_[__k]) / __p.__densities_[__k] + __p.__b_[__k]; |
244 | } |
245 | |
246 | template <class _CharT, class _Traits, class _RT> |
247 | _LIBCPP_HIDE_FROM_ABI basic_ostream<_CharT, _Traits>& |
248 | operator<<(basic_ostream<_CharT, _Traits>& __os, const piecewise_constant_distribution<_RT>& __x) { |
249 | __save_flags<_CharT, _Traits> __lx(__os); |
250 | typedef basic_ostream<_CharT, _Traits> _OStream; |
251 | __os.flags(_OStream::dec | _OStream::left | _OStream::fixed | _OStream::scientific); |
252 | _CharT __sp = __os.widen(' '); |
253 | __os.fill(__sp); |
254 | size_t __n = __x.__p_.__b_.size(); |
255 | __os << __n; |
256 | for (size_t __i = 0; __i < __n; ++__i) |
257 | __os << __sp << __x.__p_.__b_[__i]; |
258 | __n = __x.__p_.__densities_.size(); |
259 | __os << __sp << __n; |
260 | for (size_t __i = 0; __i < __n; ++__i) |
261 | __os << __sp << __x.__p_.__densities_[__i]; |
262 | __n = __x.__p_.__areas_.size(); |
263 | __os << __sp << __n; |
264 | for (size_t __i = 0; __i < __n; ++__i) |
265 | __os << __sp << __x.__p_.__areas_[__i]; |
266 | return __os; |
267 | } |
268 | |
269 | template <class _CharT, class _Traits, class _RT> |
270 | _LIBCPP_HIDE_FROM_ABI basic_istream<_CharT, _Traits>& |
271 | operator>>(basic_istream<_CharT, _Traits>& __is, piecewise_constant_distribution<_RT>& __x) { |
272 | typedef piecewise_constant_distribution<_RT> _Eng; |
273 | typedef typename _Eng::result_type result_type; |
274 | __save_flags<_CharT, _Traits> __lx(__is); |
275 | typedef basic_istream<_CharT, _Traits> _Istream; |
276 | __is.flags(_Istream::dec | _Istream::skipws); |
277 | size_t __n; |
278 | __is >> __n; |
279 | vector<result_type> __b(__n); |
280 | for (size_t __i = 0; __i < __n; ++__i) |
281 | __is >> __b[__i]; |
282 | __is >> __n; |
283 | vector<result_type> __densities(__n); |
284 | for (size_t __i = 0; __i < __n; ++__i) |
285 | __is >> __densities[__i]; |
286 | __is >> __n; |
287 | vector<result_type> __areas(__n); |
288 | for (size_t __i = 0; __i < __n; ++__i) |
289 | __is >> __areas[__i]; |
290 | if (!__is.fail()) { |
291 | swap(__x.__p_.__b_, __b); |
292 | swap(__x.__p_.__densities_, __densities); |
293 | swap(__x.__p_.__areas_, __areas); |
294 | } |
295 | return __is; |
296 | } |
297 | |
298 | _LIBCPP_END_NAMESPACE_STD |
299 | |
300 | _LIBCPP_POP_MACROS |
301 | |
302 | #endif // _LIBCPP___RANDOM_PIECEWISE_CONSTANT_DISTRIBUTION_H |
303 | |