Improve unit tests of Faure sequence

mbaudin47 · May 12, 2024 · 40ac801 · 40ac801
1 parent ea3ab35
commit 40ac801
Show file tree

Hide file tree

Showing 6 changed files with 154 additions and 2,029 deletions.
diff --git a/lib/test/CMakeLists.txt b/lib/test/CMakeLists.txt
@@ -225,8 +225,8 @@ ot_check_test (Field_draw)
 ot_check_test (Field_std)
 ot_check_test (ProcessSample_std)
 ot_check_test (RandomGenerator_std)
-ot_check_test (SobolSequence_std IGNOREOUT)
-ot_check_test (FaureSequence_std)
+ot_check_test (SobolSequence_std)
+ot_check_test (FaureSequence_std IGNOREOUT)
 ot_check_test (HaltonSequence_std IGNOREOUT)
 ot_check_test (HaselgroveSequence_std)
 ot_check_test (ReverseHaltonSequence_std)

diff --git a/lib/test/t_FaureSequence_std.cxx b/lib/test/t_FaureSequence_std.cxx
@@ -31,31 +31,74 @@ int main(int, char *[])
 
   try
   {
-    // Create a Faure sequence
+    // Create 8 points from a Faure sequence in dimension 1
+    Sample expected1D(8, 1);
+    expected1D[0] = Point({1.0 / 2.0});
+    expected1D[1] = Point({1.0 / 4.0});
+    expected1D[2] = Point({3.0 / 4.0});
+    expected1D[3] = Point({1.0 / 8.0});
+    expected1D[4] = Point({5.0 / 8.0});
+    expected1D[5] = Point({3.0 / 8.0});
+    expected1D[6] = Point({7.0 / 8.0});
+    expected1D[7] = Point({1.0 / 16.0});
+    FaureSequence sequence1D(1);
+    fullprint << sequence1D << std::endl;
+    Sample faureSample1D(sequence1D.generate(8));
+    assert_almost_equal(faureSample1D, expected1D);
+
+    // Create 8 points from a Faure sequence in dimension 2
+    Sample expected2D(8, 2);
+    expected2D[0] = Point({0.5, 0.5});
+    expected2D[1] = Point({0.25, 0.75});
+    expected2D[2] = Point({0.75, 0.25});
+    expected2D[3] = Point({0.125, 0.625});
+    expected2D[4] = Point({0.625, 0.125});
+    expected2D[5] = Point({0.375, 0.375});
+    expected2D[6] = Point({0.875, 0.875});
+    expected2D[7] = Point({0.0625, 0.9375});
+    FaureSequence sequence2D(2);
+    fullprint << sequence2D << std::endl;
+    Sample faureSample2D(sequence2D.generate(8));
+    assert_almost_equal(faureSample2D, expected2D);
+
+    // Create 8 points from a Faure sequence in dimension 3
+    Sample expected3D(8, 3);
+    expected3D[0] = Point({1.0 / 3.0, 1.0 / 3.0, 1.0 / 3.0});
+    expected3D[1] = Point({2.0 / 3.0, 2.0 / 3.0, 2.0 / 3.0});
+    expected3D[2] = Point({1.0 / 9.0, 4.0 / 9.0, 7.0 / 9.0});
+    expected3D[3] = Point({4.0 / 9.0, 7.0 / 9.0, 1.0 / 9.0});
+    expected3D[4] = Point({7.0 / 9.0, 1.0 / 9.0, 4.0 / 9.0});
+    expected3D[5] = Point({2.0 / 9.0, 8.0 / 9.0, 5.0 / 9.0});
+    expected3D[6] = Point({5.0 / 9.0, 2.0 / 9.0, 8.0 / 9.0});
+    expected3D[7] = Point({8.0 / 9.0, 5.0 / 9.0, 2.0 / 9.0});
+    FaureSequence sequence3D(3);
+    fullprint << sequence3D << std::endl;
+    Sample faureSample3D(sequence3D.generate(8));
+    assert_almost_equal(faureSample3D, expected3D);
+
+    // Create a Faure sequence in dimension 15
     FaureSequence sequence(15);
     fullprint << sequence << std::endl;
-
-    // Create a numerical sample of the sequence
-    Sample faureSample(sequence.generate(10));
-    fullprint << faureSample << std::endl;
+    sequence.generate(10);
 
     // Create another Faure' sequence of dimension 2 to estimate Pi in [0; 1)^2
     UnsignedInteger dimension = 2;
     sequence = FaureSequence(dimension);
     UnsignedInteger pointInsideCircle = 0;
-    UnsignedInteger sampleSize = 1000;
+    UnsignedInteger sampleSize = std::pow(3, 7);  // This is significant!
     for(UnsignedInteger i = 0; i < sampleSize; ++i)
     {
       Point faurePoint(sequence.generate());
-      fullprint << faurePoint << std::endl;
       if(faurePoint.norm() < 1.0)
         ++ pointInsideCircle;
     }
     Scalar probabilityEstimate = 1.0 * pointInsideCircle / sampleSize;
     Scalar probability = M_PI / 4.0;
-    Scalar relativeError = std::abs(probability - probabilityEstimate) / probability;
     fullprint << "sample size=" << sampleSize << std::endl;
-    fullprint << "relative error to Pi=" << relativeError << std::endl;
+    fullprint << "computed probability =" << probabilityEstimate << std::endl;
+    fullprint << "expected probability =" << probability << std::endl;
+    Scalar rtol = 10.0 / sampleSize;
+    assert_almost_equal(probabilityEstimate, probability, rtol);
   }
 
   catch (TestFailed & ex)