star-sp

commit 5ec1c0d0ad2dfd8913aa376663a7d5a6e13ebbaf
parent cbb39bbd4dd31a4f240b555a200fcad8a7bfe790
Author: Vincent Forest <vincent.forest@meso-star.com>
Date:   Tue, 29 Sep 2015 12:09:29 +0200

Add the gaussian and lognormal distributions

Diffstat:
M
cmake/CMakeLists.txt
 | 
4
++--
M
src/ssp.h
 | 
27
+++++++++++++++++++++++++++
M
src/ssp_rng.c
 | 
40
++++++++++++++++++++++++++++++++++++++++

3 files changed, 69 insertions(+), 2 deletions(-)
diff --git a/cmake/CMakeLists.txt b/cmake/CMakeLists.txt
@@ -66,8 +66,8 @@ rcmake_append_runtime_dirs(_runtime_dirs RSys ${Boost_LIBRARY_DIRS})
 # Configure and define targets
 ################################################################################
 set(VERSION_MAJOR 0)
-set(VERSION_MINOR 2)
-set(VERSION_PATCH 1)
+set(VERSION_MINOR 3)
+set(VERSION_PATCH 0)
 set(VERSION ${VERSION_MAJOR}.${VERSION_MINOR}.${VERSION_PATCH})
 
 set(SSP_FILES_SRC ssp_rng.c ssp_rng_proxy.c)
diff --git a/src/ssp.h b/src/ssp.h
@@ -220,6 +220,33 @@ ssp_ran_exp_pdf
   (const double x,
    const double mu);
 
+/* Random variate from the guassian (or normal) distribution with mean `mu':
+ * P(x) dx = 1 / (sigma*sqrt(2*PI)) * exp(1/2*((x-mu)/sigma)^2) dx */
+SSP_API double
+ssp_ran_gaussian
+  (struct ssp_rng* rng,
+   const double mu,
+   const double sigma);
+
+SSP_API double
+ssp_ran_gaussian_pdf
+  (const double mu,
+   const double sigma);
+
+/* Random variate from the lognormal distribution:
+ * P(x) dx = 1/(sigma*x*sqrt(2*PI)) * exp(-(ln(x)-zeta)^2 / (2*sigma^2)) dx */
+SSP_API double
+ssp_ran_lognormal
+  (struct ssp_rng* rng,
+   const double zeta,
+   const double sigma);
+
+SSP_API double
+ssp_ran_lognormal_pdf
+  (const double x,
+   const double zeta,
+   const double sigma);
+
 /*******************************************************************************
  * Sphere distribution
  ******************************************************************************/
diff --git a/src/ssp_rng.c b/src/ssp_rng.c
@@ -74,6 +74,8 @@
 #include <cstring>
 #include <limits>
 
+#define SQRT_2_PI 2.50662827463100050240
+
 /*******************************************************************************
  * KISS PRNG
  ******************************************************************************/
@@ -695,6 +697,44 @@ ssp_ran_exp_pdf(const double x, const double mu)
   return mu * exp(-x * mu);
 }
 
+double
+ssp_ran_gaussian
+  (struct ssp_rng* rng,
+   const double mu,
+   const double sigma)
+{
+  ASSERT(rng);
+  rng_cxx rng_cxx(*rng);
+  RAN_NAMESPACE::normal_distribution<double> distribution(mu, sigma);
+  return distribution(rng_cxx);
+}
+
+double
+ssp_ran_gaussian_pdf(const double x, const double mu, const double sigma)
+{
+  const double tmp = (x - mu) / sigma;
+  return 1.0 / (sigma * SQRT_2_PI) * exp(-0.5 * tmp * tmp);
+}
+
+double
+ssp_ran_lognormal
+  (struct ssp_rng* rng,
+   const double zeta,
+   const double sigma)
+{
+  ASSERT(rng);
+  rng_cxx rng_cxx(*rng);
+  RAN_NAMESPACE::normal_distribution<double> distribution(zeta, sigma);
+  return distribution(rng_cxx);
+}
+
+double
+ssp_ran_lognormal_pdf(const double x, const double zeta, const double sigma)
+{
+  const double tmp = log(x) - zeta;
+  return 1.0 / (sigma * x * SQRT_2_PI) * exp(-(tmp*tmp) / (2*sigma*sigma));
+}
+
 #ifdef COMPILER_CL
   #pragma warning(pop)
 #endif