]>
git.proxmox.com Git - ceph.git/blob - ceph/src/boost/libs/math/test/test_toms748_solve.cpp
1 // (C) Copyright John Maddock 2006.
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)
8 #define BOOST_TEST_MAIN
9 #include <boost/test/unit_test.hpp>
10 #include <boost/test/tools/floating_point_comparison.hpp>
12 #include <boost/math/tools/toms748_solve.hpp>
13 #include <boost/math/special_functions/gamma.hpp>
14 #include <boost/math/special_functions/beta.hpp>
19 // Test functor implements the same test cases as used by
20 // "Algorithm 748: Enclosing Zeros of Continuous Functions"
21 // by G.E. Alefeld, F.A. Potra and Yixun Shi.
23 // Plus two more: one for inverting the incomplete gamma,
24 // and one for inverting the incomplete beta.
29 toms748tester(unsigned i
) : k(i
), ip(0), a(0), b(0){}
30 toms748tester(unsigned i
, int ip_
) : k(i
), ip(ip_
), a(0), b(0){}
31 toms748tester(unsigned i
, T a_
, T b_
) : k(i
), ip(0), a(a_
), b(b_
){}
33 static unsigned total_calls()
50 return sin(x
) - x
/ 2;
54 for(int i
= 1; i
<= 20; ++i
)
58 r
+= p
* p
/ (q
* q
* q
);
64 return a
* x
* exp(b
* x
);
70 return 2 * x
* exp(-T(ip
)) - 2 * exp(-ip
* x
) + 1;
72 return (1 + (1 - ip
) * (1 - ip
)) * x
- (1 - ip
* x
) * (1 - ip
* x
);
74 return x
* x
- pow(1 - x
, a
);
76 return (1 + (1 - ip
) * (1 - ip
) * (1 - ip
) * (1 - ip
)) * x
- (1 - ip
* x
) * (1 - ip
* x
) * (1 - ip
* x
) * (1 - ip
* x
);
78 return exp(-ip
* x
) * (x
- 1) + pow(x
, T(ip
));
80 return (ip
* x
- 1) / ((ip
- 1) * x
);
82 return pow(x
, T(1)/ip
) - pow(T(ip
), T(1)/ip
);
84 return x
== 0 ? 0 : x
* exp(-1 / (x
* x
));
86 return x
>= 0 ? (T(ip
)/20) * (x
/ 1.5f
+ sin(x
) - 1) : -T(ip
)/20;
91 return exp(1.0) - 1.859;
93 return exp((ip
+1)*x
*1000 / 2) - 1.859;
98 return boost::math::gamma_q(x
, a
) - b
;
101 return boost::math::ibeta(x
, a
, 0.5) - b
;
106 int k
; // index of problem.
107 int ip
; // integer parameter
108 T a
, b
; // real parameter
110 static unsigned invocations
;
114 unsigned toms748tester
<T
>::invocations
= 0;
116 boost::uintmax_t total
= 0;
117 boost::uintmax_t invocations
= 0;
120 void run_test(T a
, T b
, int id
)
122 boost::uintmax_t c
= 1000;
123 std::pair
<double, double> r
= toms748_solve(toms748tester
<double>(id
),
126 boost::math::tools::eps_tolerance
<double>(std::numeric_limits
<double>::digits
),
128 BOOST_CHECK_EQUAL(c
, toms748tester
<double>::total_calls());
130 invocations
+= toms748tester
<double>::total_calls();
131 std::cout
<< "Function " << id
<< "\nresult={" << r
.first
<< ", " << r
.second
<< "} total calls=" << toms748tester
<double>::total_calls() << "\n\n";
132 toms748tester
<double>::reset();
136 void run_test(T a
, T b
, int id
, int p
)
138 boost::uintmax_t c
= 1000;
139 std::pair
<double, double> r
= toms748_solve(toms748tester
<double>(id
, p
),
142 boost::math::tools::eps_tolerance
<double>(std::numeric_limits
<double>::digits
),
144 BOOST_CHECK_EQUAL(c
, toms748tester
<double>::total_calls());
146 invocations
+= toms748tester
<double>::total_calls();
147 std::cout
<< "Function " << id
<< "\nresult={" << r
.first
<< ", " << r
.second
<< "} total calls=" << toms748tester
<double>::total_calls() << "\n\n";
148 toms748tester
<double>::reset();
152 void run_test(T a
, T b
, int id
, T p1
, T p2
)
154 boost::uintmax_t c
= 1000;
155 std::pair
<double, double> r
= toms748_solve(toms748tester
<double>(id
, p1
, p2
),
158 boost::math::tools::eps_tolerance
<double>(std::numeric_limits
<double>::digits
),
160 BOOST_CHECK_EQUAL(c
, toms748tester
<double>::total_calls());
162 invocations
+= toms748tester
<double>::total_calls();
163 std::cout
<< "Function " << id
<< "\n Result={" << r
.first
<< ", " << r
.second
<< "} total calls=" << toms748tester
<double>::total_calls() << "\n\n";
164 toms748tester
<double>::reset();
167 BOOST_AUTO_TEST_CASE( test_main
)
169 std::cout
<< std::setprecision(18);
170 run_test(3.14/2, 3.14, 1);
172 for(int i
= 1; i
<= 10; i
+= 1)
174 run_test(i
*i
+ 1e-9, (i
+1)*(i
+1) - 1e-9, 2);
177 run_test(-9.0, 31.0, 3, -40.0, -1.0);
178 run_test(-9.0, 31.0, 3, -100.0, -2.0);
179 run_test(-9.0, 31.0, 3, -200.0, -3.0);
181 for(int n
= 4; n
<= 12; n
+= 2)
183 run_test(0.0, 5.0, 4, 0.2, double(n
));
185 for(int n
= 4; n
<= 12; n
+= 2)
187 run_test(0.0, 5.0, 4, 1.0, double(n
));
189 for(int n
= 8; n
<= 14; n
+= 2)
191 run_test(-0.95, 4.05, 4, 1.0, double(n
));
193 run_test(0.0, 1.5, 5);
194 for(int n
= 1; n
<= 5; ++n
)
196 run_test(0.0, 1.0, 6, n
);
198 for(int n
= 20; n
<= 100; n
+= 20)
200 run_test(0.0, 1.0, 6, n
);
202 run_test(0.0, 1.0, 7, 5);
203 run_test(0.0, 1.0, 7, 10);
204 run_test(0.0, 1.0, 7, 20);
205 run_test(0.0, 1.0, 8, 2);
206 run_test(0.0, 1.0, 8, 5);
207 run_test(0.0, 1.0, 8, 10);
208 run_test(0.0, 1.0, 8, 15);
209 run_test(0.0, 1.0, 8, 20);
210 run_test(0.0, 1.0, 9, 1);
211 run_test(0.0, 1.0, 9, 2);
212 run_test(0.0, 1.0, 9, 4);
213 run_test(0.0, 1.0, 9, 5);
214 run_test(0.0, 1.0, 9, 8);
215 run_test(0.0, 1.0, 9, 15);
216 run_test(0.0, 1.0, 9, 20);
217 run_test(0.0, 1.0, 10, 1);
218 run_test(0.0, 1.0, 10, 5);
219 run_test(0.0, 1.0, 10, 10);
220 run_test(0.0, 1.0, 10, 15);
221 run_test(0.0, 1.0, 10, 20);
223 run_test(0.01, 1.0, 11, 2);
224 run_test(0.01, 1.0, 11, 5);
225 run_test(0.01, 1.0, 11, 15);
226 run_test(0.01, 1.0, 11, 20);
228 for(int n
= 2; n
<= 6; ++n
)
229 run_test(1.0, 100.0, 12, n
);
230 for(int n
= 7; n
<= 33; n
+=2)
231 run_test(1.0, 100.0, 12, n
);
233 run_test(-1.0, 4.0, 13);
235 for(int n
= 1; n
<= 40; ++n
)
236 run_test(-1e4
, 3.14/2, 14, n
);
238 for(int n
= 20; n
<= 40; ++n
)
239 run_test(-1e4
, 1e-4, 15, n
);
241 for(int n
= 100; n
<= 1000; n
+=100)
242 run_test(-1e4
, 1e-4, 15, n
);
244 std::cout
<< "Total iterations consumed = " << total
<< std::endl
;
245 std::cout
<< "Total function invocations consumed = " << invocations
<< std::endl
<< std::endl
;
247 BOOST_CHECK(invocations
< 3150);
249 std::cout
<< std::setprecision(18);
251 for(int n
= 5; n
<= 100; n
+= 10)
252 run_test(sqrt(double(n
)), double(n
+1), 16, (double)n
, 0.4);
254 for(int n
= 5; n
<= 100; n
+= 10)
255 run_test(double(n
/ 2), double(2*n
), 17, (double)n
, 0.4);
258 for(int n
= 4; n
<= 12; n
+= 2)
260 boost::uintmax_t c
= 1000;
261 std::pair
<double, double> r
= bracket_and_solve_root(toms748tester
<double>(4, 0.2, double(n
)),
265 boost::math::tools::eps_tolerance
<double>(std::numeric_limits
<double>::digits
),
267 std::cout
<< std::setprecision(18);
268 std::cout
<< "Function " << 4 << "\n Result={" << r
.first
<< ", " << r
.second
<< "} total calls=" << toms748tester
<double>::total_calls() << "\n\n";
269 toms748tester
<double>::reset();