[Boost-commit] svn:boost r68028 - in trunk: boost/random libs/random/test

Subject: [Boost-commit] svn:boost r68028 - in trunk: boost/random libs/random/test
From: steven_at_[hidden]
Date: 2011-01-11 23:53:09

Author: steven_watanabe
Date: 2011-01-11 23:53:08 EST (Tue, 11 Jan 2011)
New Revision: 68028
URL: http://svn.boost.org/trac/boost/changeset/68028

Log:
Implement piecewise_linear_distribution.
Added:
   trunk/boost/random/piecewise_linear_distribution.hpp   (contents, props changed)
   trunk/libs/random/test/test_piecewise_linear.cpp   (contents, props changed)
   trunk/libs/random/test/test_piecewise_linear_distribution.cpp   (contents, props changed)
Text files modified: 
   trunk/libs/random/test/Jamfile.v2 |     2 ++                                      
   1 files changed, 2 insertions(+), 0 deletions(-)

Added: trunk/boost/random/piecewise_linear_distribution.hpp
==============================================================================
--- (empty file)
+++ trunk/boost/random/piecewise_linear_distribution.hpp	2011-01-11 23:53:08 EST (Tue, 11 Jan 2011)
@@ -0,0 +1,535 @@
+/* boost random/piecewise_linear_distribution.hpp header file
+ *
+ * Copyright Steven Watanabe 2011
+ * Distributed under the Boost Software License, Version 1.0. (See
+ * accompanying file LICENSE_1_0.txt or copy at
+ * http://www.boost.org/LICENSE_1_0.txt)
+ *
+ * See http://www.boost.org for most recent version including documentation.
+ *
+ * $Id$
+ */
+
+#ifndef BOOST_RANDOM_PIECEWISE_LINEAR_DISTRIBUTION_HPP_INCLUDED
+#define BOOST_RANDOM_PIECEWISE_LINEAR_DISTRIBUTION_HPP_INCLUDED
+
+#include <vector>
+#include <algorithm>
+#include <cassert>
+#include <cmath>
+#include <cstdlib>
+#include <boost/random/uniform_real.hpp>
+#include <boost/random/discrete_distribution.hpp>
+#include <boost/random/detail/config.hpp>
+#include <boost/random/detail/operators.hpp>
+#include <boost/random/detail/vector_io.hpp>
+
+#ifndef BOOST_NO_INITIALIZER_LISTS
+#include <initializer_list>
+#endif
+
+#include <boost/range/begin.hpp>
+#include <boost/range/end.hpp>
+
+namespace boost {
+namespace random {
+
+/**
+ * The class @c piecewise_linear_distribution models a \random_distribution.
+ */
+template<class RealType = double, class WeightType = double>
+class piecewise_linear_distribution {
+public:
+    typedef RealType input_type;
+    typedef RealType result_type;
+
+    class param_type {
+    public:
+
+        typedef piecewise_linear_distribution distribution_type;
+
+        /**
+         * Constructs a @c param_type object, representing a distribution
+         * that produces values uniformly distributed in the range [0, 1).
+         */
+        param_type()
+        {
+            _weights.push_back(WeightType(1));
+            _weights.push_back(WeightType(1));
+            _intervals.push_back(RealType(0));
+            _intervals.push_back(RealType(1));
+        }
+        /**
+         * Constructs a @c param_type object from two iterator ranges
+         * containing the interval boundaries and weights at the boundaries.
+         * If there are fewer than two boundaries, then this is equivalent to
+         * the default constructor and the distribution will produce values
+         * uniformly distributed in the range [0, 1).
+         *
+         * The values of the interval boundaries must be strictly
+         * increasing, and the number of weights must be the same as
+         * the number of interval boundaries.  If there are extra
+         * weights, they are ignored.
+         */
+        template<class IntervalIter, class WeightIter>
+        param_type(IntervalIter intervals_first, IntervalIter intervals_last,
+                   WeightIter weight_first)
+          : _intervals(intervals_first, intervals_last)
+        {
+            if(_intervals.size() < 2) {
+                _intervals.clear();
+                _weights.push_back(WeightType(1));
+                _weights.push_back(WeightType(1));
+                _intervals.push_back(RealType(0));
+                _intervals.push_back(RealType(1));
+            } else {
+                _weights.reserve(_intervals.size());
+                for(std::size_t i = 0; i < _intervals.size(); ++i) {
+                    _weights.push_back(*weight_first++);
+                }
+            }
+        }
+#ifndef BOOST_NO_INITIALIZER_LISTS
+        /**
+         * Constructs a @c param_type object from an initializer_list
+         * containing the interval boundaries and a unary function
+         * specifying the weights at the boundaries.  Each weight is
+         * determined by calling the function at the corresponding point.
+         *
+         * If the initializer_list contains fewer than two elements,
+         * this is equivalent to the default constructor and the
+         * distribution will produce values uniformly distributed
+         * in the range [0, 1).
+         */
+        template<class T, class F>
+        param_type(const std::initializer_list<T>& il, F f)
+          : _intervals(il.begin(), il.end())
+        {
+            if(_intervals.size() < 2) {
+                _intervals.clear();
+                _weights.push_back(WeightType(1));
+                _weights.push_back(WeightType(1));
+                _intervals.push_back(RealType(0));
+                _intervals.push_back(RealType(1));
+            } else {
+                _weights.reserve(_intervals.size());
+                for(typename std::vector<RealType>::const_iterator
+                    iter = _intervals.begin(), end = _intervals.end();
+                    iter != end; ++iter)
+                {
+                    _weights.push_back(f(*iter));
+                }
+            }
+        }
+#endif
+        /**
+         * Constructs a @c param_type object from Boost.Range ranges holding
+         * the interval boundaries and the weights at the boundaries.  If
+         * there are fewer than two interval boundaries, this is equivalent
+         * to the default constructor and the distribution will produce
+         * values uniformly distributed in the range [0, 1).  The
+         * number of weights must be equal to the number of
+         * interval boundaries.
+         */
+        template<class IntervalRange, class WeightRange>
+        param_type(const IntervalRange& intervals_arg,
+                   const WeightRange& weights_arg)
+          : _intervals(boost::begin(intervals_arg), boost::end(intervals_arg)),
+            _weights(boost::begin(weights_arg), boost::end(weights_arg))
+        {
+            if(_intervals.size() < 2) {
+                _weights.clear();
+                _weights.push_back(WeightType(1));
+                _weights.push_back(WeightType(1));
+                _intervals.clear();
+                _intervals.push_back(RealType(0));
+                _intervals.push_back(RealType(1));
+            }
+        }
+
+        /**
+         * Constructs the parameters for a distribution that approximates a
+         * function.  The range of the distribution is [xmin, xmax).  This
+         * range is divided into nw equally sized intervals and the weights
+         * are found by calling the unary function f on the boundaries of the
+         * intervals.
+         */
+        template<class F>
+        param_type(std::size_t nw, RealType xmin, RealType xmax, F f)
+        {
+            std::size_t n = (nw == 0) ? 1 : nw;
+            double delta = (xmax - xmin) / n;
+            assert(delta > 0);
+            for(std::size_t k = 0; k < n; ++k) {
+                _weights.push_back(f(xmin + k*delta));
+                _intervals.push_back(xmin + k*delta);
+            }
+            _weights.push_back(f(xmax));
+            _intervals.push_back(xmax);
+        }
+
+        /**  Returns a vector containing the interval boundaries. */
+        std::vector<RealType> intervals() const { return _intervals; }
+
+        /**
+         * Returns a vector containing the probability densities
+         * at all the interval boundaries.
+         */
+        std::vector<WeightType> densities() const
+        {
+            WeightType sum = static_cast<WeightType>(0);
+            for(std::size_t i = 0; i < _intervals.size() - 1; ++i) {
+                RealType width = _intervals[i + 1] - _intervals[i];
+                sum += (_weights[i] + _weights[i + 1]) * width / 2;
+            }
+            std::vector<WeightType> result;
+            result.reserve(_weights.size());
+            for(typename std::vector<WeightType>::const_iterator
+                iter = _weights.begin(), end = _weights.end();
+                iter != end; ++iter)
+            {
+                result.push_back(*iter / sum);
+            }
+            return result;
+        }
+
+        /** Writes the parameters to a @c std::ostream. */
+        BOOST_RANDOM_DETAIL_OSTREAM_OPERATOR(os, param_type, parm)
+        {
+            detail::print_vector(os, parm._intervals);
+            detail::print_vector(os, parm._weights);
+            return os;
+        }
+        
+        /** Reads the parameters from a @c std::istream. */
+        BOOST_RANDOM_DETAIL_ISTREAM_OPERATOR(is, param_type, parm)
+        {
+            std::vector<RealType> new_intervals;
+            std::vector<WeightType> new_weights;
+            detail::read_vector(is, new_intervals);
+            detail::read_vector(is, new_weights);
+            if(is) {
+                parm._intervals.swap(new_intervals);
+                parm._weights.swap(new_weights);
+            }
+            return is;
+        }
+
+        /** Returns true if the two sets of parameters are the same. */
+        BOOST_RANDOM_DETAIL_EQUALITY_OPERATOR(param_type, lhs, rhs)
+        {
+            return lhs._intervals == rhs._intervals
+                && lhs._weights == rhs._weights;
+        }
+        /** Returns true if the two sets of parameters are different. */
+        BOOST_RANDOM_DETAIL_INEQUALITY_OPERATOR(param_type)
+
+    private:
+
+        /// @cond
+
+        friend class piecewise_linear_distribution;
+
+        std::vector<RealType> _intervals;
+        std::vector<WeightType> _weights;
+
+        /// @endcond
+    };
+
+    /**
+     * Creates a new @c piecewise_linear_distribution that
+     * produces values uniformly distributed in the range [0, 1).
+     */
+    piecewise_linear_distribution()
+    {
+        default_init();
+    }
+    /**
+     * Constructs a piecewise_linear_distribution from two iterator ranges
+     * containing the interval boundaries and the weights at the boundaries.
+     * If there are fewer than two boundaries, then this is equivalent to
+     * the default constructor and creates a distribution that
+     * produces values uniformly distributed in the range [0, 1).
+     *
+     * The values of the interval boundaries must be strictly
+     * increasing, and the number of weights must be equal to
+     * the number of interval boundaries.  If there are extra
+     * weights, they are ignored.
+     *
+     * For example,
+     *
+     * @code
+     * double intervals[] = { 0.0, 1.0, 2.0 };
+     * double weights[] = { 0.0, 1.0, 0.0 };
+     * piecewise_constant_distribution<> dist(
+     *     &intervals[0], &intervals[0] + 3, &weights[0]);
+     * @endcode
+     *
+     * produces a triangle distribution.
+     */
+    template<class IntervalIter, class WeightIter>
+    piecewise_linear_distribution(IntervalIter first_interval,
+                                  IntervalIter last_interval,
+                                  WeightIter first_weight)
+      : _intervals(first_interval, last_interval)
+    {
+        if(_intervals.size() < 2) {
+            default_init();
+        } else {
+            _weights.reserve(_intervals.size());
+            for(std::size_t i = 0; i < _intervals.size(); ++i) {
+                _weights.push_back(*first_weight++);
+            }
+            init();
+        }
+    }
+#ifndef BOOST_NO_INITIALIZER_LISTS
+    /**
+     * Constructs a piecewise_linear_distribution from an
+     * initializer_list containing the interval boundaries
+     * and a unary function specifying the weights.  Each
+     * weight is determined by calling the function at the
+     * corresponding interval boundary.
+     *
+     * If the initializer_list contains fewer than two elements,
+     * this is equivalent to the default constructor and the
+     * distribution will produce values uniformly distributed
+     * in the range [0, 1).
+     */
+    template<class T, class F>
+    piecewise_linear_distribution(std::initializer_list<T> il, F f)
+      : _intervals(il.begin(), il.end())
+    {
+        if(_intervals.size() < 2) {
+            default_init();
+        } else {
+            _weights.reserve(_intervals.size());
+            for(typename std::vector<RealType>::const_iterator
+                iter = _intervals.begin(), end = _intervals.end();
+                iter != end; ++iter)
+            {
+                _weights.push_back(f(*iter));
+            }
+            init();
+        }
+    }
+#endif
+    /**
+     * Constructs a piecewise_linear_distribution from Boost.Range
+     * ranges holding the interval boundaries and the weights.  If
+     * there are fewer than two interval boundaries, this is equivalent
+     * to the default constructor and the distribution will produce
+     * values uniformly distributed in the range [0, 1).  The
+     * number of weights must be equal to the number of
+     * interval boundaries.
+     */
+    template<class IntervalsRange, class WeightsRange>
+    piecewise_linear_distribution(const IntervalsRange& intervals_arg,
+                                  const WeightsRange& weights_arg)
+      : _intervals(boost::begin(intervals_arg), boost::end(intervals_arg)),
+        _weights(boost::begin(weights_arg), boost::end(weights_arg))
+    {
+        if(_intervals.size() < 2) {
+            default_init();
+        } else {
+            init();
+        }
+    }
+    /**
+     * Constructs a piecewise_linear_distribution that approximates a
+     * function.  The range of the distribution is [xmin, xmax).  This
+     * range is divided into nw equally sized intervals and the weights
+     * are found by calling the unary function f on the interval boundaries.
+     */
+    template<class F>
+    piecewise_linear_distribution(std::size_t nw,
+                                  RealType xmin,
+                                  RealType xmax,
+                                  F f)
+    {
+        if(nw == 0) { nw = 1; }
+        RealType delta = (xmax - xmin) / nw;
+        _intervals.reserve(nw + 1);
+        for(std::size_t i = 0; i < nw; ++i) {
+            RealType x = xmin + i * delta;
+            _intervals.push_back(x);
+            _weights.push_back(f(x));
+        }
+        _intervals.push_back(xmax);
+        _weights.push_back(f(xmax));
+        init();
+    }
+    /**
+     * Constructs a piecewise_linear_distribution from its parameters.
+     */
+    explicit piecewise_linear_distribution(const param_type& parm)
+      : _intervals(parm._intervals),
+        _weights(parm._weights)
+    {
+        init();
+    }
+
+    /**
+     * Returns a value distributed according to the parameters of the
+     * piecewise_linear_distribution.
+     */
+    template<class URNG>
+    RealType operator()(URNG& urng) const
+    {
+        std::size_t i = _bins(urng);
+        bool is_in_rectangle = (i % 2 == 0);
+        i = i / 2;
+        uniform_real<RealType> dist(_intervals[i], _intervals[i+1]);
+        if(is_in_rectangle) {
+            return dist(urng);
+        } else if(_weights[i] < _weights[i+1]) {
+            return std::max(dist(urng), dist(urng));
+        } else {
+            return std::min(dist(urng), dist(urng));
+        }
+    }
+    
+    /**
+     * Returns a value distributed according to the parameters
+     * specified by param.
+     */
+    template<class URNG>
+    RealType operator()(URNG& urng, const param_type& parm) const
+    {
+        return piecewise_linear_distribution(parm)(urng);
+    }
+    
+    /** Returns the smallest value that the distribution can produce. */
+    result_type min BOOST_PREVENT_MACRO_SUBSTITUTION () const
+    { return _intervals.front(); }
+    /** Returns the largest value that the distribution can produce. */
+    result_type max BOOST_PREVENT_MACRO_SUBSTITUTION () const
+    { return _intervals.back(); }
+
+    /**
+     * Returns a vector containing the probability densities
+     * at the interval boundaries.
+     */
+    std::vector<WeightType> densities() const
+    {
+        WeightType sum = static_cast<WeightType>(0);
+        for(std::size_t i = 0; i < _intervals.size() - 1; ++i) {
+            RealType width = _intervals[i + 1] - _intervals[i];
+            sum += (_weights[i] + _weights[i + 1]) * width / 2;
+        }
+        std::vector<WeightType> result;
+        result.reserve(_weights.size());
+        for(typename std::vector<WeightType>::const_iterator
+            iter = _weights.begin(), end = _weights.end();
+            iter != end; ++iter)
+        {
+            result.push_back(*iter / sum);
+        }
+        return result;
+    }
+    /**  Returns a vector containing the interval boundaries. */
+    std::vector<RealType> intervals() const { return _intervals; }
+
+    /** Returns the parameters of the distribution. */
+    param_type param() const
+    {
+        return param_type(_intervals, _weights);
+    }
+    /** Sets the parameters of the distribution. */
+    void param(const param_type& parm)
+    {
+        std::vector<RealType> new_intervals(parm._intervals);
+        std::vector<WeightType> new_weights(parm._weights);
+        init(new_intervals, new_weights);
+        _intervals.swap(new_intervals);
+        _weights.swap(new_weights);
+    }
+    
+    /**
+     * Effects: Subsequent uses of the distribution do not depend
+     * on values produced by any engine prior to invoking reset.
+     */
+    void reset() { _bins.reset(); }
+
+    /** Writes a distribution to a @c std::ostream. */
+    BOOST_RANDOM_DETAIL_OSTREAM_OPERATOR(
+        os, piecewise_linear_distribution, pld)
+    {
+        os << pld.param();
+        return os;
+    }
+
+    /** Reads a distribution from a @c std::istream */
+    BOOST_RANDOM_DETAIL_ISTREAM_OPERATOR(
+        is, piecewise_linear_distribution, pld)
+    {
+        param_type parm;
+        if(is >> parm) {
+            pld.param(parm);
+        }
+        return is;
+    }
+
+    /**
+     * Returns true if the two distributions will return the
+     * same sequence of values, when passed equal generators.
+     */
+    BOOST_RANDOM_DETAIL_EQUALITY_OPERATOR(
+        piecewise_linear_distribution, lhs,  rhs)
+    {
+        return lhs._intervals == rhs._intervals && lhs._weights == rhs._weights;
+    }
+    /**
+     * Returns true if the two distributions may return different
+     * sequences of values, when passed equal generators.
+     */
+    BOOST_RANDOM_DETAIL_INEQUALITY_OPERATOR(piecewise_linear_distribution)
+
+private:
+
+    /// @cond
+
+    void init(const std::vector<RealType>& intervals_arg,
+              const std::vector<WeightType>& weights_arg)
+    {
+        std::vector<WeightType> bin_weights;
+        bin_weights.reserve((intervals_arg.size() - 1) * 2);
+        for(std::size_t i = 0; i < intervals_arg.size() - 1; ++i) {
+            RealType width = intervals_arg[i + 1] - intervals_arg[i];
+            WeightType w1 = weights_arg[i];
+            WeightType w2 = weights_arg[i + 1];
+            bin_weights.push_back(std::min(w1, w2) * width);
+            bin_weights.push_back(std::abs(w1 - w2) * width / 2);
+        }
+        typedef discrete_distribution<std::size_t, WeightType> bins_type;
+        typename bins_type::param_type bins_param(bin_weights);
+        _bins.param(bins_param);
+    }
+
+    void init()
+    {
+        init(_intervals, _weights);
+    }
+
+    void default_init()
+    {
+        _intervals.clear();
+        _intervals.push_back(RealType(0));
+        _intervals.push_back(RealType(1));
+        _weights.clear();
+        _weights.push_back(WeightType(1));
+        _weights.push_back(WeightType(1));
+        init();
+    }
+
+    discrete_distribution<std::size_t, WeightType> _bins;
+    std::vector<RealType> _intervals;
+    std::vector<WeightType> _weights;
+
+    /// @endcond
+};
+
+}
+}
+
+#endif

