Re-enable Sobol64 distributions (#290)

CHANGELOG.md

@@ -9,10 +9,8 @@ Full documentation for rocRAND is available at [https://rocrand.readthedocs.io/e
 ### Changed
 - The `mrg_<distribution>_distribution` structures, which provided numbers based on MRG32K3A, are now replaced by `mrg_engine_<distribution>_distribution`, where `<distribution>` is `log_normal`, `normal`, `poisson`, or `uniform`. These structures provide numbers for MRG31K3P (with template type `rocrand_state_mrg31k3p`) and MRG32K3A (with template type `rocrand_state_mrg32k3a`).
 ### Fixed
-- Sobol64 now returns 64 bit instead of 32 bit random numbers.
+- Sobol64 now returns 64 bits random numbers, instead of 32 bits random numbers. As a result, the performance of this generator has regressed.
 - Fixed a bug that prevented compiling code in C++ mode (with a host compiler) when it included the rocRAND headers on Windows.
-- ### Known issues
-- uniform, normal, and lognormal distributions for `float` produce incorrect values for Sobol 64 and Scrambled Sobol 64.
 
 ## (Unreleased) rocRAND-2.10.15 for ROCm 5.3.0
 ### Changed

library/include/rocrand/rocrand_normal.h

@@ -220,6 +220,14 @@ float normal_distribution(unsigned int x)
     return v;
 }
 
+FQUALIFIERS
+float normal_distribution(unsigned long long int x)
+{
+    float p = ::rocrand_device::detail::uniform_distribution(x);
+    float v = ROCRAND_SQRT2 * ::rocrand_device::detail::roc_f_erfinv(2.0f * p - 1.0f);
+    return v;
+}
+
 FQUALIFIERS
 float2 normal_distribution2(unsigned int v1, unsigned int v2)
 {

library/include/rocrand/rocrand_uniform.h

@@ -66,6 +66,14 @@ float uniform_distribution(unsigned int v)
     return ROCRAND_2POW32_INV + (v * ROCRAND_2POW32_INV);
 }
 
+// For unsigned integer between 0 and ULLONG_MAX, returns value between
+// 0.0f and 1.0f, excluding 0.0f and including 1.0f.
+FQUALIFIERS
+float uniform_distribution(unsigned long long int v)
+{
+    return ROCRAND_2POW32_INV + (v >> 32) * ROCRAND_2POW32_INV;
+}
+
 FQUALIFIERS
 float4 uniform_distribution4(uint4 v)
 {

python/rocrand/tests/rocrand_test.py

@@ -227,9 +227,8 @@ def _test_uniform(self, dtype):
         self.assertAlmostEqual(output.mean(), 0.5, delta=0.2)
         self.assertAlmostEqual(output.std(), math.sqrt(1 / 12.0), delta=0.2 * math.sqrt(1 / 12.0))
 
-    # TODO: temporarily disable float for uniform distribution
-    # def test_uniform_float(self):
-    #     self._test_uniform(np.float32)
+    def test_uniform_float(self):
+        self._test_uniform(np.float32)
 
     def test_uniform_double(self):
         self._test_uniform(np.float64)
@@ -241,9 +240,8 @@ def _test_normal_real(self, dtype):
         self.assertAlmostEqual(output.mean(), 0.0, delta=0.2)
         self.assertAlmostEqual(output.std(), 1.0, delta=0.2)
 
-    # TODO: temporarily disable float for normal distribution
-    # def test_normal_float(self):
-    #     self._test_normal_real(np.float32)
+    def test_normal_float(self):
+        self._test_normal_real(np.float32)
 
     def test_normal_double(self):
         self._test_normal_real(np.float64)
@@ -260,9 +258,8 @@ def _test_lognormal_real(self, dtype):
         self.assertAlmostEqual(logmean, 1.6, delta=1.6 * 0.2)
         self.assertAlmostEqual(logstd, 0.25, delta=0.25 * 0.2)
 
-    # TODO: temporarily disable float for lognormal distribution
-    # def test_lognormal_float(self):
-    #     self._test_lognormal_real(np.float32)
+    def test_lognormal_float(self):
+        self._test_lognormal_real(np.float32)
 
     def test_lognormal_double(self):
         self._test_lognormal_real(np.float64)

test/test_rocrand_kernel_scrambled_sobol64.cpp

@@ -229,27 +229,26 @@ TEST(rocrand_kernel_scrambled_sobol64, rocrand)
     EXPECT_NEAR(mean, 0.5, 0.1);
 }
 
-// TODO: temporarily disable float for uniform distribution
-// TEST(rocrand_kernel_scrambled_sobol64, rocrand_uniform)
-// {
-//     using ResultType   = float;
-//     using Distribution = rocrand_uniform_f;
-
-//     // amount of generated numbers has to be a multiple of the dimensions for sobol, so size is given per dimension
-//     constexpr size_t       size_per_dimension = 8192;
-//     constexpr unsigned int dimensions         = 8;
-
-//     std::vector<ResultType> output_host;
-//     call_rocrand_kernel<ResultType, Distribution>(output_host, dimensions, size_per_dimension);
-
-//     double mean = 0;
-//     for(ResultType v : output_host)
-//     {
-//         mean += static_cast<double>(v);
-//     }
-//     mean = mean / output_host.size();
-//     EXPECT_NEAR(mean, 0.5, 0.1);
-// }
+TEST(rocrand_kernel_scrambled_sobol64, rocrand_uniform)
+{
+    using ResultType   = float;
+    using Distribution = rocrand_uniform_f;
+
+    // amount of generated numbers has to be a multiple of the dimensions for sobol, so size is given per dimension
+    constexpr size_t       size_per_dimension = 8192;
+    constexpr unsigned int dimensions         = 8;
+
+    std::vector<ResultType> output_host;
+    call_rocrand_kernel<ResultType, Distribution>(output_host, dimensions, size_per_dimension);
+
+    double mean = 0;
+    for(ResultType v : output_host)
+    {
+        mean += static_cast<double>(v);
+    }
+    mean = mean / output_host.size();
+    EXPECT_NEAR(mean, 0.5, 0.1);
+}
 
 TEST(rocrand_kernel_scrambled_sobol64, rocrand_uniform_double)
 {
@@ -272,35 +271,34 @@ TEST(rocrand_kernel_scrambled_sobol64, rocrand_uniform_double)
     EXPECT_NEAR(mean, 0.5, 0.1);
 }
 
-// TODO: temporarily disable float for normal distribution
-// TEST(rocrand_kernel_scrambled_sobol64, rocrand_normal)
-// {
-//     using ResultType   = float;
-//     using Distribution = rocrand_normal_f;
-
-//     // amount of generated numbers has to be a multiple of the dimensions for sobol, so size is given per dimension
-//     constexpr size_t       size_per_dimension = 8192;
-//     constexpr unsigned int dimensions         = 8;
-
-//     std::vector<ResultType> output_host;
-//     call_rocrand_kernel<ResultType, Distribution>(output_host, dimensions, size_per_dimension);
-
-//     double mean = 0;
-//     for(ResultType v : output_host)
-//     {
-//         mean += static_cast<double>(v);
-//     }
-//     mean = mean / output_host.size();
-//     EXPECT_NEAR(mean, 0.0, 0.2);
-
-//     double stddev = 0;
-//     for(ResultType v : output_host)
-//     {
-//         stddev += std::pow(static_cast<double>(v) - mean, 2);
-//     }
-//     stddev = stddev / output_host.size();
-//     EXPECT_NEAR(stddev, 1.0, 0.2);
-// }
+TEST(rocrand_kernel_scrambled_sobol64, rocrand_normal)
+{
+    using ResultType   = float;
+    using Distribution = rocrand_normal_f;
+
+    // amount of generated numbers has to be a multiple of the dimensions for sobol, so size is given per dimension
+    constexpr size_t       size_per_dimension = 8192;
+    constexpr unsigned int dimensions         = 8;
+
+    std::vector<ResultType> output_host;
+    call_rocrand_kernel<ResultType, Distribution>(output_host, dimensions, size_per_dimension);
+
+    double mean = 0;
+    for(ResultType v : output_host)
+    {
+        mean += static_cast<double>(v);
+    }
+    mean = mean / output_host.size();
+    EXPECT_NEAR(mean, 0.0, 0.2);
+
+    double stddev = 0;
+    for(ResultType v : output_host)
+    {
+        stddev += std::pow(static_cast<double>(v) - mean, 2);
+    }
+    stddev = stddev / output_host.size();
+    EXPECT_NEAR(stddev, 1.0, 0.2);
+}
 
 TEST(rocrand_kernel_scrambled_sobol64, rocrand_normal_double)
 {

test/test_rocrand_scrambled_sobol64_qrng.cpp

@@ -38,9 +38,7 @@ struct rocrand_scrambled_sobol64_float_tests : public ::testing::Test
     using type = T;
 };
 
-// TODO: temporarily disable float for uniform and normal distributions
-// using FloatReturnTypes = ::testing::Types<float, double>;
-using FloatReturnTypes = ::testing::Types<double>;
+using FloatReturnTypes = ::testing::Types<float, double>;
 
 TYPED_TEST_SUITE(rocrand_scrambled_sobol64_float_tests, FloatReturnTypes);