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