/** * @file NormalExp.h * @author bab2min (bab2min@gmail.com) * @brief * @version 0.3.3 * @date 2021-03-31 * * @copyright Copyright (c) 2020-2021 * */ #ifndef EIGENRAND_DISTS_NORMAL_EXP_H #define EIGENRAND_DISTS_NORMAL_EXP_H namespace Eigen { namespace Rand { /** * @brief Generator of reals on the standard normal distribution * * @tparam _Scalar */ template class StdNormalGen : OptCacheStore, public GenBase, _Scalar> { static_assert(std::is_floating_point<_Scalar>::value, "normalDist needs floating point types."); bool valid = false; StdUniformRealGen<_Scalar> ur; public: using Scalar = _Scalar; template EIGEN_STRONG_INLINE const _Scalar operator() (Rng&& rng) { using namespace Eigen::internal; if (valid) { valid = false; return OptCacheStore::get<_Scalar>(); } valid = true; _Scalar v1, v2, sx; for (int _i = 0; ; ++_i) { EIGENRAND_CHECK_INFINITY_LOOP(); v1 = 2 * ur(rng) - 1; v2 = 2 * ur(rng) - 1; sx = v1 * v1 + v2 * v2; if (sx && sx < 1) break; } _Scalar fx = std::sqrt((_Scalar)-2.0 * std::log(sx) / sx); OptCacheStore::get<_Scalar>() = fx * v2; return fx * v1; } template EIGEN_STRONG_INLINE const Packet packetOp(Rng&& rng) { using namespace Eigen::internal; if (valid) { valid = false; return OptCacheStore::template get(); } valid = true; Packet u1 = ur.template packetOp(rng), u2 = ur.template packetOp(rng); u1 = psub(pset1(1), u1); auto radius = psqrt(pmul(pset1(-2), plog(u1))); auto theta = pmul(pset1(2 * constant::pi), u2); Packet sintheta, costheta; psincos(theta, sintheta, costheta); OptCacheStore::template get() = pmul(radius, costheta); return pmul(radius, sintheta); } }; /** * @brief Generator of reals on a normal distribution * * @tparam _Scalar */ template class NormalGen : public GenBase, _Scalar> { static_assert(std::is_floating_point<_Scalar>::value, "normalDist needs floating point types."); StdNormalGen<_Scalar> stdnorm; _Scalar mean = 0, stdev = 1; public: using Scalar = _Scalar; /** * @brief Construct a new normal generator * * @param _mean mean of the distribution * @param _stdev standard deviation of the distribution */ NormalGen(_Scalar _mean = 0, _Scalar _stdev = 1) : mean{ _mean }, stdev{ _stdev } { } NormalGen(const NormalGen&) = default; NormalGen(NormalGen&&) = default; NormalGen& operator=(const NormalGen&) = default; NormalGen& operator=(NormalGen&&) = default; template EIGEN_STRONG_INLINE const _Scalar operator() (Rng&& rng) { using namespace Eigen::internal; return stdnorm(std::forward(rng)) * stdev + mean; } template EIGEN_STRONG_INLINE const Packet packetOp(Rng&& rng) { using namespace Eigen::internal; return padd(pmul( stdnorm.template packetOp(std::forward(rng)), pset1(stdev) ), pset1(mean)); } }; /** * @brief Generator of reals on a lognormal distribution * * @tparam _Scalar */ template class LognormalGen : public GenBase, _Scalar> { static_assert(std::is_floating_point<_Scalar>::value, "lognormalDist needs floating point types."); NormalGen<_Scalar> norm; public: using Scalar = _Scalar; /** * @brief Construct a new lognormal generator * * @param _mean mean of the log distribution * @param _stdev standard deviation of the log distribution */ LognormalGen(_Scalar _mean = 0, _Scalar _stdev = 1) : norm{ _mean, _stdev } { } LognormalGen(const LognormalGen&) = default; LognormalGen(LognormalGen&&) = default; LognormalGen& operator=(const LognormalGen&) = default; LognormalGen& operator=(LognormalGen&&) = default; template EIGEN_STRONG_INLINE const _Scalar operator() (Rng&& rng) { using namespace Eigen::internal; return std::exp(norm(std::forward(rng))); } template EIGEN_STRONG_INLINE const Packet packetOp(Rng&& rng) { using namespace Eigen::internal; return pexp(norm.template packetOp(std::forward(rng))); } }; /** * @brief Generator of reals on a Student's t distribution * * @tparam _Scalar */ template class StudentTGen : public GenBase, _Scalar> { static_assert(std::is_floating_point<_Scalar>::value, "studentT needs floating point types."); StdUniformRealGen<_Scalar> ur; _Scalar n; public: using Scalar = _Scalar; /** * @brief Construct a new Student's t generator * * @param _n degrees of freedom */ StudentTGen(_Scalar _n = 1) : n{ _n } { } StudentTGen(const StudentTGen&) = default; StudentTGen(StudentTGen&&) = default; StudentTGen& operator=(const StudentTGen&) = default; StudentTGen& operator=(StudentTGen&&) = default; template EIGEN_STRONG_INLINE const _Scalar operator() (Rng&& rng) { using namespace Eigen::internal; _Scalar v1, v2, sx; for (int _i = 0; ; ++_i) { EIGENRAND_CHECK_INFINITY_LOOP(); v1 = 2 * ur(rng) - 1; v2 = 2 * ur(rng) - 1; sx = v1 * v1 + v2 * v2; if (sx && sx < 1) break; } _Scalar fx = std::sqrt(n * (std::pow(sx, -2 / n) - 1) / sx); return fx * v1; } template EIGEN_STRONG_INLINE const Packet packetOp(Rng&& rng) { using namespace Eigen::internal; Packet u1 = ur.template packetOp(rng), u2 = ur.template packetOp(rng); u1 = psub(pset1(1), u1); auto pn = pset1(n); auto radius = psqrt(pmul(pn, psub(pexp(pmul(plog(u1), pset1(-2 / n))), pset1(1)) )); auto theta = pmul(pset1(2 * constant::pi), u2); //Packet sintheta, costheta; //psincos(theta, sintheta, costheta); return pmul(radius, psin(theta)); } }; template class GammaGen; /** * @brief Generator of reals on an exponential distribution * * @tparam _Scalar */ template class ExponentialGen : public GenBase, _Scalar> { friend GammaGen<_Scalar>; static_assert(std::is_floating_point<_Scalar>::value, "expDist needs floating point types."); StdUniformRealGen<_Scalar> ur; _Scalar lambda = 1; public: using Scalar = _Scalar; /** * @brief Construct a new exponential generator * * @param _lambda scale parameter of the distribution */ ExponentialGen(_Scalar _lambda = 1) : lambda{ _lambda } { } ExponentialGen(const ExponentialGen&) = default; ExponentialGen(ExponentialGen&&) = default; ExponentialGen& operator=(const ExponentialGen&) = default; ExponentialGen& operator=(ExponentialGen&&) = default; template EIGEN_STRONG_INLINE const _Scalar operator() (Rng&& rng) { using namespace Eigen::internal; return -std::log(1 - ur(std::forward(rng))) / lambda; } template EIGEN_STRONG_INLINE const Packet packetOp(Rng&& rng) { using namespace Eigen::internal; return pnegate(pdiv(plog( psub(pset1(1), ur.template packetOp(std::forward(rng))) ), pset1(lambda))); } }; template class NegativeBinomialGen; /** * @brief Generator of reals on a gamma distribution * * @tparam _Scalar */ template class GammaGen : OptCacheStore, public GenBase, _Scalar> { template friend class NegativeBinomialGen; static_assert(std::is_floating_point<_Scalar>::value, "gammaDist needs floating point types."); int cache_rest_cnt = 0; ExponentialGen<_Scalar> expon; _Scalar alpha, beta, px, sqrt; public: using Scalar = _Scalar; /** * @brief Construct a new gamma generator * * @param _alpha shape parameter of the distribution * @param _beta scale parameter of the distribution */ GammaGen(_Scalar _alpha = 1, _Scalar _beta = 1) : alpha{ _alpha }, beta{ _beta } { px = (_Scalar)(constant::e / (alpha + constant::e)); sqrt = std::sqrt(2 * alpha - 1); } GammaGen(const GammaGen&) = default; GammaGen(GammaGen&&) = default; GammaGen& operator=(const GammaGen&) = default; GammaGen& operator=(GammaGen&&) = default; template EIGEN_STRONG_INLINE const _Scalar operator() (Rng&& rng) { using namespace Eigen::internal; if (alpha < 1) { _Scalar ux, vx, xx, qx; for (int _i = 0; ; ++_i) { EIGENRAND_CHECK_INFINITY_LOOP(); ux = expon.ur(rng); vx = expon.ur.nzur_scalar(rng); if (ux < px) { xx = std::pow(vx, 1 / alpha); qx = std::exp(-xx); } else { xx = 1 - std::log(vx); qx = std::pow(xx, alpha - 1); } if (expon.ur(rng) < qx) { return beta * xx; } } } if (alpha == 1) { return beta * expon(rng); } int count; if ((count = alpha) == alpha && count < 20) { _Scalar yx; yx = expon.ur.nzur_scalar(rng); while (--count) { yx *= expon.ur.nzur_scalar(rng); } return -beta * std::log(yx); } for (int _i = 0; ; ++_i) { EIGENRAND_CHECK_INFINITY_LOOP(); _Scalar yx, xx; yx = std::tan(constant::pi * expon.ur(rng)); xx = sqrt * yx + alpha - 1; if (xx <= 0) continue; if (expon.ur(rng) <= (1 + yx * yx) * std::exp((alpha - 1) * std::log(xx / (alpha - 1)) - sqrt * yx)) { return beta * xx; } } } template EIGEN_STRONG_INLINE const Packet packetOp(Rng&& rng) { using namespace Eigen::internal; using RUtils = RandUtils; auto& cm = Rand::detail::CompressMask::get_inst(); RUtils ru; if (alpha < 1) { for (int _i = 0; ; ++_i) { EIGENRAND_CHECK_INFINITY_LOOP(); Packet ux = ru.uniform_real(rng); Packet vx = ru.nonzero_uniform_real(rng); Packet xx = pexp(pmul(pset1(1 / alpha), plog(vx))); Packet qx = pexp(pnegate(xx)); Packet xx2 = psub(pset1(1), plog(vx)); Packet qx2 = pexp(pmul(plog(xx2), pset1(alpha - 1))); auto c = pcmplt(ux, pset1(px)); xx = pblendv(c, xx, xx2); qx = pblendv(c, qx, qx2); ux = ru.uniform_real(rng); Packet cands = pmul(pset1(beta), xx); bool full = false; cache_rest_cnt = cm.compress_append(cands, pcmplt(ux, qx), OptCacheStore::template get(), cache_rest_cnt, full); if (full) return cands; } } if (alpha == 1) { return pmul(pset1(beta), expon.template packetOp(rng) ); } int count; if ((count = alpha) == alpha && count < 20) { RUtils ru; Packet ux, yx; yx = ru.nonzero_uniform_real(rng); while (--count) { yx = pmul(yx, ru.nonzero_uniform_real(rng)); } return pnegate(pmul(pset1(beta), plog(yx))); } else { for (int _i = 0; ; ++_i) { EIGENRAND_CHECK_INFINITY_LOOP(); Packet alpha_1 = pset1(alpha - 1); Packet ys, yc; psincos(pmul(pset1(constant::pi), ru.uniform_real(rng)), ys, yc); Packet yx = pdiv(ys, yc); Packet xx = padd(pmul(pset1(sqrt), yx), alpha_1); auto c = pcmplt(pset1(0), xx); Packet ux = ru.uniform_real(rng); Packet ub = pmul(padd(pmul(yx, yx), pset1(1)), pexp(psub( pmul(alpha_1, plog(pdiv(xx, alpha_1))), pmul(yx, pset1(sqrt)) )) ); c = pand(c, pcmple(ux, ub)); Packet cands = pmul(pset1(beta), xx); bool full = false; cache_rest_cnt = cm.compress_append(cands, c, OptCacheStore::template get(), cache_rest_cnt, full); if (full) return cands; } } } }; /** * @brief Generator of reals on a Weibull distribution * * @tparam _Scalar */ template class WeibullGen : public GenBase, _Scalar> { static_assert(std::is_floating_point<_Scalar>::value, "weilbullDist needs floating point types."); StdUniformRealGen<_Scalar> ur; _Scalar a = 1, b = 1; public: using Scalar = _Scalar; /** * @brief Construct a new Weibull generator * * @param _a shape parameter of the distribution * @param _b scale parameter of the distribution */ WeibullGen(_Scalar _a = 1, _Scalar _b = 1) : a{ _a }, b{ _b } { } WeibullGen(const WeibullGen&) = default; WeibullGen(WeibullGen&&) = default; WeibullGen& operator=(const WeibullGen&) = default; WeibullGen& operator=(WeibullGen&&) = default; template EIGEN_STRONG_INLINE const _Scalar operator() (Rng&& rng) { using namespace Eigen::internal; return std::pow(-std::log(1 - ur(std::forward(rng))), 1 / a) * b; } template EIGEN_STRONG_INLINE const Packet packetOp(Rng&& rng) { using namespace Eigen::internal; return pmul(pexp(pmul(plog(pnegate(plog( psub(pset1(1), ur.template packetOp(std::forward(rng))) ))), pset1(1 / a))), pset1(b)); } }; /** * @brief Generator of reals on an extreme value distribution * * @tparam _Scalar */ template class ExtremeValueGen : public GenBase, _Scalar> { static_assert(std::is_floating_point<_Scalar>::value, "extremeValueDist needs floating point types."); StdUniformRealGen<_Scalar> ur; _Scalar a = 0, b = 1; public: using Scalar = _Scalar; /** * @brief Construct a new extreme value generator * * @param _a location parameter of the distribution * @param _b scale parameter of the distribution */ ExtremeValueGen(_Scalar _a = 0, _Scalar _b = 1) : a{ _a }, b{ _b } { } ExtremeValueGen(const ExtremeValueGen&) = default; ExtremeValueGen(ExtremeValueGen&&) = default; ExtremeValueGen& operator=(const ExtremeValueGen&) = default; ExtremeValueGen& operator=(ExtremeValueGen&&) = default; template EIGEN_STRONG_INLINE const _Scalar operator() (Rng&& rng) { using namespace Eigen::internal; return (a - b * std::log(-std::log(ur.nzur_scalar(std::forward(rng))))); } template EIGEN_STRONG_INLINE const Packet packetOp(Rng&& rng) { using namespace Eigen::internal; using RUtils = RandUtils; return psub(pset1(a), pmul(plog(pnegate(plog(RUtils{}.nonzero_uniform_real(std::forward(rng))))), pset1(b)) ); } }; /** * @brief Generator of reals on a chi-squared distribution * * @tparam _Scalar */ template class ChiSquaredGen : public GenBase, _Scalar> { static_assert(std::is_floating_point<_Scalar>::value, "chiSquaredDist needs floating point types."); GammaGen<_Scalar> gamma; public: using Scalar = _Scalar; /** * @brief Construct a new chi-squared generator * * @param n degrees of freedom */ ChiSquaredGen(_Scalar n = 1) : gamma{ n * _Scalar(0.5), 2 } { } ChiSquaredGen(const ChiSquaredGen&) = default; ChiSquaredGen(ChiSquaredGen&&) = default; ChiSquaredGen& operator=(const ChiSquaredGen&) = default; ChiSquaredGen& operator=(ChiSquaredGen&&) = default; template EIGEN_STRONG_INLINE const _Scalar operator() (Rng&& rng) { return gamma(rng); } template EIGEN_STRONG_INLINE const Packet packetOp(Rng&& rng) { return gamma.template packetOp(rng); } }; /** * @brief Generator of reals on a Cauchy distribution * * @tparam _Scalar */ template class CauchyGen : public GenBase, _Scalar> { static_assert(std::is_floating_point<_Scalar>::value, "cauchyDist needs floating point types."); StdUniformRealGen<_Scalar> ur; _Scalar a = 0, b = 1; public: using Scalar = _Scalar; /** * @brief Construct a new Cauchy generator * * @param _a location parameter of the distribution * @param _b scale parameter of the distribution */ CauchyGen(_Scalar _a = 0, _Scalar _b = 1) : a{ _a }, b{ _b } { } CauchyGen(const CauchyGen&) = default; CauchyGen(CauchyGen&&) = default; CauchyGen& operator=(const CauchyGen&) = default; CauchyGen& operator=(CauchyGen&&) = default; template EIGEN_STRONG_INLINE const _Scalar operator() (Rng&& rng) { using namespace Eigen::internal; return a + b * std::tan(constant::pi * (ur(std::forward(rng)) - 0.5)); } template EIGEN_STRONG_INLINE const Packet packetOp(Rng&& rng) { using namespace Eigen::internal; Packet s, c; psincos(pmul(pset1(constant::pi), psub(ur.template packetOp(std::forward(rng)), pset1(0.5)) ), s, c); return padd(pset1(a), pmul(pset1(b), pdiv(s, c)) ); } }; template class FisherFGen; /** * @brief Generator of reals on a beta distribution * * @tparam _Scalar */ template class BetaGen : OptCacheStore, public GenBase, _Scalar> { friend FisherFGen<_Scalar>; static_assert(std::is_floating_point<_Scalar>::value, "betaDist needs floating point types."); int cache_rest_cnt = 0; StdUniformRealGen<_Scalar> ur; _Scalar a, b; GammaGen<_Scalar> gd1, gd2; public: using Scalar = _Scalar; /** * @brief Construct a new beta generator * * @param _a, _b shape parameter */ BetaGen(_Scalar _a = 1, _Scalar _b = 1) : a{ _a }, b{ _b }, gd1{ _a }, gd2{ _b } { } BetaGen(const BetaGen&) = default; BetaGen(BetaGen&&) = default; BetaGen& operator=(const BetaGen&) = default; BetaGen& operator=(BetaGen&&) = default; template EIGEN_STRONG_INLINE const _Scalar operator() (Rng&& rng) { using namespace Eigen::internal; if (a < 1 && b < 1) { _Scalar x, p1, p2; for (int _i = 0; ; ++_i) { EIGENRAND_CHECK_INFINITY_LOOP(); p1 = std::pow(ur(rng), 1 / a); p2 = std::pow(ur(rng), 1 / b); x = p1 + p2; if (x <= 1) break; } return p1 / x; } else { _Scalar p1 = gd1(rng), p2 = gd2(rng); return p1 / (p1 + p2); } } template EIGEN_STRONG_INLINE const Packet packetOp(Rng&& rng) { using namespace Eigen::internal; if (a < 1 && b < 1) { auto& cm = Rand::detail::CompressMask::get_inst(); Packet x, p1, p2; for (int _i = 0; ; ++_i) { EIGENRAND_CHECK_INFINITY_LOOP(); p1 = pexp(pmul(plog(ur.template packetOp(rng)), pset1(1 / a))); p2 = pexp(pmul(plog(ur.template packetOp(rng)), pset1(1 / b))); x = padd(p1, p2); Packet cands = pdiv(p1, x); bool full = false; cache_rest_cnt = cm.compress_append(cands, pcmple(x, pset1(1)), OptCacheStore::template get(), cache_rest_cnt, full); if (full) return cands; } } else { auto p1 = gd1.template packetOp(rng), p2 = gd2.template packetOp(rng); return pdiv(p1, padd(p1, p2)); } } }; /** * @brief Generator of reals on a Fisher's f distribution * * @tparam _Scalar */ template class FisherFGen : public GenBase, _Scalar> { static_assert(std::is_floating_point<_Scalar>::value, "fisherF needs floating point types."); BetaGen<_Scalar> beta; public: using Scalar = _Scalar; /** * @brief Construct a new Fisher's f generator * * @param m, n degrees of freedom */ FisherFGen(_Scalar m = 1, _Scalar n = 1) : beta{ m * _Scalar(0.5), n * _Scalar(0.5) } { } FisherFGen(const FisherFGen&) = default; FisherFGen(FisherFGen&&) = default; FisherFGen& operator=(const FisherFGen&) = default; FisherFGen& operator=(FisherFGen&&) = default; template EIGEN_STRONG_INLINE const _Scalar operator() (Rng&& rng) { using namespace Eigen::internal; auto x = beta(std::forward(rng)); return beta.b / beta.a * x / (1 - x); } template EIGEN_STRONG_INLINE const Packet packetOp(Rng&& rng) { using namespace Eigen::internal; auto x = beta.template packetOp(std::forward(rng)); return pdiv(pmul(pset1(beta.b / beta.a), x), psub(pset1(1), x)); } }; template using BetaType = CwiseNullaryOp, typename Derived::Scalar, Urng, true>, const Derived>; /** * @brief generates reals on the beta distribution. * * @tparam Derived * @tparam Urng * @param rows the number of rows being generated * @param cols the number of columns being generated * @param urng c++11-style random number generator * @param a,b shape parameter * @return a random matrix expression with a shape (`rows`, `cols`) * * @see Eigen::Rand::BetaGen */ template inline const BetaType beta(Index rows, Index cols, Urng&& urng, typename Derived::Scalar a = 1, typename Derived::Scalar b = 1) { return { rows, cols, { std::forward(urng), BetaGen{a, b} } }; } /** * @brief generates reals on the beta distribution. * * @tparam Derived * @tparam Urng * @param o an instance of any type of Eigen::DenseBase * @param urng c++11-style random number generator * @param a,b shape parameter * @return a random matrix expression of the same shape as `o` * * @see Eigen::Rand::BetaGen */ template inline const BetaType betaLike(Derived& o, Urng&& urng, typename Derived::Scalar a = 1, typename Derived::Scalar b = 1) { return { o.rows(), o.cols(), { std::forward(urng), BetaGen{a, b} } }; } template using CauchyType = CwiseNullaryOp, typename Derived::Scalar, Urng, true>, const Derived>; /** * @brief generates reals on the Cauchy distribution. * * @tparam Derived * @tparam Urng * @param rows the number of rows being generated * @param cols the number of columns being generated * @param urng c++11-style random number generator * @param a a location parameter of the distribution * @param b a scale parameter of the distribution * @return a random matrix expression with a shape (`rows`, `cols`) * * @see Eigen::Rand::CauchyGen */ template inline const CauchyType cauchy(Index rows, Index cols, Urng&& urng, typename Derived::Scalar a = 0, typename Derived::Scalar b = 1) { return { rows, cols, { std::forward(urng), CauchyGen{a, b} } }; } /** * @brief generates reals on the Cauchy distribution. * * @tparam Derived * @tparam Urng * @param o an instance of any type of Eigen::DenseBase * @param urng c++11-style random number generator * @param a a location parameter of the distribution * @param b a scale parameter of the distribution * @return a random matrix expression of the same shape as `o` * * @see Eigen::Rand::CauchyGen */ template inline const CauchyType cauchyLike(Derived& o, Urng&& urng, typename Derived::Scalar a = 0, typename Derived::Scalar b = 1) { return { o.rows(), o.cols(), { std::forward(urng), CauchyGen{a, b} } }; } template using NormalType = CwiseNullaryOp, typename Derived::Scalar, Urng, true>, const Derived>; /** * @brief generates reals on a standard normal distribution (`mean` = 0, `stdev`=1) * * @tparam Derived a type of Eigen::DenseBase * @tparam Urng * @param rows the number of rows being generated * @param cols the number of columns being generated * @param urng c++11-style random number generator * @return a random matrix expression with a shape (`rows`, `cols`) * * @see Eigen::Rand::StdNormalGen */ template inline const NormalType normal(Index rows, Index cols, Urng&& urng) { return { rows, cols, { std::forward(urng) } }; } /** * @brief generates reals on a standard normal distribution (`mean` = 0, `stdev`=1) * * @tparam Derived * @tparam Urng * @param o an instance of any type of Eigen::DenseBase * @param urng c++11-style random number generator * @return a random matrix expression of the same shape as `o` * * @see Eigen::Rand::StdNormalGen */ template inline const NormalType normalLike(Derived& o, Urng&& urng) { return { o.rows(), o.cols(), { std::forward(urng) } }; } template using Normal2Type = CwiseNullaryOp, typename Derived::Scalar, Urng, true>, const Derived>; /** * @brief generates reals on a normal distribution with arbitrary `mean` and `stdev`. * * @tparam Derived * @tparam Urng * @param rows the number of rows being generated * @param cols the number of columns being generated * @param urng c++11-style random number generator * @param mean a mean value of the distribution * @param stdev a standard deviation value of the distribution * @return a random matrix expression with a shape (`rows`, `cols`) * * @see Eigen::Rand::NormalGen */ template inline const Normal2Type normal(Index rows, Index cols, Urng&& urng, typename Derived::Scalar mean, typename Derived::Scalar stdev = 1) { return { rows, cols, { std::forward(urng), NormalGen{mean, stdev} } }; } /** * @brief generates reals on a normal distribution with arbitrary `mean` and `stdev`. * * @tparam Derived * @tparam Urng * @param o an instance of any type of Eigen::DenseBase * @param urng c++11-style random number generator * @param mean a mean value of the distribution * @param stdev a standard deviation value of the distribution * @return a random matrix expression of the same shape as `o` * * @see Eigen::Rand::NormalGen */ template inline const Normal2Type normalLike(Derived& o, Urng&& urng, typename Derived::Scalar mean, typename Derived::Scalar stdev = 1) { return { o.rows(), o.cols(), { std::forward(urng), NormalGen{mean, stdev} } }; } template using LognormalType = CwiseNullaryOp, typename Derived::Scalar, Urng, true>, const Derived>; /** * @brief generates reals on a lognormal distribution with arbitrary `mean` and `stdev`. * * @tparam Derived * @tparam Urng * @param rows the number of rows being generated * @param cols the number of columns being generated * @param urng c++11-style random number generator * @param mean a mean value of the distribution * @param stdev a standard deviation value of the distribution * @return a random matrix expression with a shape (`rows`, `cols`) * * @see Eigen::Rand::LognormalGen */ template inline const LognormalType lognormal(Index rows, Index cols, Urng&& urng, typename Derived::Scalar mean = 0, typename Derived::Scalar stdev = 1) { return { rows, cols, { std::forward(urng), LognormalGen{mean, stdev} } }; } /** * @brief generates reals on a lognormal distribution with arbitrary `mean` and `stdev`. * * @tparam Derived * @tparam Urng * @param o an instance of any type of Eigen::DenseBase * @param urng c++11-style random number generator * @param mean a mean value of the distribution * @param stdev a standard deviation value of the distribution * @return a random matrix expression of the same shape as `o` * * @see Eigen::Rand::LognormalGen */ template inline const LognormalType lognormalLike(Derived& o, Urng&& urng, typename Derived::Scalar mean = 0, typename Derived::Scalar stdev = 1) { return { o.rows(), o.cols(), { std::forward(urng), LognormalGen{mean, stdev} } }; } template using StudentTType = CwiseNullaryOp, typename Derived::Scalar, Urng, true>, const Derived>; /** * @brief generates reals on the Student's t distribution with arbirtrary degress of freedom. * * @tparam Derived a type of Eigen::DenseBase * @tparam Urng * @param rows the number of rows being generated * @param cols the number of columns being generated * @param urng c++11-style random number generator * @param n degrees of freedom * @return a random matrix expression with a shape (`rows`, `cols`) * * @see Eigen::Rand::StudentTGen */ template inline const StudentTType studentT(Index rows, Index cols, Urng&& urng, typename Derived::Scalar n = 1) { return { rows, cols, { std::forward(urng), StudentTGen{n} } }; } /** * @brief generates reals on the Student's t distribution with arbirtrary degress of freedom. * * @tparam Derived * @tparam Urng * @param o an instance of any type of Eigen::DenseBase * @param urng c++11-style random number generator * @param n degrees of freedom * @return a random matrix expression of the same shape as `o` * * @see Eigen::Rand::StudentTGen */ template inline const StudentTType studentTLike(Derived& o, Urng&& urng, typename Derived::Scalar n = 1) { return { o.rows(), o.cols(), { std::forward(urng), StudentTGen{n} } }; } template using ExponentialType = CwiseNullaryOp, typename Derived::Scalar, Urng, true>, const Derived>; /** * @brief generates reals on an exponential distribution with arbitrary scale parameter. * * @tparam Derived * @tparam Urng * @param rows the number of rows being generated * @param cols the number of columns being generated * @param urng c++11-style random number generator * @param lambda a scale parameter of the distribution * @return a random matrix expression with a shape (`rows`, `cols`) * * @see Eigen::Rand::ExponentialGen */ template inline const ExponentialType exponential(Index rows, Index cols, Urng&& urng, typename Derived::Scalar lambda = 1) { return { rows, cols, { std::forward(urng), ExponentialGen{lambda} } }; } /** * @brief generates reals on an exponential distribution with arbitrary scale parameter. * * @tparam Derived * @tparam Urng * @param o an instance of any type of Eigen::DenseBase * @param urng c++11-style random number generator * @param lambda a scale parameter of the distribution * @return a random matrix expression of the same shape as `o` * * @see Eigen::Rand::ExponentialGen */ template inline const ExponentialType exponentialLike(Derived& o, Urng&& urng, typename Derived::Scalar lambda = 1) { return { o.rows(), o.cols(), { std::forward(urng), ExponentialGen{lambda} } }; } template using GammaType = CwiseNullaryOp, typename Derived::Scalar, Urng, true>, const Derived>; /** * @brief generates reals on a gamma distribution with arbitrary shape and scale parameter. * * @tparam Derived * @tparam Urng * @param rows the number of rows being generated * @param cols the number of columns being generated * @param urng c++11-style random number generator * @param alpha a shape parameter of the distribution * @param beta a scale parameter of the distribution * @return a random matrix expression with a shape (`rows`, `cols`) * * @see Eigen::Rand::GammaGen */ template inline const GammaType gamma(Index rows, Index cols, Urng&& urng, typename Derived::Scalar alpha = 1, typename Derived::Scalar beta = 1) { return { rows, cols, { std::forward(urng), GammaGen{alpha, beta} } }; } /** * @brief generates reals on a gamma distribution with arbitrary shape and scale parameter. * * @tparam Derived * @tparam Urng * @param o an instance of any type of Eigen::DenseBase * @param urng c++11-style random number generator * @param alpha a shape parameter of the distribution * @param beta a scale parameter of the distribution * @return a random matrix expression of the same shape as `o` * * @see Eigen::Rand::GammaGen */ template inline const GammaType gammaLike(Derived& o, Urng&& urng, typename Derived::Scalar alpha = 1, typename Derived::Scalar beta = 1) { return { o.rows(), o.cols(), { std::forward(urng), GammaGen{alpha, beta} } }; } template using WeibullType = CwiseNullaryOp, typename Derived::Scalar, Urng, true>, const Derived>; /** * @brief generates reals on a Weibull distribution with arbitrary shape and scale parameter. * * @tparam Derived * @tparam Urng * @param rows the number of rows being generated * @param cols the number of columns being generated * @param urng c++11-style random number generator * @param a a shape parameter of the distribution * @param b a scale parameter of the distribution * @return a random matrix expression with a shape (`rows`, `cols`) * * @see Eigen::Rand::WeibullGen */ template inline const WeibullType weibull(Index rows, Index cols, Urng&& urng, typename Derived::Scalar a = 1, typename Derived::Scalar b = 1) { return { rows, cols, { std::forward(urng), WeibullGen{a, b} } }; } /** * @brief generates reals on a Weibull distribution with arbitrary shape and scale parameter. * * @tparam Derived * @tparam Urng * @param o an instance of any type of Eigen::DenseBase * @param urng c++11-style random number generator * @param a a shape parameter of the distribution * @param b a scale parameter of the distribution * @return a random matrix expression of the same shape as `o` * * @see Eigen::Rand::WeibullGen */ template inline const WeibullType weibullLike(Derived& o, Urng&& urng, typename Derived::Scalar a = 1, typename Derived::Scalar b = 1) { return { o.rows(), o.cols(), { std::forward(urng), WeibullGen{a, b} } }; } template using ExtremeValueType = CwiseNullaryOp, typename Derived::Scalar, Urng, true>, const Derived>; /** * @brief generates reals on an extreme value distribution * (a.k.a Gumbel Type I, log-Weibull, Fisher-Tippett Type I) with arbitrary shape and scale parameter. * * @tparam Derived * @tparam Urng * @param rows the number of rows being generated * @param cols the number of columns being generated * @param urng c++11-style random number generator * @param a a location parameter of the distribution * @param b a scale parameter of the distribution * @return a random matrix expression with a shape (`rows`, `cols`) * * @see Eigen::Rand::ExtremeValueGen */ template inline const ExtremeValueType extremeValue(Index rows, Index cols, Urng&& urng, typename Derived::Scalar a = 0, typename Derived::Scalar b = 1) { return { rows, cols, { std::forward(urng), ExtremeValueGen{a, b} } }; } /** * @brief generates reals on an extreme value distribution * (a.k.a Gumbel Type I, log-Weibull, Fisher-Tippett Type I) with arbitrary shape and scale parameter. * * @tparam Derived * @tparam Urng * @param o an instance of any type of Eigen::DenseBase * @param urng c++11-style random number generator * @param a a location parameter of the distribution * @param b a scale parameter of the distribution * @return a random matrix expression of the same shape as `o` * * @see Eigen::Rand::ExtremeValueGen */ template inline const ExtremeValueType extremeValueLike(Derived& o, Urng&& urng, typename Derived::Scalar a = 0, typename Derived::Scalar b = 1) { return { o.rows(), o.cols(), { std::forward(urng), ExtremeValueGen{a, b} } }; } template using ChiSquaredType = CwiseNullaryOp, typename Derived::Scalar, Urng, true>, const Derived>; /** * @brief generates reals on the Chi-squared distribution with arbitrary degrees of freedom. * * @tparam Derived * @tparam Urng * @param rows the number of rows being generated * @param cols the number of columns being generated * @param urng c++11-style random number generator * @param n the degrees of freedom of the distribution * @return a random matrix expression with a shape (`rows`, `cols`) * * @see Eigen::Rand::ChiSquaredGen */ template inline const ChiSquaredType chiSquared(Index rows, Index cols, Urng&& urng, typename Derived::Scalar n = 1) { return { rows, cols, { std::forward(urng), ChiSquaredGen{n} } }; } /** * @brief generates reals on the Chi-squared distribution with arbitrary degrees of freedom. * * @tparam Derived * @tparam Urng * @param o an instance of any type of Eigen::DenseBase * @param urng c++11-style random number generator * @param n the degrees of freedom of the distribution * @return a random matrix expression of the same shape as `o` * * @see Eigen::Rand::ChiSquaredGen */ template inline const ChiSquaredType chiSquaredLike(Derived& o, Urng&& urng, typename Derived::Scalar n = 1) { return { o.rows(), o.cols(), { std::forward(urng), ChiSquaredGen{n} } }; } template using FisherFType = CwiseNullaryOp, typename Derived::Scalar, Urng, true>, const Derived>; /** * @brief generates reals on the Fisher's F distribution. * * @tparam Derived * @tparam Urng * @param rows the number of rows being generated * @param cols the number of columns being generated * @param urng c++11-style random number generator * @param m degrees of freedom * @param n degrees of freedom * @return a random matrix expression with a shape (`rows`, `cols`) * * @see Eigen::Rand::FisherFGen */ template inline const FisherFType fisherF(Index rows, Index cols, Urng&& urng, typename Derived::Scalar m = 1, typename Derived::Scalar n = 1) { return { rows, cols, { std::forward(urng), FisherFGen{m, n} } }; } /** * @brief generates reals on the Fisher's F distribution. * * @tparam Derived * @tparam Urng * @param o an instance of any type of Eigen::DenseBase * @param urng c++11-style random number generator * @param m degrees of freedom * @param n degrees of freedom * @return a random matrix expression of the same shape as `o` * * @see Eigen::Rand::FisherFGen */ template inline const FisherFType fisherFLike(Derived& o, Urng&& urng, typename Derived::Scalar m = 1, typename Derived::Scalar n = 1) { return { o.rows(), o.cols(), { std::forward(urng), FisherFGen{m, n} } }; } } } #endif