mirror of
https://github.com/saitohirga/WSJT-X.git
synced 2025-10-30 12:30:23 -04:00
274 lines
7.7 KiB
C++
274 lines
7.7 KiB
C++
// (C) Copyright John Maddock 2006.
|
|
// Use, modification and distribution are subject to the
|
|
// Boost Software License, Version 1.0. (See accompanying file
|
|
// LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
|
|
|
|
#include <pch.hpp>
|
|
|
|
#define BOOST_TEST_MAIN
|
|
#include <boost/test/unit_test.hpp>
|
|
#include <boost/test/floating_point_comparison.hpp>
|
|
|
|
#include <boost/math/tools/toms748_solve.hpp>
|
|
#include <boost/math/special_functions/gamma.hpp>
|
|
#include <boost/math/special_functions/beta.hpp>
|
|
#include <iostream>
|
|
#include <iomanip>
|
|
|
|
//
|
|
// Test functor implements the same test cases as used by
|
|
// "Algorithm 748: Enclosing Zeros of Continuous Functions"
|
|
// by G.E. Alefeld, F.A. Potra and Yixun Shi.
|
|
//
|
|
// Plus two more: one for inverting the incomplete gamma,
|
|
// and one for inverting the incomplete beta.
|
|
//
|
|
template <class T>
|
|
struct toms748tester
|
|
{
|
|
toms748tester(unsigned i) : k(i), ip(0), a(0), b(0){}
|
|
toms748tester(unsigned i, int ip_) : k(i), ip(ip_), a(0), b(0){}
|
|
toms748tester(unsigned i, T a_, T b_) : k(i), ip(0), a(a_), b(b_){}
|
|
|
|
static unsigned total_calls()
|
|
{
|
|
return invocations;
|
|
}
|
|
static void reset()
|
|
{
|
|
invocations = 0;
|
|
}
|
|
|
|
T operator()(T x)
|
|
{
|
|
using namespace std;
|
|
|
|
++invocations;
|
|
switch(k)
|
|
{
|
|
case 1:
|
|
return sin(x) - x / 2;
|
|
case 2:
|
|
{
|
|
T r = 0;
|
|
for(int i = 1; i <= 20; ++i)
|
|
{
|
|
T p = (2 * i - 5);
|
|
T q = x - i * i;
|
|
r += p * p / (q * q * q);
|
|
}
|
|
r *= -2;
|
|
return r;
|
|
}
|
|
case 3:
|
|
return a * x * exp(b * x);
|
|
case 4:
|
|
return pow(x, b) - a;
|
|
case 5:
|
|
return sin(x) - 0.5;
|
|
case 6:
|
|
return 2 * x * exp(-T(ip)) - 2 * exp(-ip * x) + 1;
|
|
case 7:
|
|
return (1 + (1 - ip) * (1 - ip)) * x - (1 - ip * x) * (1 - ip * x);
|
|
case 8:
|
|
return x * x - pow(1 - x, a);
|
|
case 9:
|
|
return (1 + (1 - ip) * (1 - ip) * (1 - ip) * (1 - ip)) * x - (1 - ip * x) * (1 - ip * x) * (1 - ip * x) * (1 - ip * x);
|
|
case 10:
|
|
return exp(-ip * x) * (x - 1) + pow(x, T(ip));
|
|
case 11:
|
|
return (ip * x - 1) / ((ip - 1) * x);
|
|
case 12:
|
|
return pow(x, T(1)/ip) - pow(T(ip), T(1)/ip);
|
|
case 13:
|
|
return x == 0 ? 0 : x * exp(-1 / (x * x));
|
|
case 14:
|
|
return x >= 0 ? (T(ip)/20) * (x / 1.5f + sin(x) - 1) : -T(ip)/20;
|
|
case 15:
|
|
{
|
|
T d = 2e-3f / (1+ip);
|
|
if(x > d)
|
|
return exp(1.0) - 1.859;
|
|
else if(x > 0)
|
|
return exp((ip+1)*x*1000 / 2) - 1.859;
|
|
return -0.859f;
|
|
}
|
|
case 16:
|
|
{
|
|
return boost::math::gamma_q(x, a) - b;
|
|
}
|
|
case 17:
|
|
return boost::math::ibeta(x, a, 0.5) - b;
|
|
}
|
|
return 0;
|
|
}
|
|
private:
|
|
int k; // index of problem.
|
|
int ip; // integer parameter
|
|
T a, b; // real parameter
|
|
|
|
static unsigned invocations;
|
|
};
|
|
|
|
template <class T>
|
|
unsigned toms748tester<T>::invocations = 0;
|
|
|
|
boost::uintmax_t total = 0;
|
|
boost::uintmax_t invocations = 0;
|
|
|
|
template <class T>
|
|
void run_test(T a, T b, int id)
|
|
{
|
|
boost::uintmax_t c = 1000;
|
|
std::pair<double, double> r = toms748_solve(toms748tester<double>(id),
|
|
a,
|
|
b,
|
|
boost::math::tools::eps_tolerance<double>(std::numeric_limits<double>::digits),
|
|
c);
|
|
BOOST_CHECK_EQUAL(c, toms748tester<double>::total_calls());
|
|
total += c;
|
|
invocations += toms748tester<double>::total_calls();
|
|
std::cout << "Function " << id << "\nresult={" << r.first << ", " << r.second << "} total calls=" << toms748tester<double>::total_calls() << "\n\n";
|
|
toms748tester<double>::reset();
|
|
}
|
|
|
|
template <class T>
|
|
void run_test(T a, T b, int id, int p)
|
|
{
|
|
boost::uintmax_t c = 1000;
|
|
std::pair<double, double> r = toms748_solve(toms748tester<double>(id, p),
|
|
a,
|
|
b,
|
|
boost::math::tools::eps_tolerance<double>(std::numeric_limits<double>::digits),
|
|
c);
|
|
BOOST_CHECK_EQUAL(c, toms748tester<double>::total_calls());
|
|
total += c;
|
|
invocations += toms748tester<double>::total_calls();
|
|
std::cout << "Function " << id << "\nresult={" << r.first << ", " << r.second << "} total calls=" << toms748tester<double>::total_calls() << "\n\n";
|
|
toms748tester<double>::reset();
|
|
}
|
|
|
|
template <class T>
|
|
void run_test(T a, T b, int id, T p1, T p2)
|
|
{
|
|
boost::uintmax_t c = 1000;
|
|
std::pair<double, double> r = toms748_solve(toms748tester<double>(id, p1, p2),
|
|
a,
|
|
b,
|
|
boost::math::tools::eps_tolerance<double>(std::numeric_limits<double>::digits),
|
|
c);
|
|
BOOST_CHECK_EQUAL(c, toms748tester<double>::total_calls());
|
|
total += c;
|
|
invocations += toms748tester<double>::total_calls();
|
|
std::cout << "Function " << id << "\n Result={" << r.first << ", " << r.second << "} total calls=" << toms748tester<double>::total_calls() << "\n\n";
|
|
toms748tester<double>::reset();
|
|
}
|
|
|
|
BOOST_AUTO_TEST_CASE( test_main )
|
|
{
|
|
std::cout << std::setprecision(18);
|
|
run_test(3.14/2, 3.14, 1);
|
|
|
|
for(int i = 1; i <= 10; i += 1)
|
|
{
|
|
run_test(i*i + 1e-9, (i+1)*(i+1) - 1e-9, 2);
|
|
}
|
|
|
|
run_test(-9.0, 31.0, 3, -40.0, -1.0);
|
|
run_test(-9.0, 31.0, 3, -100.0, -2.0);
|
|
run_test(-9.0, 31.0, 3, -200.0, -3.0);
|
|
|
|
for(int n = 4; n <= 12; n += 2)
|
|
{
|
|
run_test(0.0, 5.0, 4, 0.2, double(n));
|
|
}
|
|
for(int n = 4; n <= 12; n += 2)
|
|
{
|
|
run_test(0.0, 5.0, 4, 1.0, double(n));
|
|
}
|
|
for(int n = 8; n <= 14; n += 2)
|
|
{
|
|
run_test(-0.95, 4.05, 4, 1.0, double(n));
|
|
}
|
|
run_test(0.0, 1.5, 5);
|
|
for(int n = 1; n <= 5; ++n)
|
|
{
|
|
run_test(0.0, 1.0, 6, n);
|
|
}
|
|
for(int n = 20; n <= 100; n += 20)
|
|
{
|
|
run_test(0.0, 1.0, 6, n);
|
|
}
|
|
run_test(0.0, 1.0, 7, 5);
|
|
run_test(0.0, 1.0, 7, 10);
|
|
run_test(0.0, 1.0, 7, 20);
|
|
run_test(0.0, 1.0, 8, 2);
|
|
run_test(0.0, 1.0, 8, 5);
|
|
run_test(0.0, 1.0, 8, 10);
|
|
run_test(0.0, 1.0, 8, 15);
|
|
run_test(0.0, 1.0, 8, 20);
|
|
run_test(0.0, 1.0, 9, 1);
|
|
run_test(0.0, 1.0, 9, 2);
|
|
run_test(0.0, 1.0, 9, 4);
|
|
run_test(0.0, 1.0, 9, 5);
|
|
run_test(0.0, 1.0, 9, 8);
|
|
run_test(0.0, 1.0, 9, 15);
|
|
run_test(0.0, 1.0, 9, 20);
|
|
run_test(0.0, 1.0, 10, 1);
|
|
run_test(0.0, 1.0, 10, 5);
|
|
run_test(0.0, 1.0, 10, 10);
|
|
run_test(0.0, 1.0, 10, 15);
|
|
run_test(0.0, 1.0, 10, 20);
|
|
|
|
run_test(0.01, 1.0, 11, 2);
|
|
run_test(0.01, 1.0, 11, 5);
|
|
run_test(0.01, 1.0, 11, 15);
|
|
run_test(0.01, 1.0, 11, 20);
|
|
|
|
for(int n = 2; n <= 6; ++n)
|
|
run_test(1.0, 100.0, 12, n);
|
|
for(int n = 7; n <= 33; n+=2)
|
|
run_test(1.0, 100.0, 12, n);
|
|
|
|
run_test(-1.0, 4.0, 13);
|
|
|
|
for(int n = 1; n <= 40; ++n)
|
|
run_test(-1e4, 3.14/2, 14, n);
|
|
|
|
for(int n = 20; n <= 40; ++n)
|
|
run_test(-1e4, 1e-4, 15, n);
|
|
|
|
for(int n = 100; n <= 1000; n+=100)
|
|
run_test(-1e4, 1e-4, 15, n);
|
|
|
|
std::cout << "Total iterations consumed = " << total << std::endl;
|
|
std::cout << "Total function invocations consumed = " << invocations << std::endl << std::endl;
|
|
|
|
BOOST_CHECK(invocations < 3150);
|
|
|
|
std::cout << std::setprecision(18);
|
|
|
|
for(int n = 5; n <= 100; n += 10)
|
|
run_test(sqrt(double(n)), double(n+1), 16, (double)n, 0.4);
|
|
|
|
for(int n = 5; n <= 100; n += 10)
|
|
run_test(double(n / 2), double(2*n), 17, (double)n, 0.4);
|
|
|
|
|
|
for(int n = 4; n <= 12; n += 2)
|
|
{
|
|
boost::uintmax_t c = 1000;
|
|
std::pair<double, double> r = bracket_and_solve_root(toms748tester<double>(4, 0.2, double(n)),
|
|
2.0,
|
|
2.0,
|
|
true,
|
|
boost::math::tools::eps_tolerance<double>(std::numeric_limits<double>::digits),
|
|
c);
|
|
std::cout << std::setprecision(18);
|
|
std::cout << "Function " << 4 << "\n Result={" << r.first << ", " << r.second << "} total calls=" << toms748tester<double>::total_calls() << "\n\n";
|
|
toms748tester<double>::reset();
|
|
BOOST_CHECK(c < 20);
|
|
}
|
|
}
|
|
|