[lmi-commits] [lmi] master 8f2f355 19/23: Augment unit tests; record some observations

[Greg Chicares]

branch: master
commit 8f2f355bec2d4e8a82183a0a09b28a5a47e9e84c
Author: Gregory W. Chicares <gchicares@sbcglobal.net>
Commit: Gregory W. Chicares <gchicares@sbcglobal.net>

    Augment unit tests; record some observations
    
    Statistics for pure bisection, using two different means, with
    trivial patches inlined in this commit message:
    
      evaluations  max   mean  sample-std-dev    commit        date
          7212      63   18.6       6.94       d6bd8029e  20210718T1630Z
          7501      75   19.3       7.36       d2dcbf860  20210723T2039Z
          7331      64   18.9       7.12       dc63e62cd  20210724T2119Z
         14987      39   38.6       0.76       (1) below: arithmetic mean
         10543      35   27.2       2.31       (2) below: binary64 mean
    
    The observations added in the unit test agree closely with the results
    tabulated above for the arithmetic mean. The binary64 mean performs
    better than the unit-test discussion might suggest, presumably because
    the objective function performs currency rounding: e.g., there are many
    distinct binary64 values in [$0.00, $0.01], but only two rounded values.
    
    (1) Using arithmetic-mean bisection only:
    
    diff --git a/zero.hpp b/zero.hpp
    index 1b89b0f92..6c1a5fade 100644
    --- a/zero.hpp
    +++ b/zero.hpp
    @@ -492,7 +492,7 @@ root_type lmi_root
                     ; // Do nothing.
                     }
                 }
    -        if(std::fabs(e) < tol)
    +        if(true || std::fabs(e) < tol)
                 {
                 impetus = dithering_near_root;
                 d = e = n - b;
    
    (2) Using binary64 bisection only:
    
    ------------------------ >8 ------------------------
    diff --git a/zero.hpp b/zero.hpp
    index 1b89b0f92..2e13df69b 100644
    --- a/zero.hpp
    +++ b/zero.hpp
    @@ -492,7 +492,13 @@ root_type lmi_root
                     ; // Do nothing.
                     }
                 }
    -        if(std::fabs(e) < tol)
    +        if(true)
    +            {
    +            impetus = dithering_near_root;
    +            n = binary64_midpoint(b, c); // "next" iterate
    +            d = e = n - b;
    +            }
    +        else if(std::fabs(e) < tol)
                 {
                 impetus = dithering_near_root;
                 d = e = n - b;
    ------------------------ >8 ------------------------
    [This is the end of the commit message.]
---
 zero_test.cpp | 26 ++++++++++++++++++++++++++
 1 file changed, 26 insertions(+)

diff --git a/zero_test.cpp b/zero_test.cpp
index 28696cc..29a2333 100644
--- a/zero_test.cpp
+++ b/zero_test.cpp
@@ -474,6 +474,32 @@ void test_binary64_midpoint()
 
     LMI_TEST_EQUAL( 1.0e100, binary64_midpoint( 1.0e100,  1.0e100));
     LMI_TEST_EQUAL(-1.0e100, binary64_midpoint(-1.0e100, -1.0e100));
+
+    // Illustration solves typically search in [0, 1.0e9]. For
+    // binary64 bisection (only):
+    //  - using 1.0e9 instead of DBL_MAX saves only about one function
+    //    evaluation (but risks overflow), so an even more implausible
+    //    upper limit like 1.0e18 would cost little;
+    //  - sometimes 0.0 is the correct answer, but the next higher
+    //    currency amount is $0.01, which is very far from zero: the
+    //    interval [1.0e-2, 1.0e9] can be searched exhaustively in
+    //    about fifty-seven function evaluations.
+    std::uint64_t bignum        = 0x7FEFFFFFFFFFFFFF;
+    std::uint64_t one_e_300     = 0x7E37E43C8800759C;
+    std::uint64_t one_billion   = 0x41CDCD6500000000;
+    std::uint64_t one_hundredth = 0x3F847AE147AE147B;
+    LMI_TEST_EQUAL(4741671816366391296, one_billion);
+    LMI_TEST_EQUAL(4576918229304087675, one_hundredth);
+    LMI_TEST(materially_equal(62.9993, std::log2(bignum       ), 1.0e-4));
+    LMI_TEST(materially_equal(62.9798, std::log2(one_e_300    ), 1.0e-4));
+    LMI_TEST(materially_equal(62.0401, std::log2(one_billion  ), 1.0e-4));
+    LMI_TEST(materially_equal(61.9891, std::log2(one_hundredth), 1.0e-4));
+    LMI_TEST(materially_equal(57.1931, std::log2(one_billion - one_hundredth), 
1.0e-4));
+    // The same [0, 1.0e9] interval could be searched exhaustively for
+    // integral cents in fewer iterations using the arithmetic mean:
+    LMI_TEST(materially_equal(3.49808e-150, binary64_midpoint(0.0, 1.0e9), 
1.0e-5));
+    LMI_TEST_EQUAL(39, max_n_eval_bolzano(0.0, 1.0e9, 0.005, 1.0e9));
+    LMI_TEST_EQUAL(39, max_n_eval_bolzano(0.0, 1.0e9, 0.005, 0.0));
 }
 
 /// A function whose value almost everywhere in (-1.0e100, 1.0e100)