summaryrefslogtreecommitdiff
path: root/boost/random
diff options
context:
space:
mode:
Diffstat (limited to 'boost/random')
-rw-r--r--boost/random/binomial_distribution.hpp54
-rw-r--r--boost/random/detail/vector_io.hpp3
-rw-r--r--boost/random/discrete_distribution.hpp4
-rw-r--r--boost/random/hyperexponential_distribution.hpp11
-rw-r--r--boost/random/independent_bits.hpp12
-rw-r--r--boost/random/linear_congruential.hpp6
-rw-r--r--boost/random/poisson_distribution.hpp34
-rw-r--r--boost/random/seed_seq.hpp2
-rw-r--r--boost/random/uniform_on_sphere.hpp8
9 files changed, 73 insertions, 61 deletions
diff --git a/boost/random/binomial_distribution.hpp b/boost/random/binomial_distribution.hpp
index 3efc905746..78d1a123a4 100644
--- a/boost/random/binomial_distribution.hpp
+++ b/boost/random/binomial_distribution.hpp
@@ -278,18 +278,18 @@ private:
m = static_cast<IntType>((t+1)*p);
if(use_inversion()) {
- q_n = pow((1 - p), static_cast<RealType>(t));
+ _u.q_n = pow((1 - p), static_cast<RealType>(t));
} else {
- btrd.r = p/(1-p);
- btrd.nr = (t+1)*btrd.r;
- btrd.npq = t*p*(1-p);
- RealType sqrt_npq = sqrt(btrd.npq);
- btrd.b = 1.15 + 2.53 * sqrt_npq;
- btrd.a = -0.0873 + 0.0248*btrd.b + 0.01*p;
- btrd.c = t*p + 0.5;
- btrd.alpha = (2.83 + 5.1/btrd.b) * sqrt_npq;
- btrd.v_r = 0.92 - 4.2/btrd.b;
- btrd.u_rv_r = 0.86*btrd.v_r;
+ _u.btrd.r = p/(1-p);
+ _u.btrd.nr = (t+1)*_u.btrd.r;
+ _u.btrd.npq = t*p*(1-p);
+ RealType sqrt_npq = sqrt(_u.btrd.npq);
+ _u.btrd.b = 1.15 + 2.53 * sqrt_npq;
+ _u.btrd.a = -0.0873 + 0.0248*_u.btrd.b + 0.01*p;
+ _u.btrd.c = t*p + 0.5;
+ _u.btrd.alpha = (2.83 + 5.1/_u.btrd.b) * sqrt_npq;
+ _u.btrd.v_r = 0.92 - 4.2/_u.btrd.b;
+ _u.btrd.u_rv_r = 0.86*_u.btrd.v_r;
}
}
@@ -303,24 +303,24 @@ private:
while(true) {
RealType u;
RealType v = uniform_01<RealType>()(urng);
- if(v <= btrd.u_rv_r) {
- u = v/btrd.v_r - 0.43;
+ if(v <= _u.btrd.u_rv_r) {
+ u = v/_u.btrd.v_r - 0.43;
return static_cast<IntType>(floor(
- (2*btrd.a/(0.5 - abs(u)) + btrd.b)*u + btrd.c));
+ (2*_u.btrd.a/(0.5 - abs(u)) + _u.btrd.b)*u + _u.btrd.c));
}
- if(v >= btrd.v_r) {
+ if(v >= _u.btrd.v_r) {
u = uniform_01<RealType>()(urng) - 0.5;
} else {
- u = v/btrd.v_r - 0.93;
+ u = v/_u.btrd.v_r - 0.93;
u = ((u < 0)? -0.5 : 0.5) - u;
- v = uniform_01<RealType>()(urng) * btrd.v_r;
+ v = uniform_01<RealType>()(urng) * _u.btrd.v_r;
}
RealType us = 0.5 - abs(u);
- IntType k = static_cast<IntType>(floor((2*btrd.a/us + btrd.b)*u + btrd.c));
+ IntType k = static_cast<IntType>(floor((2*_u.btrd.a/us + _u.btrd.b)*u + _u.btrd.c));
if(k < 0 || k > _t) continue;
- v = v*btrd.alpha/(btrd.a/(us*us) + btrd.b);
+ v = v*_u.btrd.alpha/(_u.btrd.a/(us*us) + _u.btrd.b);
RealType km = abs(k - m);
if(km <= 15) {
RealType f = 1;
@@ -328,13 +328,13 @@ private:
IntType i = m;
do {
++i;
- f = f*(btrd.nr/i - btrd.r);
+ f = f*(_u.btrd.nr/i - _u.btrd.r);
} while(i != k);
} else if(m > k) {
IntType i = k;
do {
++i;
- v = v*(btrd.nr/i - btrd.r);
+ v = v*(_u.btrd.nr/i - _u.btrd.r);
} while(i != m);
}
if(v <= f) return k;
@@ -343,18 +343,18 @@ private:
// final acceptance/rejection
v = log(v);
RealType rho =
- (km/btrd.npq)*(((km/3. + 0.625)*km + 1./6)/btrd.npq + 0.5);
- RealType t = -km*km/(2*btrd.npq);
+ (km/_u.btrd.npq)*(((km/3. + 0.625)*km + 1./6)/_u.btrd.npq + 0.5);
+ RealType t = -km*km/(2*_u.btrd.npq);
if(v < t - rho) return k;
if(v > t + rho) continue;
IntType nm = _t - m + 1;
- RealType h = (m + 0.5)*log((m + 1)/(btrd.r*nm))
+ RealType h = (m + 0.5)*log((m + 1)/(_u.btrd.r*nm))
+ fc(m) + fc(_t - m);
IntType nk = _t - k + 1;
if(v <= h + (_t+1)*log(static_cast<RealType>(nm)/nk)
- + (k + 0.5)*log(nk*btrd.r/(k+1))
+ + (k + 0.5)*log(nk*_u.btrd.r/(k+1))
- fc(k)
- fc(_t - k))
{
@@ -372,7 +372,7 @@ private:
RealType q = 1 - p;
RealType s = p / q;
RealType a = (t + 1) * s;
- RealType r = q_n;
+ RealType r = _u.q_n;
RealType u = uniform_01<RealType>()(urng);
IntType x = 0;
while(u > r) {
@@ -417,7 +417,7 @@ private:
} btrd;
// for inversion
RealType q_n;
- };
+ } _u;
/// @endcond
};
diff --git a/boost/random/detail/vector_io.hpp b/boost/random/detail/vector_io.hpp
index 24508c210e..fe3869abb5 100644
--- a/boost/random/detail/vector_io.hpp
+++ b/boost/random/detail/vector_io.hpp
@@ -16,6 +16,7 @@
#include <vector>
#include <iosfwd>
#include <istream>
+#include <boost/io/ios_state.hpp>
namespace boost {
namespace random {
@@ -52,12 +53,14 @@ void read_vector(std::basic_istream<CharT, Traits>& is, std::vector<T>& vec)
is.setstate(std::ios_base::failbit);
return;
}
+ boost::io::basic_ios_exception_saver<CharT, Traits> e(is, std::ios_base::goodbit);
T val;
while(is >> std::ws >> val) {
vec.push_back(val);
}
if(is.fail()) {
is.clear();
+ e.restore();
if(!(is >> ch)) {
return;
}
diff --git a/boost/random/discrete_distribution.hpp b/boost/random/discrete_distribution.hpp
index 3d14d4a02c..9574cbcae5 100644
--- a/boost/random/discrete_distribution.hpp
+++ b/boost/random/discrete_distribution.hpp
@@ -92,7 +92,7 @@ struct integer_alias_table {
return _alias_table == other._alias_table &&
_average == other._average && _excess == other._excess;
}
- static WeightType normalize(WeightType val, WeightType average)
+ static WeightType normalize(WeightType val, WeightType /* average */)
{
return val;
}
@@ -183,7 +183,7 @@ struct real_alias_table {
{
return true;
}
- static WeightType try_get_sum(const std::vector<WeightType>& weights)
+ static WeightType try_get_sum(const std::vector<WeightType>& /* weights */)
{
return static_cast<WeightType>(1);
}
diff --git a/boost/random/hyperexponential_distribution.hpp b/boost/random/hyperexponential_distribution.hpp
index 046d0fc401..a2833e8a9a 100644
--- a/boost/random/hyperexponential_distribution.hpp
+++ b/boost/random/hyperexponential_distribution.hpp
@@ -490,14 +490,10 @@ class hyperexponential_distribution
rates.clear();
}
- bool fail = false;
-
// Adjust vector sizes (if needed)
if (param.probs_.size() != param.rates_.size()
|| param.probs_.size() == 0)
{
- fail = true;
-
const std::size_t np = param.probs_.size();
const std::size_t nr = param.rates_.size();
@@ -521,13 +517,6 @@ class hyperexponential_distribution
// can be changed due to size conformance
hyperexp_detail::normalize(param.probs_);
- // Set the error state in the underlying stream in case of invalid input
- if (fail)
- {
- // This throws an exception if ios_base::exception(failbit) is enabled
- is.setstate(std::ios_base::failbit);
- }
-
//post: vector size conformance
assert(param.probs_.size() == param.rates_.size());
diff --git a/boost/random/independent_bits.hpp b/boost/random/independent_bits.hpp
index dec63b3f46..cd2dd1f225 100644
--- a/boost/random/independent_bits.hpp
+++ b/boost/random/independent_bits.hpp
@@ -158,6 +158,18 @@ public:
BOOST_ASSERT(n0*w0 + (n - n0)*(w0 + 1) == w);
+ BOOST_ASSERT((n == 1) == (w0 == w));
+
+ // special case to avoid undefined behavior from shifting
+ if(n == 1) {
+ BOOST_ASSERT(n0 == 1);
+ base_unsigned u;
+ do {
+ u = detail::subtract<base_result_type>()(_base(), (_base.min)());
+ } while(u > base_unsigned(y0 - 1));
+ return u & y0_mask;
+ }
+
result_type S = 0;
for(std::size_t k = 0; k < n0; ++k) {
base_unsigned u;
diff --git a/boost/random/linear_congruential.hpp b/boost/random/linear_congruential.hpp
index de3a1d0749..22ade3f3db 100644
--- a/boost/random/linear_congruential.hpp
+++ b/boost/random/linear_congruential.hpp
@@ -124,8 +124,12 @@ public:
* distinct seeds in the range [1,m) will leave the generator in distinct
* states. If c is not zero, the range is [0,m).
*/
- BOOST_RANDOM_DETAIL_ARITHMETIC_SEED(linear_congruential_engine, IntType, x0)
+ BOOST_RANDOM_DETAIL_ARITHMETIC_SEED(linear_congruential_engine, IntType, x0_)
{
+ // Work around a msvc 12/14 optimizer bug, which causes
+ // the line _x = 1 to run unconditionally sometimes.
+ // Creating a local copy seems to make it work.
+ IntType x0 = x0_;
// wrap _x if it doesn't fit in the destination
if(modulus == 0) {
_x = x0;
diff --git a/boost/random/poisson_distribution.hpp b/boost/random/poisson_distribution.hpp
index 759f206e42..1281a7d859 100644
--- a/boost/random/poisson_distribution.hpp
+++ b/boost/random/poisson_distribution.hpp
@@ -255,13 +255,13 @@ private:
using std::exp;
if(use_inversion()) {
- _exp_mean = exp(-_mean);
+ _u._exp_mean = exp(-_mean);
} else {
- _ptrd.smu = sqrt(_mean);
- _ptrd.b = 0.931 + 2.53 * _ptrd.smu;
- _ptrd.a = -0.059 + 0.02483 * _ptrd.b;
- _ptrd.inv_alpha = 1.1239 + 1.1328 / (_ptrd.b - 3.4);
- _ptrd.v_r = 0.9277 - 3.6224 / (_ptrd.b - 2);
+ _u._ptrd.smu = sqrt(_mean);
+ _u._ptrd.b = 0.931 + 2.53 * _u._ptrd.smu;
+ _u._ptrd.a = -0.059 + 0.02483 * _u._ptrd.b;
+ _u._ptrd.inv_alpha = 1.1239 + 1.1328 / (_u._ptrd.b - 3.4);
+ _u._ptrd.v_r = 0.9277 - 3.6224 / (_u._ptrd.b - 2);
}
}
@@ -275,18 +275,18 @@ private:
while(true) {
RealType u;
RealType v = uniform_01<RealType>()(urng);
- if(v <= 0.86 * _ptrd.v_r) {
- u = v / _ptrd.v_r - 0.43;
+ if(v <= 0.86 * _u._ptrd.v_r) {
+ u = v / _u._ptrd.v_r - 0.43;
return static_cast<IntType>(floor(
- (2*_ptrd.a/(0.5-abs(u)) + _ptrd.b)*u + _mean + 0.445));
+ (2*_u._ptrd.a/(0.5-abs(u)) + _u._ptrd.b)*u + _mean + 0.445));
}
- if(v >= _ptrd.v_r) {
+ if(v >= _u._ptrd.v_r) {
u = uniform_01<RealType>()(urng) - 0.5;
} else {
- u = v/_ptrd.v_r - 0.93;
+ u = v/_u._ptrd.v_r - 0.93;
u = ((u < 0)? -0.5 : 0.5) - u;
- v = uniform_01<RealType>()(urng) * _ptrd.v_r;
+ v = uniform_01<RealType>()(urng) * _u._ptrd.v_r;
}
RealType us = 0.5 - abs(u);
@@ -294,13 +294,13 @@ private:
continue;
}
- RealType k = floor((2*_ptrd.a/us + _ptrd.b)*u+_mean+0.445);
- v = v*_ptrd.inv_alpha/(_ptrd.a/(us*us) + _ptrd.b);
+ RealType k = floor((2*_u._ptrd.a/us + _u._ptrd.b)*u+_mean+0.445);
+ v = v*_u._ptrd.inv_alpha/(_u._ptrd.a/(us*us) + _u._ptrd.b);
RealType log_sqrt_2pi = 0.91893853320467267;
if(k >= 10) {
- if(log(v*_ptrd.smu) <= (k + 0.5)*log(_mean/k)
+ if(log(v*_u._ptrd.smu) <= (k + 0.5)*log(_mean/k)
- _mean
- log_sqrt_2pi
+ k
@@ -320,7 +320,7 @@ private:
template<class URNG>
IntType invert(URNG& urng) const
{
- RealType p = _exp_mean;
+ RealType p = _u._exp_mean;
IntType x = 0;
RealType u = uniform_01<RealType>()(urng);
while(u > p) {
@@ -344,7 +344,7 @@ private:
} _ptrd;
// for inversion
RealType _exp_mean;
- };
+ } _u;
/// @endcond
};
diff --git a/boost/random/seed_seq.hpp b/boost/random/seed_seq.hpp
index d76aef4f5e..b0b1de24e3 100644
--- a/boost/random/seed_seq.hpp
+++ b/boost/random/seed_seq.hpp
@@ -97,7 +97,7 @@ public:
& mask);
r3 = r3 ^ (r3 >> 27);
r3 = (r3 * 1566083941u) & mask;
- value_type r4 = static_cast<value_type>(r3 - k%m);
+ value_type r4 = static_cast<value_type>(r3 - k%n);
*(first + (k+p)%n) ^= r3;
*(first + (k+q)%n) ^= r4;
*(first + k%n) = r4;
diff --git a/boost/random/uniform_on_sphere.hpp b/boost/random/uniform_on_sphere.hpp
index 72c25ef781..ce2e35237e 100644
--- a/boost/random/uniform_on_sphere.hpp
+++ b/boost/random/uniform_on_sphere.hpp
@@ -225,8 +225,12 @@ public:
}
} while(sqsum == 0);
// for all i: result[i] /= sqrt(sqsum)
- std::transform(_container.begin(), _container.end(), _container.begin(),
- std::bind2nd(std::multiplies<RealType>(), 1/sqrt(sqsum)));
+ RealType inverse_distance = 1 / sqrt(sqsum);
+ for(typename Cont::iterator it = _container.begin();
+ it != _container.end();
+ ++it) {
+ *it *= inverse_distance;
+ }
}
}
return _container;