Modified: trunk/libs/random/test/Jamfile.v2
==============================================================================
--- trunk/libs/random/test/Jamfile.v2	(original)
+++ trunk/libs/random/test/Jamfile.v2	2011-01-11 23:53:08 EST (Tue, 11 Jan 2011)
@@ -71,6 +71,8 @@
 run test_normal_distribution.cpp /boost//unit_test_framework ;
 run test_piecewise_constant.cpp ;
 run test_piecewise_constant_distribution.cpp /boost//unit_test_framework ;
+run test_piecewise_linear.cpp ;
+run test_piecewise_linear_distribution.cpp /boost//unit_test_framework ;
 
 # run nondet_random_speed.cpp ;
 # run random_device.cpp ;

Added: trunk/libs/random/test/test_piecewise_linear.cpp
==============================================================================
--- (empty file)
+++ trunk/libs/random/test/test_piecewise_linear.cpp	2011-01-11 23:53:08 EST (Tue, 11 Jan 2011)
@@ -0,0 +1,131 @@
+/* test_piecewise_linear.cpp
+ *
+ * Copyright Steven Watanabe 2011
+ * Distributed under the Boost Software License, Version 1.0. (See
+ * accompanying file LICENSE_1_0.txt or copy at
+ * http://www.boost.org/LICENSE_1_0.txt)
+ *
+ * $Id$
+ *
+ */
+
+#include <boost/random/piecewise_linear_distribution.hpp>
+#include <boost/random/uniform_int.hpp>
+#include <boost/random/mersenne_twister.hpp>
+#include <boost/lexical_cast.hpp>
+#include <boost/exception/diagnostic_information.hpp>
+#include <vector>
+#include <iostream>
+#include <numeric>
+
+#include "chi_squared_test.hpp"
+
+bool do_test(int n, long long max) {
+    std::cout << "running piecewise_linear(p0, p1, ..., p" << n-1 << ")" << " " << max << " times: " << std::flush;
+
+    std::vector<double> weights;
+    std::vector<double> expected;
+    {
+        boost::mt19937 egen;
+        for(int i = 0; i < n; ++i) {
+            weights.push_back(egen());
+        }
+        for(int i = 0; i < n - 1; ++i) {
+            expected.push_back((weights[i] + weights[i + 1]) / 2);
+        }
+        double sum = std::accumulate(expected.begin(), expected.end(), 0.0);
+        for(std::vector<double>::iterator iter = expected.begin(), end = expected.end(); iter != end; ++iter) {
+            *iter /= sum;
+        }
+    }
+    std::vector<double> intervals;
+    for(int i = 0; i < n; ++i) {
+        intervals.push_back(i);
+    }
+    
+    boost::random::piecewise_linear_distribution<> dist(intervals, weights);
+    boost::mt19937 gen;
+    std::vector<long long> results(expected.size());
+    for(long long i = 0; i < max; ++i) {
+        ++results[static_cast<std::size_t>(dist(gen))];
+    }
+
+    long long sum = std::accumulate(results.begin(), results.end(), 0ll);
+    if(sum != max) {
+        std::cout << "*** Failed: incorrect total: " << sum << " ***" << std::endl;
+        return false;
+    }
+    double chsqr = chi_squared_test(results, expected, max);
+
+    bool result = chsqr < 0.99;
+    const char* err = result? "" : "*";
+    std::cout << std::setprecision(17) << chsqr << err << std::endl;
+
+    std::cout << std::setprecision(6);
+
+    return result;
+}
+
+bool do_tests(int repeat, int max_n, long long trials) {
+    boost::mt19937 gen;
+    boost::uniform_int<> idist(1, max_n);
+    int errors = 0;
+    for(int i = 0; i < repeat; ++i) {
+        if(!do_test(idist(gen), trials)) {
+            ++errors;
+        }
+    }
+    if(errors != 0) {
+        std::cout << "*** " << errors << " errors detected ***" << std::endl;
+    }
+    return errors == 0;
+}
+
+int usage() {
+    std::cerr << "Usage: test_piecewise_linear -r <repeat> -n <max n> -t <trials>" << std::endl;
+    return 2;
+}
+
+template<class T>
+bool handle_option(int& argc, char**& argv, char opt, T& value) {
+    if(argv[0][1] == opt && argc > 1) {
+        --argc;
+        ++argv;
+        value = boost::lexical_cast<T>(argv[0]);
+        return true;
+    } else {
+        return false;
+    }
+}
+
+int main(int argc, char** argv) {
+    int repeat = 10;
+    int max_n = 100000;
+    long long trials = 1000000ll;
+
+    if(argc > 0) {
+        --argc;
+        ++argv;
+    }
+    while(argc > 0) {
+        if(argv[0][0] != '-') return usage();
+        else if(!handle_option(argc, argv, 'r', repeat)
+             && !handle_option(argc, argv, 'n', max_n)
+             && !handle_option(argc, argv, 't', trials)) {
+            return usage();
+        }
+        --argc;
+        ++argv;
+    }
+
+    try {
+        if(do_tests(repeat, max_n, trials)) {
+            return 0;
+        } else {
+            return EXIT_FAILURE;
+        }
+    } catch(...) {
+        std::cerr << boost::current_exception_diagnostic_information() << std::endl;
+        return EXIT_FAILURE;
+    }
+}

