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_LINEAR_DISTRIBUTION_H 10 #define _LIBCPP___RANDOM_PIECEWISE_LINEAR_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_linear_distribution 31 { 32 public: 33 // types 34 typedef _RealType result_type; 35 36 class _LIBCPP_TEMPLATE_VIS param_type 37 { 38 vector<result_type> __b_; 39 vector<result_type> __densities_; 40 vector<result_type> __areas_; 41 public: 42 typedef piecewise_linear_distribution distribution_type; 43 44 param_type(); 45 template<class _InputIteratorB, class _InputIteratorW> 46 param_type(_InputIteratorB __f_b, _InputIteratorB __l_b, 47 _InputIteratorW __f_w); 48 #ifndef _LIBCPP_CXX03_LANG 49 template<class _UnaryOperation> 50 param_type(initializer_list<result_type> __bl, _UnaryOperation __fw); 51 #endif // _LIBCPP_CXX03_LANG 52 template<class _UnaryOperation> 53 param_type(size_t __nw, result_type __xmin, result_type __xmax, 54 _UnaryOperation __fw); 55 param_type(param_type const&) = default; 56 param_type & operator=(const param_type& __rhs); 57 58 _LIBCPP_INLINE_VISIBILITY 59 vector<result_type> intervals() const {return __b_;} 60 _LIBCPP_INLINE_VISIBILITY 61 vector<result_type> densities() const {return __densities_;} 62 63 friend _LIBCPP_INLINE_VISIBILITY 64 bool operator==(const param_type& __x, const param_type& __y) 65 {return __x.__densities_ == __y.__densities_ && __x.__b_ == __y.__b_;} 66 friend _LIBCPP_INLINE_VISIBILITY 67 bool operator!=(const param_type& __x, const param_type& __y) 68 {return !(__x == __y);} 69 70 private: 71 void __init(); 72 73 friend class piecewise_linear_distribution; 74 75 template <class _CharT, class _Traits, class _RT> 76 friend 77 basic_ostream<_CharT, _Traits>& 78 operator<<(basic_ostream<_CharT, _Traits>& __os, 79 const piecewise_linear_distribution<_RT>& __x); 80 81 template <class _CharT, class _Traits, class _RT> 82 friend 83 basic_istream<_CharT, _Traits>& 84 operator>>(basic_istream<_CharT, _Traits>& __is, 85 piecewise_linear_distribution<_RT>& __x); 86 }; 87 88 private: 89 param_type __p_; 90 91 public: 92 // constructor and reset functions 93 _LIBCPP_INLINE_VISIBILITY 94 piecewise_linear_distribution() {} 95 template<class _InputIteratorB, class _InputIteratorW> 96 _LIBCPP_INLINE_VISIBILITY 97 piecewise_linear_distribution(_InputIteratorB __f_b, 98 _InputIteratorB __l_b, 99 _InputIteratorW __f_w) 100 : __p_(__f_b, __l_b, __f_w) {} 101 102 #ifndef _LIBCPP_CXX03_LANG 103 template<class _UnaryOperation> 104 _LIBCPP_INLINE_VISIBILITY 105 piecewise_linear_distribution(initializer_list<result_type> __bl, 106 _UnaryOperation __fw) 107 : __p_(__bl, __fw) {} 108 #endif // _LIBCPP_CXX03_LANG 109 110 template<class _UnaryOperation> 111 _LIBCPP_INLINE_VISIBILITY 112 piecewise_linear_distribution(size_t __nw, result_type __xmin, 113 result_type __xmax, _UnaryOperation __fw) 114 : __p_(__nw, __xmin, __xmax, __fw) {} 115 116 _LIBCPP_INLINE_VISIBILITY 117 explicit piecewise_linear_distribution(const param_type& __p) 118 : __p_(__p) {} 119 120 _LIBCPP_INLINE_VISIBILITY 121 void reset() {} 122 123 // generating functions 124 template<class _URNG> 125 _LIBCPP_INLINE_VISIBILITY 126 result_type operator()(_URNG& __g) 127 {return (*this)(__g, __p_);} 128 template<class _URNG> result_type operator()(_URNG& __g, const param_type& __p); 129 130 // property functions 131 _LIBCPP_INLINE_VISIBILITY 132 vector<result_type> intervals() const {return __p_.intervals();} 133 _LIBCPP_INLINE_VISIBILITY 134 vector<result_type> densities() const {return __p_.densities();} 135 136 _LIBCPP_INLINE_VISIBILITY 137 param_type param() const {return __p_;} 138 _LIBCPP_INLINE_VISIBILITY 139 void param(const param_type& __p) {__p_ = __p;} 140 141 _LIBCPP_INLINE_VISIBILITY 142 result_type min() const {return __p_.__b_.front();} 143 _LIBCPP_INLINE_VISIBILITY 144 result_type max() const {return __p_.__b_.back();} 145 146 friend _LIBCPP_INLINE_VISIBILITY 147 bool operator==(const piecewise_linear_distribution& __x, 148 const piecewise_linear_distribution& __y) 149 {return __x.__p_ == __y.__p_;} 150 friend _LIBCPP_INLINE_VISIBILITY 151 bool operator!=(const piecewise_linear_distribution& __x, 152 const piecewise_linear_distribution& __y) 153 {return !(__x == __y);} 154 155 template <class _CharT, class _Traits, class _RT> 156 friend 157 basic_ostream<_CharT, _Traits>& 158 operator<<(basic_ostream<_CharT, _Traits>& __os, 159 const piecewise_linear_distribution<_RT>& __x); 160 161 template <class _CharT, class _Traits, class _RT> 162 friend 163 basic_istream<_CharT, _Traits>& 164 operator>>(basic_istream<_CharT, _Traits>& __is, 165 piecewise_linear_distribution<_RT>& __x); 166 }; 167 168 template<class _RealType> 169 typename piecewise_linear_distribution<_RealType>::param_type & 170 piecewise_linear_distribution<_RealType>::param_type::operator= 171 (const param_type& __rhs) 172 { 173 // These can throw 174 __b_.reserve (__rhs.__b_.size ()); 175 __densities_.reserve(__rhs.__densities_.size()); 176 __areas_.reserve (__rhs.__areas_.size()); 177 178 // These can not throw 179 __b_ = __rhs.__b_; 180 __densities_ = __rhs.__densities_; 181 __areas_ = __rhs.__areas_; 182 return *this; 183 } 184 185 186 template<class _RealType> 187 void 188 piecewise_linear_distribution<_RealType>::param_type::__init() 189 { 190 __areas_.assign(__densities_.size() - 1, result_type()); 191 result_type _Sp = 0; 192 for (size_t __i = 0; __i < __areas_.size(); ++__i) 193 { 194 __areas_[__i] = (__densities_[__i+1] + __densities_[__i]) * 195 (__b_[__i+1] - __b_[__i]) * .5; 196 _Sp += __areas_[__i]; 197 } 198 for (size_t __i = __areas_.size(); __i > 1;) 199 { 200 --__i; 201 __areas_[__i] = __areas_[__i-1] / _Sp; 202 } 203 __areas_[0] = 0; 204 for (size_t __i = 1; __i < __areas_.size(); ++__i) 205 __areas_[__i] += __areas_[__i-1]; 206 for (size_t __i = 0; __i < __densities_.size(); ++__i) 207 __densities_[__i] /= _Sp; 208 } 209 210 template<class _RealType> 211 piecewise_linear_distribution<_RealType>::param_type::param_type() 212 : __b_(2), 213 __densities_(2, 1.0), 214 __areas_(1, 0.0) 215 { 216 __b_[1] = 1; 217 } 218 219 template<class _RealType> 220 template<class _InputIteratorB, class _InputIteratorW> 221 piecewise_linear_distribution<_RealType>::param_type::param_type( 222 _InputIteratorB __f_b, _InputIteratorB __l_b, _InputIteratorW __f_w) 223 : __b_(__f_b, __l_b) 224 { 225 if (__b_.size() < 2) 226 { 227 __b_.resize(2); 228 __b_[0] = 0; 229 __b_[1] = 1; 230 __densities_.assign(2, 1.0); 231 __areas_.assign(1, 0.0); 232 } 233 else 234 { 235 __densities_.reserve(__b_.size()); 236 for (size_t __i = 0; __i < __b_.size(); ++__i, ++__f_w) 237 __densities_.push_back(*__f_w); 238 __init(); 239 } 240 } 241 242 #ifndef _LIBCPP_CXX03_LANG 243 244 template<class _RealType> 245 template<class _UnaryOperation> 246 piecewise_linear_distribution<_RealType>::param_type::param_type( 247 initializer_list<result_type> __bl, _UnaryOperation __fw) 248 : __b_(__bl.begin(), __bl.end()) 249 { 250 if (__b_.size() < 2) 251 { 252 __b_.resize(2); 253 __b_[0] = 0; 254 __b_[1] = 1; 255 __densities_.assign(2, 1.0); 256 __areas_.assign(1, 0.0); 257 } 258 else 259 { 260 __densities_.reserve(__b_.size()); 261 for (size_t __i = 0; __i < __b_.size(); ++__i) 262 __densities_.push_back(__fw(__b_[__i])); 263 __init(); 264 } 265 } 266 267 #endif // _LIBCPP_CXX03_LANG 268 269 template<class _RealType> 270 template<class _UnaryOperation> 271 piecewise_linear_distribution<_RealType>::param_type::param_type( 272 size_t __nw, result_type __xmin, result_type __xmax, _UnaryOperation __fw) 273 : __b_(__nw == 0 ? 2 : __nw + 1) 274 { 275 size_t __n = __b_.size() - 1; 276 result_type __d = (__xmax - __xmin) / __n; 277 __densities_.reserve(__b_.size()); 278 for (size_t __i = 0; __i < __n; ++__i) 279 { 280 __b_[__i] = __xmin + __i * __d; 281 __densities_.push_back(__fw(__b_[__i])); 282 } 283 __b_[__n] = __xmax; 284 __densities_.push_back(__fw(__b_[__n])); 285 __init(); 286 } 287 288 template<class _RealType> 289 template<class _URNG> 290 _RealType 291 piecewise_linear_distribution<_RealType>::operator()(_URNG& __g, const param_type& __p) 292 { 293 static_assert(__libcpp_random_is_valid_urng<_URNG>::value, ""); 294 typedef uniform_real_distribution<result_type> _Gen; 295 result_type __u = _Gen()(__g); 296 ptrdiff_t __k = _VSTD::upper_bound(__p.__areas_.begin(), __p.__areas_.end(), 297 __u) - __p.__areas_.begin() - 1; 298 __u -= __p.__areas_[__k]; 299 const result_type __dk = __p.__densities_[__k]; 300 const result_type __dk1 = __p.__densities_[__k+1]; 301 const result_type __deltad = __dk1 - __dk; 302 const result_type __bk = __p.__b_[__k]; 303 if (__deltad == 0) 304 return __u / __dk + __bk; 305 const result_type __bk1 = __p.__b_[__k+1]; 306 const result_type __deltab = __bk1 - __bk; 307 return (__bk * __dk1 - __bk1 * __dk + 308 _VSTD::sqrt(__deltab * (__deltab * __dk * __dk + 2 * __deltad * __u))) / 309 __deltad; 310 } 311 312 template <class _CharT, class _Traits, class _RT> 313 basic_ostream<_CharT, _Traits>& 314 operator<<(basic_ostream<_CharT, _Traits>& __os, 315 const piecewise_linear_distribution<_RT>& __x) 316 { 317 __save_flags<_CharT, _Traits> __lx(__os); 318 typedef basic_ostream<_CharT, _Traits> _OStream; 319 __os.flags(_OStream::dec | _OStream::left | _OStream::fixed | 320 _OStream::scientific); 321 _CharT __sp = __os.widen(' '); 322 __os.fill(__sp); 323 size_t __n = __x.__p_.__b_.size(); 324 __os << __n; 325 for (size_t __i = 0; __i < __n; ++__i) 326 __os << __sp << __x.__p_.__b_[__i]; 327 __n = __x.__p_.__densities_.size(); 328 __os << __sp << __n; 329 for (size_t __i = 0; __i < __n; ++__i) 330 __os << __sp << __x.__p_.__densities_[__i]; 331 __n = __x.__p_.__areas_.size(); 332 __os << __sp << __n; 333 for (size_t __i = 0; __i < __n; ++__i) 334 __os << __sp << __x.__p_.__areas_[__i]; 335 return __os; 336 } 337 338 template <class _CharT, class _Traits, class _RT> 339 basic_istream<_CharT, _Traits>& 340 operator>>(basic_istream<_CharT, _Traits>& __is, 341 piecewise_linear_distribution<_RT>& __x) 342 { 343 typedef piecewise_linear_distribution<_RT> _Eng; 344 typedef typename _Eng::result_type result_type; 345 __save_flags<_CharT, _Traits> __lx(__is); 346 typedef basic_istream<_CharT, _Traits> _Istream; 347 __is.flags(_Istream::dec | _Istream::skipws); 348 size_t __n; 349 __is >> __n; 350 vector<result_type> __b(__n); 351 for (size_t __i = 0; __i < __n; ++__i) 352 __is >> __b[__i]; 353 __is >> __n; 354 vector<result_type> __densities(__n); 355 for (size_t __i = 0; __i < __n; ++__i) 356 __is >> __densities[__i]; 357 __is >> __n; 358 vector<result_type> __areas(__n); 359 for (size_t __i = 0; __i < __n; ++__i) 360 __is >> __areas[__i]; 361 if (!__is.fail()) 362 { 363 swap(__x.__p_.__b_, __b); 364 swap(__x.__p_.__densities_, __densities); 365 swap(__x.__p_.__areas_, __areas); 366 } 367 return __is; 368 } 369 370 _LIBCPP_END_NAMESPACE_STD 371 372 _LIBCPP_POP_MACROS 373 374 #endif // _LIBCPP___RANDOM_PIECEWISE_LINEAR_DISTRIBUTION_H 375