Line data Source code
1 : /** \file randomT.hpp
2 : * \author Jared R. Males
3 : * \brief Defines a random number type
4 : * \ingroup gen_math_files
5 : *
6 : */
7 :
8 : //***********************************************************************//
9 : // Copyright 2015, 2016, 2017, 2020 Jared R. Males (jaredmales@gmail.com)
10 : //
11 : // This file is part of mxlib.
12 : //
13 : // mxlib is free software: you can redistribute it and/or modify
14 : // it under the terms of the GNU General Public License as published by
15 : // the Free Software Foundation, either version 3 of the License, or
16 : // (at your option) any later version.
17 : //
18 : // mxlib is distributed in the hope that it will be useful,
19 : // but WITHOUT ANY WARRANTY; without even the implied warranty of
20 : // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
21 : // GNU General Public License for more details.
22 : //
23 : // You should have received a copy of the GNU General Public License
24 : // along with mxlib. If not, see <http://www.gnu.org/licenses/>.
25 : //***********************************************************************//
26 :
27 : #ifndef mx_math_randomT_hpp
28 : #define mx_math_randomT_hpp
29 :
30 : #include <random>
31 :
32 : #include "randomSeed.hpp"
33 :
34 : namespace mx
35 : {
36 : namespace math
37 : {
38 :
39 : /// A random number type, which functions like any other arithmetic type.
40 : /**
41 : * Combines a random engine, and a random distribution. Using the type conversion operator
42 : * randomT returns the next random deviate whenever it is referenced.
43 : *
44 : * Example:
45 : *
46 : \code
47 : //This can also be done using the alias definition mx::math::normDistT.
48 : randomT<double, std::mt19937_64, std::normal_distribution<double> > norm_distd;
49 : norm_distd.seed(); //
50 :
51 : double d1 = norm_distd; //get a normally distributed value
52 : double d2 = norm_distd; //get the next normally distributed value
53 : \endcode
54 : *
55 : * \ingroup random
56 : */
57 : template <class typeT, class _ranengT, class _randistT>
58 : class randomT
59 : {
60 :
61 : public:
62 : /// Typedef for the distribution type
63 : typedef _randistT randistT;
64 :
65 : /// Typedef for the engine type
66 : typedef _ranengT ranengT;
67 :
68 : /// Constructor
69 : /** By default this calls the seed method, which will use /dev/random to seed the generator on linux, and time(0) on
70 : * other systems. Set to false to suppress seeding, and/or set a seed with seed(x).
71 : */
72 30 : randomT( bool doSeed = true /**< [in] [optional] if true then the seed method is called upon construction.*/ )
73 35 : {
74 30 : if( doSeed )
75 : {
76 9 : seed();
77 : }
78 30 : }
79 :
80 : /// The random distribution
81 : _randistT distribution;
82 :
83 : /// The random engine
84 : ranengT engine;
85 :
86 : /// The conversion operator, returns the next value in the sequence, according to the distribution.
87 14578298 : operator typeT()
88 : {
89 14578298 : return distribution( engine );
90 : }
91 :
92 : /// Set the seed of the random engine.
93 : /** Calls the engines seed member function.
94 : *
95 : */
96 42 : void seed( typename ranengT::result_type seedval /**< [in] the argument to pass to ranengT::seed() */ )
97 : {
98 42 : engine.seed( seedval );
99 42 : }
100 :
101 : /// Seed the random engine with a good value
102 : /** Calls \ref mx::math::randomSeed to get the value. On Linux this uses `/dev/urandom`. On other systems, this uses
103 : * `time(0)`.
104 : */
105 18 : void seed()
106 : {
107 : typename ranengT::result_type seedval;
108 18 : randomSeed( seedval );
109 18 : seed( seedval );
110 18 : }
111 : };
112 :
113 : /**
114 : * @brief The Laplace (double exponential) continuous distribution for random numbers.
115 : *
116 : * The formula for the exponential probability density function is
117 : * @f$p(x|\lambda) = \frac{\lambda}{2} e^{-\lambda x}@f$.
118 : *
119 : * <table border=1 cellpadding=10 cellspacing=0>
120 : * <caption align=top>Distribution Statistics</caption>
121 : * <tr><td>Mean</td><td>@f$0@f$</td></tr>
122 : * <tr><td>Median</td><td>@f$0@f$</td></tr>
123 : * <tr><td>Mode</td><td>@f$0@f$</td></tr>
124 : * <tr><td>Range</td><td>@f$[-\infty, \infty]@f$</td></tr>
125 : * <tr><td>Standard Deviation</td><td>@f$\frac{\sqrt{2}}{\lambda}@f$</td></tr>
126 : * </table>
127 : *
128 : * This is based on the implementation of the exponential distribution in the GNU ISO C++ Library version 4.6.
129 : *
130 : * \ingroup random
131 : */
132 : template <typename _RealType = double>
133 : class laplace_distribution
134 : {
135 : static_assert( std::is_floating_point<_RealType>::value, "template argument not a floating point type" );
136 :
137 : public:
138 : /** The type of the range of the distribution. */
139 : typedef _RealType result_type;
140 :
141 : /** Parameter type. */
142 : struct param_type
143 : {
144 : typedef laplace_distribution<_RealType> distribution_type;
145 :
146 5 : explicit param_type( _RealType __lambda = _RealType( 1 ) ) : _M_lambda( __lambda )
147 : {
148 5 : }
149 :
150 6 : _RealType lambda() const
151 : {
152 6 : return _M_lambda;
153 : }
154 :
155 : friend bool operator==( const param_type &__p1, const param_type &__p2 )
156 : {
157 : return __p1._M_lambda == __p2._M_lambda;
158 : }
159 :
160 : private:
161 : _RealType _M_lambda;
162 : };
163 :
164 : public:
165 : /**
166 : * @brief Constructs a Laplace (double exponential) distribution with inverse scale
167 : * parameter @f$\lambda@f$.
168 : */
169 10 : explicit laplace_distribution( const result_type &__lambda = result_type( 1 ) ) : _M_param( __lambda )
170 : {
171 5 : }
172 :
173 : explicit laplace_distribution( const param_type &__p ) : _M_param( __p )
174 : {
175 : }
176 :
177 : /**
178 : * @brief Resets the distribution state.
179 : *
180 : * Has no effect on Laplace distributions.
181 : */
182 : void reset()
183 : {
184 : }
185 :
186 : /**
187 : * @brief Returns the inverse scale parameter of the distribution.
188 : */
189 : _RealType lambda() const
190 : {
191 : return _M_param.lambda();
192 : }
193 :
194 : /**
195 : * @brief Returns the parameter set of the distribution.
196 : */
197 6 : param_type param() const
198 : {
199 6 : return _M_param;
200 : }
201 :
202 : /**
203 : * @brief Sets the parameter set of the distribution.
204 : * @param __param The new parameter set of the distribution.
205 : */
206 : void param( const param_type &__param )
207 : {
208 : _M_param = __param;
209 : }
210 :
211 : /**
212 : * @brief Returns the greatest lower bound value of the distribution.
213 : */
214 : result_type min() const
215 : {
216 : return std::numeric_limits<result_type>::min();
217 : }
218 :
219 : /**
220 : * @brief Returns the least upper bound value of the distribution.
221 : */
222 : result_type max() const
223 : {
224 : return std::numeric_limits<result_type>::max();
225 : }
226 :
227 : /**
228 : * @brief Generating functions.
229 : */
230 : template <typename _UniformRandomNumberGenerator>
231 6 : result_type operator()( _UniformRandomNumberGenerator &__urng )
232 : {
233 6 : return this->operator()( __urng, this->param() );
234 : }
235 :
236 : template <typename _UniformRandomNumberGenerator>
237 6 : result_type operator()( _UniformRandomNumberGenerator &__urng, const param_type &__p )
238 : {
239 6 : int sgnx = 1;
240 :
241 6 : result_type __aurng = 0.5 - std::generate_canonical<result_type,
242 : std::numeric_limits<result_type>::digits,
243 6 : _UniformRandomNumberGenerator>( __urng );
244 :
245 6 : if( __aurng < 0 )
246 : {
247 5 : sgnx = -1;
248 5 : __aurng *= -1;
249 : }
250 :
251 6 : return 0. - __p.lambda() * sgnx * std::log( 1. - 2. * __aurng );
252 : }
253 :
254 : private:
255 : param_type _M_param;
256 : };
257 :
258 : /**
259 : * @brief Return true if two exponential distributions have the same
260 : * parameters.
261 : */
262 : template <typename _RealType>
263 : bool operator==( const laplace_distribution<_RealType> &__d1, const laplace_distribution<_RealType> &__d2 )
264 : {
265 : return __d1.param() == __d2.param();
266 : }
267 :
268 : /**
269 : * @brief Return true if two exponential distributions have different
270 : * parameters.
271 : */
272 : template <typename _RealType>
273 : bool operator!=( const laplace_distribution<_RealType> &__d1, const laplace_distribution<_RealType> &__d2 )
274 : {
275 : return !( __d1 == __d2 );
276 : }
277 :
278 : /**
279 : * @brief Inserts a %laplace_distribution random number distribution
280 : * @p __x into the output stream @p __os.
281 : *
282 : * @param __os An output stream.
283 : * @param __x A %laplace_distribution random number distribution.
284 : *
285 : * @returns The output stream with the state of @p __x inserted or in
286 : * an error state.
287 : */
288 : template <typename _RealType, typename _CharT, typename _Traits>
289 : std::basic_ostream<_CharT, _Traits> &operator<<( std::basic_ostream<_CharT, _Traits> &,
290 : const laplace_distribution<_RealType> & );
291 :
292 : /**
293 : * @brief Extracts a %laplace_distribution random number distribution
294 : * @p __x from the input stream @p __is.
295 : *
296 : * @param __is An input stream.
297 : * @param __x A %laplace_distribution random number
298 : * generator engine.
299 : *
300 : * @returns The input stream with @p __x extracted or in an error state.
301 : */
302 : template <typename _RealType, typename _CharT, typename _Traits>
303 : std::basic_istream<_CharT, _Traits> &operator>>( std::basic_istream<_CharT, _Traits> &,
304 : laplace_distribution<_RealType> & );
305 :
306 : /** \ingroup random
307 : * @{
308 : */
309 :
310 : /// Alias for a uniform random variate
311 : template <typename realT>
312 : using uniDistT = randomT<realT, std::mt19937_64, std::uniform_real_distribution<realT>>;
313 :
314 : /// Alias for a standard normal random variate
315 : template <typename realT>
316 : using normDistT = randomT<realT, std::mt19937_64, std::normal_distribution<realT>>;
317 :
318 : /// Alias for an exponential random variate
319 : template <typename realT>
320 : using expDistT = randomT<realT, std::mt19937_64, std::exponential_distribution<realT>>;
321 :
322 : /// Alias for a laplace random variate
323 : template <typename realT>
324 : using lapDistT = randomT<realT, std::mt19937_64, laplace_distribution<realT>>;
325 :
326 : /// Alias for a poisson random variate
327 : template <typename intT>
328 : using poissonDistT = randomT<intT, std::mt19937_64, std::poisson_distribution<intT>>;
329 :
330 : /// Alias for a log normal variate
331 : template <typename realT>
332 : using lognormDistT = randomT<realT, std::mt19937_64, std::lognormal_distribution<realT>>;
333 :
334 : /// @}
335 :
336 : } // namespace math
337 : } // namespace mx
338 :
339 : #endif // mx_math_randomT_hpp
340 :
341 : /*5/22/06: added 64 bit mersenne twister support
342 : */
|
