From 7bdc076a460a5085c978c45a6a5d4cf76be093b2 Mon Sep 17 00:00:00 2001
From: Maximilian Reininghaus <maximilian.reininghaus@kit.edu>
Date: Mon, 17 Dec 2018 17:54:04 +0100
Subject: [PATCH] uniform real distribution with units
---
Framework/Random/CMakeLists.txt | 1 +
Framework/Random/UniformRealDistribution.h | 32 +++++++++++++++
Framework/Random/testRandom.cc | 48 +++++++++++++++++++++-
3 files changed, 80 insertions(+), 1 deletion(-)
create mode 100644 Framework/Random/UniformRealDistribution.h
diff --git a/Framework/Random/CMakeLists.txt b/Framework/Random/CMakeLists.txt
index 3c1947d2..230f84f6 100644
--- a/Framework/Random/CMakeLists.txt
+++ b/Framework/Random/CMakeLists.txt
@@ -7,6 +7,7 @@ set (
set (
CORSIKArandom_HEADERS
RNGManager.h
+ UniformRealDistribution.h
)
set (
diff --git a/Framework/Random/UniformRealDistribution.h b/Framework/Random/UniformRealDistribution.h
new file mode 100644
index 00000000..c90ab674
--- /dev/null
+++ b/Framework/Random/UniformRealDistribution.h
@@ -0,0 +1,32 @@
+#ifndef _include_random_distributions_h
+#define _include_random_distributions_h
+
+#include <corsika/units/PhysicalUnits.h>
+#include <random>
+
+namespace corsika::random {
+
+ template <class TQuantity>
+ class UniformRealDistribution {
+ using RealType = typename TQuantity::value_type;
+ std::uniform_real_distribution<RealType> dist{RealType(0.), RealType(1.)};
+
+ TQuantity const a, b;
+
+ public:
+ UniformRealDistribution(TQuantity b)
+ : a{TQuantity(phys::units::detail::magnitude_tag, 0)}
+ , b(b) {}
+ UniformRealDistribution(TQuantity a, TQuantity b)
+ : a(a)
+ , b(b) {}
+
+ template <class Generator>
+ TQuantity operator()(Generator& g) {
+ return a + dist(g) * (b - a);
+ }
+ };
+
+} // namespace corsika::random
+
+#endif
diff --git a/Framework/Random/testRandom.cc b/Framework/Random/testRandom.cc
index 041a96c2..56d8157a 100644
--- a/Framework/Random/testRandom.cc
+++ b/Framework/Random/testRandom.cc
@@ -1,4 +1,3 @@
-
/**
* (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
*
@@ -14,7 +13,11 @@
#include <catch2/catch.hpp>
#include <corsika/random/RNGManager.h>
+#include <corsika/random/UniformRealDistribution.h>
+#include <corsika/units/PhysicalUnits.h>
#include <iostream>
+#include <limits>
+#include <random>
using namespace corsika::random;
@@ -37,3 +40,46 @@ SCENARIO("random-number streams can be registered and retrieved") {
}
}
}
+
+TEST_CASE("UniformRealDistribution") {
+ using namespace corsika::units::si;
+ std::mt19937 rng;
+
+ corsika::random::UniformRealDistribution<LengthType> dist(1_m, 2_m);
+
+ SECTION("range") {
+ corsika::random::UniformRealDistribution<LengthType> dist(1_m, 2_m);
+
+ LengthType min =
+ +1_m * std::numeric_limits<typename LengthType::value_type>::infinity();
+ LengthType max =
+ -1_m * std::numeric_limits<typename LengthType::value_type>::infinity();
+
+ for (int i{0}; i < 1'000'000; ++i) {
+ LengthType x = dist(rng);
+ min = std::min(min, x);
+ max = std::max(max, x);
+ }
+
+ CHECK(min / 1_m == Approx(1.));
+ CHECK(max / 2_m == Approx(1.));
+ }
+
+ SECTION("range") {
+ corsika::random::UniformRealDistribution<LengthType> dist(18_cm);
+
+ LengthType min =
+ +1_m * std::numeric_limits<typename LengthType::value_type>::infinity();
+ LengthType max =
+ -1_m * std::numeric_limits<typename LengthType::value_type>::infinity();
+
+ for (int i{0}; i < 1'000'000; ++i) {
+ LengthType x = dist(rng);
+ min = std::min(min, x);
+ max = std::max(max, x);
+ }
+
+ CHECK(min / 1_m == Approx(0.).margin(1e-3));
+ CHECK(max / 18_cm == Approx(1.));
+ }
+}
--
GitLab