1 /* boost random/piecewise_linear_distribution.hpp header file
3 * Copyright Steven Watanabe 2011
4 * Distributed under the Boost Software License, Version 1.0. (See
5 * accompanying file LICENSE_1_0.txt or copy at
6 * http://www.boost.org/LICENSE_1_0.txt)
8 * See http://www.boost.org for most recent version including documentation.
13 #ifndef BOOST_RANDOM_PIECEWISE_LINEAR_DISTRIBUTION_HPP_INCLUDED
14 #define BOOST_RANDOM_PIECEWISE_LINEAR_DISTRIBUTION_HPP_INCLUDED
20 #include <boost/assert.hpp>
21 #include <boost/random/uniform_real.hpp>
22 #include <boost/random/discrete_distribution.hpp>
23 #include <boost/random/detail/config.hpp>
24 #include <boost/random/detail/operators.hpp>
25 #include <boost/random/detail/vector_io.hpp>
27 #ifndef BOOST_NO_CXX11_HDR_INITIALIZER_LIST
28 #include <initializer_list>
31 #include <boost/range/begin.hpp>
32 #include <boost/range/end.hpp>
38 * The class @c piecewise_linear_distribution models a \random_distribution.
40 template<class RealType = double>
41 class piecewise_linear_distribution {
43 typedef std::size_t input_type;
44 typedef RealType result_type;
49 typedef piecewise_linear_distribution distribution_type;
52 * Constructs a @c param_type object, representing a distribution
53 * that produces values uniformly distributed in the range [0, 1).
57 _weights.push_back(RealType(1));
58 _weights.push_back(RealType(1));
59 _intervals.push_back(RealType(0));
60 _intervals.push_back(RealType(1));
63 * Constructs a @c param_type object from two iterator ranges
64 * containing the interval boundaries and weights at the boundaries.
65 * If there are fewer than two boundaries, then this is equivalent to
66 * the default constructor and the distribution will produce values
67 * uniformly distributed in the range [0, 1).
69 * The values of the interval boundaries must be strictly
70 * increasing, and the number of weights must be the same as
71 * the number of interval boundaries. If there are extra
72 * weights, they are ignored.
74 template<class IntervalIter, class WeightIter>
75 param_type(IntervalIter intervals_first, IntervalIter intervals_last,
76 WeightIter weight_first)
77 : _intervals(intervals_first, intervals_last)
79 if(_intervals.size() < 2) {
81 _weights.push_back(RealType(1));
82 _weights.push_back(RealType(1));
83 _intervals.push_back(RealType(0));
84 _intervals.push_back(RealType(1));
86 _weights.reserve(_intervals.size());
87 for(std::size_t i = 0; i < _intervals.size(); ++i) {
88 _weights.push_back(*weight_first++);
92 #ifndef BOOST_NO_CXX11_HDR_INITIALIZER_LIST
94 * Constructs a @c param_type object from an initializer_list
95 * containing the interval boundaries and a unary function
96 * specifying the weights at the boundaries. Each weight is
97 * determined by calling the function at the corresponding point.
99 * If the initializer_list contains fewer than two elements,
100 * this is equivalent to the default constructor and the
101 * distribution will produce values uniformly distributed
102 * in the range [0, 1).
104 template<class T, class F>
105 param_type(const std::initializer_list<T>& il, F f)
106 : _intervals(il.begin(), il.end())
108 if(_intervals.size() < 2) {
110 _weights.push_back(RealType(1));
111 _weights.push_back(RealType(1));
112 _intervals.push_back(RealType(0));
113 _intervals.push_back(RealType(1));
115 _weights.reserve(_intervals.size());
116 for(typename std::vector<RealType>::const_iterator
117 iter = _intervals.begin(), end = _intervals.end();
120 _weights.push_back(f(*iter));
126 * Constructs a @c param_type object from Boost.Range ranges holding
127 * the interval boundaries and the weights at the boundaries. If
128 * there are fewer than two interval boundaries, this is equivalent
129 * to the default constructor and the distribution will produce
130 * values uniformly distributed in the range [0, 1). The
131 * number of weights must be equal to the number of
132 * interval boundaries.
134 template<class IntervalRange, class WeightRange>
135 param_type(const IntervalRange& intervals_arg,
136 const WeightRange& weights_arg)
137 : _intervals(boost::begin(intervals_arg), boost::end(intervals_arg)),
138 _weights(boost::begin(weights_arg), boost::end(weights_arg))
140 if(_intervals.size() < 2) {
142 _weights.push_back(RealType(1));
143 _weights.push_back(RealType(1));
145 _intervals.push_back(RealType(0));
146 _intervals.push_back(RealType(1));
151 * Constructs the parameters for a distribution that approximates a
152 * function. The range of the distribution is [xmin, xmax). This
153 * range is divided into nw equally sized intervals and the weights
154 * are found by calling the unary function f on the boundaries of the
158 param_type(std::size_t nw, RealType xmin, RealType xmax, F f)
160 std::size_t n = (nw == 0) ? 1 : nw;
161 double delta = (xmax - xmin) / n;
162 BOOST_ASSERT(delta > 0);
163 for(std::size_t k = 0; k < n; ++k) {
164 _weights.push_back(f(xmin + k*delta));
165 _intervals.push_back(xmin + k*delta);
167 _weights.push_back(f(xmax));
168 _intervals.push_back(xmax);
171 /** Returns a vector containing the interval boundaries. */
172 std::vector<RealType> intervals() const { return _intervals; }
175 * Returns a vector containing the probability densities
176 * at all the interval boundaries.
178 std::vector<RealType> densities() const
180 RealType sum = static_cast<RealType>(0);
181 for(std::size_t i = 0; i < _intervals.size() - 1; ++i) {
182 RealType width = _intervals[i + 1] - _intervals[i];
183 sum += (_weights[i] + _weights[i + 1]) * width / 2;
185 std::vector<RealType> result;
186 result.reserve(_weights.size());
187 for(typename std::vector<RealType>::const_iterator
188 iter = _weights.begin(), end = _weights.end();
191 result.push_back(*iter / sum);
196 /** Writes the parameters to a @c std::ostream. */
197 BOOST_RANDOM_DETAIL_OSTREAM_OPERATOR(os, param_type, parm)
199 detail::print_vector(os, parm._intervals);
200 detail::print_vector(os, parm._weights);
204 /** Reads the parameters from a @c std::istream. */
205 BOOST_RANDOM_DETAIL_ISTREAM_OPERATOR(is, param_type, parm)
207 std::vector<RealType> new_intervals;
208 std::vector<RealType> new_weights;
209 detail::read_vector(is, new_intervals);
210 detail::read_vector(is, new_weights);
212 parm._intervals.swap(new_intervals);
213 parm._weights.swap(new_weights);
218 /** Returns true if the two sets of parameters are the same. */
219 BOOST_RANDOM_DETAIL_EQUALITY_OPERATOR(param_type, lhs, rhs)
221 return lhs._intervals == rhs._intervals
222 && lhs._weights == rhs._weights;
224 /** Returns true if the two sets of parameters are different. */
225 BOOST_RANDOM_DETAIL_INEQUALITY_OPERATOR(param_type)
228 friend class piecewise_linear_distribution;
230 std::vector<RealType> _intervals;
231 std::vector<RealType> _weights;
235 * Creates a new @c piecewise_linear_distribution that
236 * produces values uniformly distributed in the range [0, 1).
238 piecewise_linear_distribution()
243 * Constructs a piecewise_linear_distribution from two iterator ranges
244 * containing the interval boundaries and the weights at the boundaries.
245 * If there are fewer than two boundaries, then this is equivalent to
246 * the default constructor and creates a distribution that
247 * produces values uniformly distributed in the range [0, 1).
249 * The values of the interval boundaries must be strictly
250 * increasing, and the number of weights must be equal to
251 * the number of interval boundaries. If there are extra
252 * weights, they are ignored.
257 * double intervals[] = { 0.0, 1.0, 2.0 };
258 * double weights[] = { 0.0, 1.0, 0.0 };
259 * piecewise_constant_distribution<> dist(
260 * &intervals[0], &intervals[0] + 3, &weights[0]);
263 * produces a triangle distribution.
265 template<class IntervalIter, class WeightIter>
266 piecewise_linear_distribution(IntervalIter first_interval,
267 IntervalIter last_interval,
268 WeightIter first_weight)
269 : _intervals(first_interval, last_interval)
271 if(_intervals.size() < 2) {
274 _weights.reserve(_intervals.size());
275 for(std::size_t i = 0; i < _intervals.size(); ++i) {
276 _weights.push_back(*first_weight++);
281 #ifndef BOOST_NO_CXX11_HDR_INITIALIZER_LIST
283 * Constructs a piecewise_linear_distribution from an
284 * initializer_list containing the interval boundaries
285 * and a unary function specifying the weights. Each
286 * weight is determined by calling the function at the
287 * corresponding interval boundary.
289 * If the initializer_list contains fewer than two elements,
290 * this is equivalent to the default constructor and the
291 * distribution will produce values uniformly distributed
292 * in the range [0, 1).
294 template<class T, class F>
295 piecewise_linear_distribution(std::initializer_list<T> il, F f)
296 : _intervals(il.begin(), il.end())
298 if(_intervals.size() < 2) {
301 _weights.reserve(_intervals.size());
302 for(typename std::vector<RealType>::const_iterator
303 iter = _intervals.begin(), end = _intervals.end();
306 _weights.push_back(f(*iter));
313 * Constructs a piecewise_linear_distribution from Boost.Range
314 * ranges holding the interval boundaries and the weights. If
315 * there are fewer than two interval boundaries, this is equivalent
316 * to the default constructor and the distribution will produce
317 * values uniformly distributed in the range [0, 1). The
318 * number of weights must be equal to the number of
319 * interval boundaries.
321 template<class IntervalsRange, class WeightsRange>
322 piecewise_linear_distribution(const IntervalsRange& intervals_arg,
323 const WeightsRange& weights_arg)
324 : _intervals(boost::begin(intervals_arg), boost::end(intervals_arg)),
325 _weights(boost::begin(weights_arg), boost::end(weights_arg))
327 if(_intervals.size() < 2) {
334 * Constructs a piecewise_linear_distribution that approximates a
335 * function. The range of the distribution is [xmin, xmax). This
336 * range is divided into nw equally sized intervals and the weights
337 * are found by calling the unary function f on the interval boundaries.
340 piecewise_linear_distribution(std::size_t nw,
345 if(nw == 0) { nw = 1; }
346 RealType delta = (xmax - xmin) / nw;
347 _intervals.reserve(nw + 1);
348 for(std::size_t i = 0; i < nw; ++i) {
349 RealType x = xmin + i * delta;
350 _intervals.push_back(x);
351 _weights.push_back(f(x));
353 _intervals.push_back(xmax);
354 _weights.push_back(f(xmax));
358 * Constructs a piecewise_linear_distribution from its parameters.
360 explicit piecewise_linear_distribution(const param_type& parm)
361 : _intervals(parm._intervals),
362 _weights(parm._weights)
368 * Returns a value distributed according to the parameters of the
369 * piecewise_linear_distribution.
372 RealType operator()(URNG& urng) const
374 std::size_t i = _bins(urng);
375 bool is_in_rectangle = (i % 2 == 0);
377 uniform_real<RealType> dist(_intervals[i], _intervals[i+1]);
378 if(is_in_rectangle) {
380 } else if(_weights[i] < _weights[i+1]) {
381 return (std::max)(dist(urng), dist(urng));
383 return (std::min)(dist(urng), dist(urng));
388 * Returns a value distributed according to the parameters
389 * specified by param.
392 RealType operator()(URNG& urng, const param_type& parm) const
394 return piecewise_linear_distribution(parm)(urng);
397 /** Returns the smallest value that the distribution can produce. */
398 result_type min BOOST_PREVENT_MACRO_SUBSTITUTION () const
399 { return _intervals.front(); }
400 /** Returns the largest value that the distribution can produce. */
401 result_type max BOOST_PREVENT_MACRO_SUBSTITUTION () const
402 { return _intervals.back(); }
405 * Returns a vector containing the probability densities
406 * at the interval boundaries.
408 std::vector<RealType> densities() const
410 RealType sum = static_cast<RealType>(0);
411 for(std::size_t i = 0; i < _intervals.size() - 1; ++i) {
412 RealType width = _intervals[i + 1] - _intervals[i];
413 sum += (_weights[i] + _weights[i + 1]) * width / 2;
415 std::vector<RealType> result;
416 result.reserve(_weights.size());
417 for(typename std::vector<RealType>::const_iterator
418 iter = _weights.begin(), end = _weights.end();
421 result.push_back(*iter / sum);
425 /** Returns a vector containing the interval boundaries. */
426 std::vector<RealType> intervals() const { return _intervals; }
428 /** Returns the parameters of the distribution. */
429 param_type param() const
431 return param_type(_intervals, _weights);
433 /** Sets the parameters of the distribution. */
434 void param(const param_type& parm)
436 std::vector<RealType> new_intervals(parm._intervals);
437 std::vector<RealType> new_weights(parm._weights);
438 init(new_intervals, new_weights);
439 _intervals.swap(new_intervals);
440 _weights.swap(new_weights);
444 * Effects: Subsequent uses of the distribution do not depend
445 * on values produced by any engine prior to invoking reset.
447 void reset() { _bins.reset(); }
449 /** Writes a distribution to a @c std::ostream. */
450 BOOST_RANDOM_DETAIL_OSTREAM_OPERATOR(
451 os, piecewise_linear_distribution, pld)
457 /** Reads a distribution from a @c std::istream */
458 BOOST_RANDOM_DETAIL_ISTREAM_OPERATOR(
459 is, piecewise_linear_distribution, pld)
469 * Returns true if the two distributions will return the
470 * same sequence of values, when passed equal generators.
472 BOOST_RANDOM_DETAIL_EQUALITY_OPERATOR(
473 piecewise_linear_distribution, lhs, rhs)
475 return lhs._intervals == rhs._intervals && lhs._weights == rhs._weights;
478 * Returns true if the two distributions may return different
479 * sequences of values, when passed equal generators.
481 BOOST_RANDOM_DETAIL_INEQUALITY_OPERATOR(piecewise_linear_distribution)
485 /// @cond \show_private
487 void init(const std::vector<RealType>& intervals_arg,
488 const std::vector<RealType>& weights_arg)
491 std::vector<RealType> bin_weights;
492 bin_weights.reserve((intervals_arg.size() - 1) * 2);
493 for(std::size_t i = 0; i < intervals_arg.size() - 1; ++i) {
494 RealType width = intervals_arg[i + 1] - intervals_arg[i];
495 RealType w1 = weights_arg[i];
496 RealType w2 = weights_arg[i + 1];
497 bin_weights.push_back((std::min)(w1, w2) * width);
498 bin_weights.push_back(abs(w1 - w2) * width / 2);
500 typedef discrete_distribution<std::size_t, RealType> bins_type;
501 typename bins_type::param_type bins_param(bin_weights);
502 _bins.param(bins_param);
507 init(_intervals, _weights);
513 _intervals.push_back(RealType(0));
514 _intervals.push_back(RealType(1));
516 _weights.push_back(RealType(1));
517 _weights.push_back(RealType(1));
521 discrete_distribution<std::size_t, RealType> _bins;
522 std::vector<RealType> _intervals;
523 std::vector<RealType> _weights;