[flang-commits] [flang] e111404 - [flang] RANDOM_NUMBER, RANDOM_SEED, RANDOM_INIT in runtime
peter klausler via flang-commits
flang-commits at lists.llvm.org
Fri Apr 9 17:03:21 PDT 2021
Author: peter klausler
Date: 2021-04-09T17:03:11-07:00
New Revision: e11140451b5f8b4d0a34297ba8aac77ded4dbc2a
URL: https://github.com/llvm/llvm-project/commit/e11140451b5f8b4d0a34297ba8aac77ded4dbc2a
DIFF: https://github.com/llvm/llvm-project/commit/e11140451b5f8b4d0a34297ba8aac77ded4dbc2a.diff
LOG: [flang] RANDOM_NUMBER, RANDOM_SEED, RANDOM_INIT in runtime
Add APIs, initial non-coarray implementations, and unit
tests for the intrinsic subroutines for pseudo-random
number generation.
Differential Revision: https://reviews.llvm.org/D100064
Added:
flang/runtime/random.cpp
flang/runtime/random.h
flang/unittests/RuntimeGTest/Random.cpp
Modified:
flang/runtime/CMakeLists.txt
flang/unittests/RuntimeGTest/CMakeLists.txt
Removed:
################################################################################
diff --git a/flang/runtime/CMakeLists.txt b/flang/runtime/CMakeLists.txt
index 6f95ff89d0411..c391cb2c9bd39 100644
--- a/flang/runtime/CMakeLists.txt
+++ b/flang/runtime/CMakeLists.txt
@@ -53,6 +53,7 @@ add_flang_library(FortranRuntime
memory.cpp
misc-intrinsic.cpp
numeric.cpp
+ random.cpp
reduction.cpp
stat.cpp
stop.cpp
diff --git a/flang/runtime/random.cpp b/flang/runtime/random.cpp
new file mode 100644
index 0000000000000..b5158d07c816c
--- /dev/null
+++ b/flang/runtime/random.cpp
@@ -0,0 +1,193 @@
+//===-- runtime/random.cpp ------------------------------------------------===//
+//
+// Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions.
+// See https://llvm.org/LICENSE.txt for license information.
+// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
+//
+//===----------------------------------------------------------------------===//
+
+// Implements the intrinsic subroutines RANDOM_INIT, RANDOM_NUMBER, and
+// RANDOM_SEED.
+
+#include "random.h"
+#include "cpp-type.h"
+#include "descriptor.h"
+#include "lock.h"
+#include "flang/Common/leading-zero-bit-count.h"
+#include "flang/Common/uint128.h"
+#include <algorithm>
+#include <cmath>
+#include <cstdint>
+#include <ctime>
+#include <limits>
+#include <memory>
+#include <random>
+
+namespace Fortran::runtime {
+
+// Newer "Minimum standard", recommended by Park, Miller, and Stockmeyer in
+// 1993. Same as C++17 std::minstd_rand, but explicitly instantiated for
+// permanence.
+using Generator =
+ std::linear_congruential_engine<std::uint_fast32_t, 48271, 0, 2147483647>;
+
+using GeneratedWord = typename Generator::result_type;
+static constexpr std::uint64_t range{
+ static_cast<std::uint64_t>(Generator::max() - Generator::min() + 1)};
+static constexpr bool rangeIsPowerOfTwo{(range & (range - 1)) == 0};
+static constexpr int rangeBits{
+ 64 - common::LeadingZeroBitCount(range) - !rangeIsPowerOfTwo};
+
+static Lock lock;
+static Generator generator;
+
+template <typename REAL, int PREC>
+inline void Generate(const Descriptor &harvest) {
+ static constexpr std::size_t minBits{
+ std::max<std::size_t>(PREC, 8 * sizeof(GeneratedWord))};
+ using Int = common::HostUnsignedIntType<minBits>;
+ static constexpr std::size_t words{
+ static_cast<std::size_t>(PREC + rangeBits - 1) / rangeBits};
+ std::size_t elements{harvest.Elements()};
+ SubscriptValue at[maxRank];
+ harvest.GetLowerBounds(at);
+ {
+ CriticalSection critical{lock};
+ for (std::size_t j{0}; j < elements; ++j) {
+ Int fraction{generator()};
+ if constexpr (words > 1) {
+ for (std::size_t k{1}; k < words; ++k) {
+ static constexpr Int rangeMask{(Int{1} << rangeBits) - 1};
+ GeneratedWord word{(generator() - generator.min()) & rangeMask};
+ fraction = (fraction << rangeBits) | word;
+ }
+ }
+ fraction >>= words * rangeBits - PREC;
+ *harvest.Element<REAL>(at) =
+ std::ldexp(static_cast<REAL>(fraction), -(PREC + 1));
+ harvest.IncrementSubscripts(at);
+ }
+ }
+}
+
+extern "C" {
+
+void RTNAME(RandomInit)(bool repeatable, bool /*image_distinct*/) {
+ // TODO: multiple images and image_distinct: add image number
+ {
+ CriticalSection critical{lock};
+ if (repeatable) {
+ generator.seed(0);
+ } else {
+ generator.seed(std::time(nullptr));
+ }
+ }
+}
+
+void RTNAME(RandomNumber)(
+ const Descriptor &harvest, const char *source, int line) {
+ Terminator terminator{source, line};
+ auto typeCode{harvest.type().GetCategoryAndKind()};
+ RUNTIME_CHECK(terminator, typeCode && typeCode->first == TypeCategory::Real);
+ int kind{typeCode->second};
+ switch (kind) {
+ // TODO: REAL (2 & 3)
+ case 4:
+ Generate<CppTypeFor<TypeCategory::Real, 4>, 24>(harvest);
+ break;
+ case 8:
+ Generate<CppTypeFor<TypeCategory::Real, 8>, 53>(harvest);
+ break;
+#if LONG_DOUBLE == 80
+ case 10:
+ Generate<CppTypeFor<TypeCategory::Real, 10>, 64>(harvest);
+ break;
+#elif LONG_DOUBLE == 128
+ case 4:
+ Generate<CppTypeFor<TypeCategory::Real, 16>, 113>(harvest);
+ break;
+#endif
+ default:
+ terminator.Crash("RANDOM_NUMBER(): unsupported REAL kind %d", kind);
+ }
+}
+
+void RTNAME(RandomSeedSize)(
+ const Descriptor &size, const char *source, int line) {
+ Terminator terminator{source, line};
+ auto typeCode{size.type().GetCategoryAndKind()};
+ RUNTIME_CHECK(terminator,
+ size.rank() == 0 && typeCode && typeCode->first == TypeCategory::Integer);
+ int kind{typeCode->second};
+ switch (kind) {
+ case 4:
+ *size.OffsetElement<CppTypeFor<TypeCategory::Integer, 4>>() = 1;
+ break;
+ case 8:
+ *size.OffsetElement<CppTypeFor<TypeCategory::Integer, 8>>() = 1;
+ break;
+ default:
+ terminator.Crash("RANDOM_SEED(SIZE=): bad kind %d\n", kind);
+ }
+}
+
+void RTNAME(RandomSeedPut)(
+ const Descriptor &put, const char *source, int line) {
+ Terminator terminator{source, line};
+ auto typeCode{put.type().GetCategoryAndKind()};
+ RUNTIME_CHECK(terminator,
+ put.rank() == 1 && typeCode && typeCode->first == TypeCategory::Integer &&
+ put.GetDimension(0).Extent() >= 1);
+ int kind{typeCode->second};
+ GeneratedWord seed;
+ switch (kind) {
+ case 4:
+ seed = *put.OffsetElement<CppTypeFor<TypeCategory::Integer, 4>>();
+ break;
+ case 8:
+ seed = *put.OffsetElement<CppTypeFor<TypeCategory::Integer, 8>>();
+ break;
+ default:
+ terminator.Crash("RANDOM_SEED(PUT=): bad kind %d\n", kind);
+ }
+ {
+ CriticalSection critical{lock};
+ generator.seed(seed);
+ }
+}
+
+void RTNAME(RandomSeedDefaultPut)() {
+ // TODO: should this be time &/or image dependent?
+ {
+ CriticalSection critical{lock};
+ generator.seed(0);
+ }
+}
+
+void RTNAME(RandomSeedGet)(
+ const Descriptor &got, const char *source, int line) {
+ Terminator terminator{source, line};
+ auto typeCode{got.type().GetCategoryAndKind()};
+ RUNTIME_CHECK(terminator,
+ got.rank() == 1 && typeCode && typeCode->first == TypeCategory::Integer &&
+ got.GetDimension(0).Extent() >= 1);
+ int kind{typeCode->second};
+ GeneratedWord seed;
+ {
+ CriticalSection critical{lock};
+ seed = generator();
+ generator.seed(seed);
+ }
+ switch (kind) {
+ case 4:
+ *got.OffsetElement<CppTypeFor<TypeCategory::Integer, 4>>() = seed;
+ break;
+ case 8:
+ *got.OffsetElement<CppTypeFor<TypeCategory::Integer, 8>>() = seed;
+ break;
+ default:
+ terminator.Crash("RANDOM_SEED(GET=): bad kind %d\n", kind);
+ }
+}
+} // extern "C"
+} // namespace Fortran::runtime
diff --git a/flang/runtime/random.h b/flang/runtime/random.h
new file mode 100644
index 0000000000000..5a8488f8059be
--- /dev/null
+++ b/flang/runtime/random.h
@@ -0,0 +1,30 @@
+//===-- runtime/random.h --------------------------------------------------===//
+//
+// Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions.
+// See https://llvm.org/LICENSE.txt for license information.
+// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
+//
+//===----------------------------------------------------------------------===//
+
+// Intrinsic subroutines RANDOM_INIT, RANDOM_NUMBER, and RANDOM_SEED.
+
+#include "entry-names.h"
+#include <cstdint>
+
+namespace Fortran::runtime {
+class Descriptor;
+extern "C" {
+
+void RTNAME(RandomInit)(bool repeatable, bool image_distinct);
+
+void RTNAME(RandomNumber)(
+ const Descriptor &harvest, const char *source, int line);
+
+// Subroutine RANDOM_SEED can be called with at most one of its optional
+// arguments; they each (plus the default case) map to these entry points.
+void RTNAME(RandomSeedSize)(const Descriptor &, const char *source, int line);
+void RTNAME(RandomSeedPut)(const Descriptor &, const char *source, int line);
+void RTNAME(RandomSeedGet)(const Descriptor &, const char *source, int line);
+void RTNAME(RandomSeedDefaultPut)();
+} // extern "C"
+} // namespace Fortran::runtime
diff --git a/flang/unittests/RuntimeGTest/CMakeLists.txt b/flang/unittests/RuntimeGTest/CMakeLists.txt
index 38f08d7912bad..80314c124831b 100644
--- a/flang/unittests/RuntimeGTest/CMakeLists.txt
+++ b/flang/unittests/RuntimeGTest/CMakeLists.txt
@@ -4,6 +4,7 @@ add_flang_unittest(FlangRuntimeTests
MiscIntrinsic.cpp
Numeric.cpp
NumericalFormatTest.cpp
+ Random.cpp
Reduction.cpp
RuntimeCrashTest.cpp
)
diff --git a/flang/unittests/RuntimeGTest/Random.cpp b/flang/unittests/RuntimeGTest/Random.cpp
new file mode 100644
index 0000000000000..f16f3dc1a47fa
--- /dev/null
+++ b/flang/unittests/RuntimeGTest/Random.cpp
@@ -0,0 +1,63 @@
+//===-- flang/unittests/RuntimeGTest/Random.cpp -----------------*- C++ -*-===//
+//
+// Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions.
+// See https://llvm.org/LICENSE.txt for license information.
+// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
+//
+//===----------------------------------------------------------------------===//
+
+#include "../../runtime/random.h"
+#include "gtest/gtest.h"
+#include "../../runtime/descriptor.h"
+#include "../../runtime/type-code.h"
+#include <cmath>
+
+using namespace Fortran::runtime;
+
+TEST(RandomNumber, Real4) {
+ StaticDescriptor<1> statDesc;
+ Descriptor &harvest{statDesc.descriptor()};
+ static constexpr int n{10000};
+ float xs[n]{};
+ SubscriptValue extent[1]{n};
+ harvest.Establish(TypeCategory::Real, 4, xs, 1, extent);
+ RTNAME(RandomNumber)(harvest, __FILE__, __LINE__);
+ double sum{0};
+ for (int j{0}; j < n; ++j) {
+ sum += xs[j];
+ }
+ double mean{sum / n};
+ std::fprintf(stderr, "mean of %d random numbers: %g\n", n, mean);
+ EXPECT_GE(mean, 0.95 * 0.5); // mean of uniform dist [0..1] is of course 0.5
+ EXPECT_LE(mean, 1.05 * 0.5);
+ double sumsq{0};
+ for (int j{0}; j < n; ++j) {
+ double
diff {xs[j] - mean};
+ sumsq +=
diff *
diff ;
+ }
+ double sdev{std::sqrt(sumsq / n)};
+ std::fprintf(stderr, "stddev of %d random numbers: %g\n", n, sdev);
+ double expect{1.0 / std::sqrt(12.0)}; // stddev of uniform dist [0..1]
+ EXPECT_GE(sdev, 0.95 * expect);
+ EXPECT_LT(sdev, 1.05 * expect);
+}
+
+TEST(RandomNumber, RandomSeed) {
+ StaticDescriptor<1> statDesc[2];
+ Descriptor &desc{statDesc[0].descriptor()};
+ std::int32_t n;
+ desc.Establish(TypeCategory::Integer, 4, &n, 0, nullptr);
+ RTNAME(RandomSeedSize)(desc, __FILE__, __LINE__);
+ EXPECT_EQ(n, 1);
+ SubscriptValue extent[1]{1};
+ desc.Establish(TypeCategory::Integer, 4, &n, 1, extent);
+ RTNAME(RandomSeedGet)(desc, __FILE__, __LINE__);
+ Descriptor &harvest{statDesc[1].descriptor()};
+ float x;
+ harvest.Establish(TypeCategory::Real, 4, &x, 1, extent);
+ RTNAME(RandomNumber)(harvest, __FILE__, __LINE__);
+ float got{x};
+ RTNAME(RandomSeedPut)(desc, __FILE__, __LINE__); // n from RandomSeedGet()
+ RTNAME(RandomNumber)(harvest, __FILE__, __LINE__);
+ EXPECT_EQ(x, got);
+}
More information about the flang-commits
mailing list