Added: trunk/libs/random/test/test_piecewise_linear_distribution.cpp
==============================================================================
--- (empty file)
+++ trunk/libs/random/test/test_piecewise_linear_distribution.cpp	2011-01-11 23:53:08 EST (Tue, 11 Jan 2011)
@@ -0,0 +1,247 @@
+/* test_piecewise_linear_distribution.cpp
+ *
+ * Copyright Steven Watanabe 2011
+ * Distributed under the Boost Software License, Version 1.0. (See
+ * accompanying file LICENSE_1_0.txt or copy at
+ * http://www.boost.org/LICENSE_1_0.txt)
+ *
+ * $Id$
+ *
+ */
+
+#include <boost/random/piecewise_linear_distribution.hpp>
+#include <boost/random/linear_congruential.hpp>
+#include <boost/assign/list_of.hpp>
+#include <sstream>
+#include <vector>
+
+#define BOOST_TEST_MAIN
+#include <boost/test/unit_test.hpp>
+
+struct gen {
+    double operator()(double arg) {
+        if(arg < 97) return 100;
+        else if(arg < 101) return 3;
+        else if(arg < 105) return 1;
+        else if(arg < 109) return 2;
+        else if(arg < 113) return 1;
+        else if(arg < 117) return 5;
+        else return 100;
+    }
+};
+
+#define CHECK_SEQUENCE(actual, expected)            \
+    do {                                            \
+        std::vector<double> _actual = (actual);     \
+        std::vector<double> _expected = (expected); \
+        BOOST_CHECK_EQUAL_COLLECTIONS(              \
+            _actual.begin(), _actual.end(),         \
+            _expected.begin(), _expected.end());    \
+    } while(false)
+
+using boost::assign::list_of;
+
+BOOST_AUTO_TEST_CASE(test_constructors) {
+    boost::random::piecewise_linear_distribution<> dist;
+    CHECK_SEQUENCE(dist.intervals(), list_of(0.0)(1.0));
+    CHECK_SEQUENCE(dist.densities(), list_of(1.0)(1.0));
+
+#ifndef BOOST_NO_INITIALIZER_LISTS
+    boost::random::piecewise_linear_distribution<> dist_il = {
+        { 99, 103, 107, 111, 115 },
+        gen()
+    };
+    CHECK_SEQUENCE(dist_il.intervals(), list_of(99)(103)(107)(111)(115));
+    CHECK_SEQUENCE(dist_il.densities(),
+        list_of(.09375)(.03125)(0.0625)(.03125)(.15625));
+
+    boost::random::piecewise_linear_distribution<> dist_il2 = {
+        { 99 },
+        gen()
+    };
+    CHECK_SEQUENCE(dist_il2.intervals(), list_of(0.0)(1.0));
+    CHECK_SEQUENCE(dist_il2.densities(), list_of(1.0)(1.0));
+#endif
+    std::vector<double> intervals = boost::assign::list_of(0)(1)(2)(3)(5);
+    std::vector<double> weights = boost::assign::list_of(3)(1)(2)(1)(2);
+    std::vector<double> intervals2 = boost::assign::list_of(99);
+    std::vector<double> weights2 = boost::assign::list_of(2);
+
+    boost::random::piecewise_linear_distribution<> dist_r(intervals, weights);
+    CHECK_SEQUENCE(dist_r.intervals(), list_of(0)(1)(2)(3)(5));
+    CHECK_SEQUENCE(dist_r.densities(), list_of(.375)(.125)(.25)(.125)(.25));
+
+    boost::random::piecewise_linear_distribution<>
+        dist_r2(intervals2, weights2);
+    CHECK_SEQUENCE(dist_r2.intervals(), list_of(0.0)(1.0));
+    CHECK_SEQUENCE(dist_r2.densities(), list_of(1.0)(1.0));
+    
+    boost::random::piecewise_linear_distribution<> dist_it(
+        intervals.begin(), intervals.end(), weights.begin());
+    CHECK_SEQUENCE(dist_it.intervals(), list_of(0)(1)(2)(3)(5));
+    CHECK_SEQUENCE(dist_it.densities(), list_of(.375)(.125)(.25)(.125)(.25));
+    
+    boost::random::piecewise_linear_distribution<> dist_it2(
+        intervals2.begin(), intervals2.end(), weights2.begin());
+    CHECK_SEQUENCE(dist_it2.intervals(), list_of(0.0)(1.0));
+    CHECK_SEQUENCE(dist_it2.densities(), list_of(1.0)(1.0));
+    
+    boost::random::piecewise_linear_distribution<> dist_fun(4, 99,115, gen());
+    CHECK_SEQUENCE(dist_fun.intervals(), list_of(99)(103)(107)(111)(115));
+    CHECK_SEQUENCE(dist_fun.densities(),
+        list_of(.09375)(.03125)(0.0625)(.03125)(.15625));
+    
+    boost::random::piecewise_linear_distribution<>
+        dist_fun2(1, 99, 115, gen());
+    CHECK_SEQUENCE(dist_fun2.intervals(), list_of(99)(115));
+    CHECK_SEQUENCE(dist_fun2.densities(), list_of(0.046875)(0.078125));
+
+    boost::random::piecewise_linear_distribution<> copy(dist);
+    BOOST_CHECK_EQUAL(dist, copy);
+    boost::random::piecewise_linear_distribution<> copy_r(dist_r);
+    BOOST_CHECK_EQUAL(dist_r, copy_r);
+
+    boost::random::piecewise_linear_distribution<> notpow2(3, 99, 111, gen());
+    BOOST_REQUIRE_EQUAL(notpow2.densities().size(), 4u);
+    BOOST_CHECK_CLOSE_FRACTION(notpow2.densities()[0], 0.15, 1e-12);
+    BOOST_CHECK_CLOSE_FRACTION(notpow2.densities()[1], 0.05, 1e-12);
+    BOOST_CHECK_CLOSE_FRACTION(notpow2.densities()[2], 0.1, 1e-12);
+    BOOST_CHECK_CLOSE_FRACTION(notpow2.densities()[3], 0.05, 1e-12);
+    boost::random::piecewise_linear_distribution<> copy_notpow2(notpow2);
+    BOOST_CHECK_EQUAL(notpow2, copy_notpow2);
+}
+
+BOOST_AUTO_TEST_CASE(test_param) {
+    std::vector<double> intervals = boost::assign::list_of(0)(1)(2)(3)(5);
+    std::vector<double> weights = boost::assign::list_of(3)(1)(2)(1)(2);
+    std::vector<double> intervals2 = boost::assign::list_of(99);
+    std::vector<double> weights2 = boost::assign::list_of(2);
+    boost::random::piecewise_linear_distribution<> dist(intervals, weights);
+    boost::random::piecewise_linear_distribution<>::param_type
+        param = dist.param();
+    CHECK_SEQUENCE(param.intervals(), list_of(0)(1)(2)(3)(5));
+    CHECK_SEQUENCE(param.densities(), list_of(.375)(.125)(.25)(.125)(.25));
+    boost::random::piecewise_linear_distribution<> copy1(param);
+    BOOST_CHECK_EQUAL(dist, copy1);
+    boost::random::piecewise_linear_distribution<> copy2;
+    copy2.param(param);
+    BOOST_CHECK_EQUAL(dist, copy2);
+
+    boost::random::piecewise_linear_distribution<>::param_type
+        param_copy = param;
+    BOOST_CHECK_EQUAL(param, param_copy);
+    BOOST_CHECK(param == param_copy);
+    BOOST_CHECK(!(param != param_copy));
+    boost::random::piecewise_linear_distribution<>::param_type param_default;
+    CHECK_SEQUENCE(param_default.intervals(), list_of(0.0)(1.0));
+    CHECK_SEQUENCE(param_default.densities(), list_of(1.0)(1.0));
+    BOOST_CHECK(param != param_default);
+    BOOST_CHECK(!(param == param_default));
+    
+#ifndef BOOST_NO_INITIALIZER_LISTS
+    boost::random::piecewise_linear_distribution<>::param_type parm_il = {
+        { 99, 103, 107, 111, 115 },
+        gen()
+    };
+    CHECK_SEQUENCE(parm_il.intervals(), list_of(99)(103)(107)(111)(115));
+    CHECK_SEQUENCE(parm_il.densities(),
+        list_of(.09375)(.03125)(0.0625)(.03125)(.15625));
+
+    boost::random::piecewise_linear_distribution<>::param_type parm_il2 = {
+        { 99 },
+        gen()
+    };
+    CHECK_SEQUENCE(parm_il2.intervals(), list_of(0.0)(1.0));
+    CHECK_SEQUENCE(parm_il2.densities(), list_of(1.0)(1.0));
+#endif
+
+    boost::random::piecewise_linear_distribution<>::param_type
+        parm_r(intervals, weights);
+    CHECK_SEQUENCE(parm_r.intervals(), list_of(0)(1)(2)(3)(5));
+    CHECK_SEQUENCE(parm_r.densities(), list_of(.375)(.125)(.25)(.125)(.25));
+
+    boost::random::piecewise_linear_distribution<>::param_type
+        parm_r2(intervals2, weights2);
+    CHECK_SEQUENCE(parm_r2.intervals(), list_of(0.0)(1.0));
+    CHECK_SEQUENCE(parm_r2.densities(), list_of(1.0)(1.0));
+    
+    boost::random::piecewise_linear_distribution<>::param_type
+        parm_it(intervals.begin(), intervals.end(), weights.begin());
+    CHECK_SEQUENCE(parm_it.intervals(), list_of(0)(1)(2)(3)(5));
+    CHECK_SEQUENCE(parm_it.densities(), list_of(.375)(.125)(.25)(.125)(.25));
+    
+    boost::random::piecewise_linear_distribution<>::param_type
+        parm_it2(intervals2.begin(), intervals2.end(), weights2.begin());
+    CHECK_SEQUENCE(parm_it2.intervals(), list_of(0.0)(1.0));
+    CHECK_SEQUENCE(parm_it2.densities(), list_of(1.0)(1.0));
+    
+    boost::random::piecewise_linear_distribution<>::param_type
+        parm_fun(4, 99, 115, gen());
+    CHECK_SEQUENCE(parm_fun.intervals(), list_of(99)(103)(107)(111)(115));
+    CHECK_SEQUENCE(parm_fun.densities(),
+        list_of(.09375)(.03125)(0.0625)(.03125)(.15625));
+    
+    boost::random::piecewise_linear_distribution<>::param_type
+        parm_fun2(1, 99, 115, gen());
+    CHECK_SEQUENCE(parm_fun2.intervals(), list_of(99)(115));
+    CHECK_SEQUENCE(parm_fun2.densities(), list_of(0.046875)(0.078125));
+}
+
+BOOST_AUTO_TEST_CASE(test_min_max) {
+    std::vector<double> intervals = boost::assign::list_of(0)(1)(2)(3)(5);
+    std::vector<double> weights = boost::assign::list_of(3)(1)(2)(1)(2);
+    boost::random::piecewise_linear_distribution<> dist;
+    BOOST_CHECK_EQUAL((dist.min)(), 0.0);
+    BOOST_CHECK_EQUAL((dist.max)(), 1.0);
+    boost::random::piecewise_linear_distribution<> dist_r(intervals, weights);
+    BOOST_CHECK_EQUAL((dist_r.min)(), 0.0);
+    BOOST_CHECK_EQUAL((dist_r.max)(), 5.0);
+}
+
+BOOST_AUTO_TEST_CASE(test_comparison) {
+    std::vector<double> intervals = boost::assign::list_of(0)(1)(2)(3)(5);
+    std::vector<double> weights = boost::assign::list_of(3)(1)(2)(1)(2);
+    boost::random::piecewise_linear_distribution<> dist;
+    boost::random::piecewise_linear_distribution<> dist_copy(dist);
+    boost::random::piecewise_linear_distribution<> dist_r(intervals, weights);
+    boost::random::piecewise_linear_distribution<> dist_r_copy(dist_r);
+    BOOST_CHECK(dist == dist_copy);
+    BOOST_CHECK(!(dist != dist_copy));
+    BOOST_CHECK(dist_r == dist_r_copy);
+    BOOST_CHECK(!(dist_r != dist_r_copy));
+    BOOST_CHECK(dist != dist_r);
+    BOOST_CHECK(!(dist == dist_r));
+}
+
+BOOST_AUTO_TEST_CASE(test_streaming) {
+    std::vector<double> intervals = boost::assign::list_of(0)(1)(2)(3)(5);
+    std::vector<double> weights = boost::assign::list_of(3)(1)(2)(1)(2);
+    boost::random::piecewise_linear_distribution<> dist(intervals, weights);
+    std::stringstream stream;
+    stream << dist;
+    boost::random::piecewise_linear_distribution<> restored_dist;
+    stream >> restored_dist;
+    BOOST_CHECK_EQUAL(dist, restored_dist);
+}
+
+BOOST_AUTO_TEST_CASE(test_generation) {
+    std::vector<double> intervals = boost::assign::list_of(1)(2);
+    std::vector<double> weights = boost::assign::list_of(1)(1);
+    boost::minstd_rand0 gen;
+    boost::random::piecewise_linear_distribution<> dist;
+    boost::random::piecewise_linear_distribution<> dist_r(intervals, weights);
+    for(int i = 0; i < 10; ++i) {
+        double value = dist(gen);
+        BOOST_CHECK_GE(value, 0.0);
+        BOOST_CHECK_LT(value, 1.0);
+        double value_r = dist_r(gen);
+        BOOST_CHECK_GE(value_r, 1.0);
+        BOOST_CHECK_LT(value_r, 2.0);
+        double value_param = dist_r(gen, dist.param());
+        BOOST_CHECK_GE(value_param, 0.0);
+        BOOST_CHECK_LT(value_param, 1.0);
+        double value_r_param = dist(gen, dist_r.param());
+        BOOST_CHECK_GE(value_r_param, 1.0);
+        BOOST_CHECK_LT(value_r_param, 2.0);
+    }
+}