ceph/src/boost/libs/math/reporting/performance/cohen_acceleration_performance.cpp

   1 //  (C) Copyright Nick Thompson 2020.
   2 //  Use, modification and distribution are subject to the
   3 //  Boost Software License, Version 1.0. (See accompanying file
   4 //  LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
   5
   6 #include <random>
   7 #include <iostream>
   8 #include <iomanip>
   9 #include <benchmark/benchmark.h>
  10 #include <boost/math/tools/cohen_acceleration.hpp>
  11 #include <boost/multiprecision/float128.hpp>
  12 #include <boost/multiprecision/mpfr.hpp>
  13 #include <boost/multiprecision/cpp_bin_float.hpp>
  14 #include <boost/core/demangle.hpp>
  15
  16 using boost::multiprecision::number;
  17 using boost::multiprecision::mpfr_float_backend;
  18 using boost::multiprecision::float128;
  19 using boost::multiprecision::cpp_bin_float_50;
  20 using boost::multiprecision::cpp_bin_float_100;
  21 using boost::math::tools::cohen_acceleration;
  22 using boost::math::constants::pi;
  23
  24 template<typename Real>
  25 class G {
  26 public:
  27     G(){
  28         k_ = 0;
  29     }
  30
  31     Real operator()() {
  32         k_ += 1;
  33         return 1/(k_*k_);
  34     }
  35
  36 private:
  37     Real k_;
  38 };
  39
  40
  41 template<class Real>
  42 void CohenAcceleration(benchmark::State& state)
  43 {
  44     using std::abs;
  45     Real x = pi<Real>()*pi<Real>()/12;
  46     for (auto _ : state)
  47     {
  48         auto g = G<Real>();
  49         x = cohen_acceleration(g);
  50         benchmark::DoNotOptimize(x);
  51     }
  52     if (abs(x - pi<Real>()*pi<Real>()/12) > 16*std::numeric_limits<Real>::epsilon())
  53     {
  54         std::cerr << std::setprecision(std::numeric_limits<Real>::max_digits10);
  55         std::cerr << "Cohen acceleration computed " << x << " on type " << boost::core::demangle(typeid(Real).name()) << "\n";
  56         std::cerr << "But expected value is       " << pi<Real>()*pi<Real>()/12 << "\n";
  57     }
  58 }
  59
  60 BENCHMARK_TEMPLATE(CohenAcceleration, float);
  61 BENCHMARK_TEMPLATE(CohenAcceleration, double);
  62 BENCHMARK_TEMPLATE(CohenAcceleration, long double);
  63 BENCHMARK_TEMPLATE(CohenAcceleration, float128);
  64 BENCHMARK_TEMPLATE(CohenAcceleration, cpp_bin_float_50);
  65 BENCHMARK_TEMPLATE(CohenAcceleration, cpp_bin_float_100);
  66 BENCHMARK_TEMPLATE(CohenAcceleration, number<mpfr_float_backend<100>>);
  67 BENCHMARK_TEMPLATE(CohenAcceleration, number<mpfr_float_backend<200>>);
  68 BENCHMARK_TEMPLATE(CohenAcceleration, number<mpfr_float_backend<300>>);
  69 BENCHMARK_TEMPLATE(CohenAcceleration, number<mpfr_float_backend<400>>);
  70 BENCHMARK_TEMPLATE(CohenAcceleration, number<mpfr_float_backend<1000>>);
  71
  72
  73 template<class Real>
  74 void NaiveSum(benchmark::State& state)
  75 {
  76     using std::abs;
  77     Real x = pi<Real>()*pi<Real>()/12;
  78     for (auto _ : state)
  79     {
  80         auto g = G<Real>();
  81         Real term = g();
  82         x = term;
  83         bool even = false;
  84         while (term > std::numeric_limits<Real>::epsilon()/2) {
  85             term = g();
  86             if (even) {
  87                 x += term;
  88                 even = false;
  89             } else {
  90                 x -= term;
  91                 even = true;
  92             }
  93         }
  94         benchmark::DoNotOptimize(x);
  95     }
  96     // The accuracy tests don't pass because the sum is ill-conditioned:
  97     /*if (abs(x - pi<Real>()*pi<Real>()/12) > 16*std::numeric_limits<Real>::epsilon())
  98     {
  99         std::cerr << std::setprecision(std::numeric_limits<Real>::max_digits10);
 100         std::cerr << "Cohen acceleration computed " << x << " on type " << boost::core::demangle(typeid(Real).name()) << "\n";
 101         std::cerr << "But expected value is       " << pi<Real>()*pi<Real>()/12 << "\n";
 102     }*/
 103 }
 104
 105 BENCHMARK_TEMPLATE(NaiveSum, float);
 106 BENCHMARK_TEMPLATE(NaiveSum, double);
 107 BENCHMARK_TEMPLATE(NaiveSum, long double);
 108 BENCHMARK_TEMPLATE(NaiveSum, float128);
 109 BENCHMARK_TEMPLATE(NaiveSum, cpp_bin_float_50);
 110 BENCHMARK_TEMPLATE(NaiveSum, cpp_bin_float_100);
 111 BENCHMARK_TEMPLATE(NaiveSum, number<mpfr_float_backend<100>>);
 112 BENCHMARK_TEMPLATE(NaiveSum, number<mpfr_float_backend<200>>);
 113 BENCHMARK_TEMPLATE(NaiveSum, number<mpfr_float_backend<300>>);
 114 BENCHMARK_TEMPLATE(NaiveSum, number<mpfr_float_backend<400>>);
 115 BENCHMARK_TEMPLATE(NaiveSum, number<mpfr_float_backend<1000>>);
 116
 117 BENCHMARK_MAIN();