random/distribution/ExponentialDistribution.hh

0001 //------------------------------- -*- C++ -*- -------------------------------//
0002 // Copyright Celeritas contributors: see top-level COPYRIGHT file for details
0003 // SPDX-License-Identifier: (Apache-2.0 OR MIT)
0004 //---------------------------------------------------------------------------//
0005 //! \file corecel/random/distribution/ExponentialDistribution.hh
0006 //---------------------------------------------------------------------------//
0007 #pragma once
0008
0009 #include <cmath>
0010
0011 #include "corecel/Assert.hh"
0012 #include "corecel/Macros.hh"
0013 #include "corecel/Types.hh"
0014
0015 #include "GenerateCanonical.hh"
0016
0017 namespace celeritas
0018 {
0019 //---------------------------------------------------------------------------//
0020 /*!
0021  * Sample from an exponential distribution.
0022  *
0023  * Sample from a probability distribution function with the normalized PDF:
0024  * \f[
0025    f(x; \lambda) = \lambda e^{-\lambda x} \quad \mathrm{for} \  x >= 0
0026    \f]
0027  * which integrated into a CDF and inverted gives a sample:
0028  * \f[
0029   x = -\frac{\log \xi}{ \lambda}
0030    \f]
0031  *
0032  * This is simply an implementation of std::exponential_distribution that uses
0033  * the Celeritas canonical generator and is independent of library
0034  * implementation.
0035  *
0036  * Note (for performance-critical sections of code) that if this class is
0037  * constructed locally with the default value of lambda = 1.0, the
0038  * inversion and multiplication will be optimized out.
0039  */
0040 template<class RealType = ::celeritas::real_type>
0041 class ExponentialDistribution
0042 {
0043   public:
0044     //!@{
0045     //! \name Type aliases
0046     using real_type = RealType;
0047     using result_type = RealType;
0048     //!@}
0049
0050   public:
0051     // Construct with defaults
0052     explicit inline CELER_FUNCTION
0053     ExponentialDistribution(real_type lambda = 1);
0054
0055     // Sample using the given random number generator
0056     template<class Generator>
0057     inline CELER_FUNCTION result_type operator()(Generator& g);
0058
0059   private:
0060     result_type neg_inv_lambda_;
0061 };
0062
0063 //---------------------------------------------------------------------------//
0064 // INLINE DEFINITIONS
0065 //---------------------------------------------------------------------------//
0066 /*!
0067  * Construct from the mean of the exponential distribution.
0068  */
0069 template<class RT>
0070 CELER_FUNCTION
0071 ExponentialDistribution<RT>::ExponentialDistribution(real_type lambda)
0072     : neg_inv_lambda_(real_type{-1} / lambda)
0073 {
0074     CELER_EXPECT(lambda > real_type{0});
0075 }
0076
0077 //---------------------------------------------------------------------------//
0078 /*!
0079  * Sample a random number according to the distribution.
0080  */
0081 template<class RT>
0082 template<class Generator>
0083 CELER_FUNCTION auto
0084 ExponentialDistribution<RT>::operator()(Generator& rng) -> result_type
0085 {
0086     return std::log(generate_canonical<RT>(rng)) * neg_inv_lambda_;
0087 }
0088
0089 //---------------------------------------------------------------------------//
0090 }  // namespace celeritas