]>
Commit | Line | Data |
---|---|---|
1 | // Copyright John Maddock 2006. | |
2 | // Copyright Paul A. Bristow 2007, 2009 | |
3 | // Use, modification and distribution are subject to the | |
4 | // Boost Software License, Version 1.0. (See accompanying file | |
5 | // LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt) | |
6 | ||
7 | #include <boost/math/concepts/real_concept.hpp> | |
8 | #define BOOST_TEST_MAIN | |
9 | #include <boost/test/unit_test.hpp> | |
10 | #include <boost/test/floating_point_comparison.hpp> | |
11 | #include <boost/math/special_functions/math_fwd.hpp> | |
12 | #include <boost/math/tools/stats.hpp> | |
13 | #include <boost/math/tools/test.hpp> | |
14 | #include <boost/math/constants/constants.hpp> | |
15 | #include <boost/type_traits/is_floating_point.hpp> | |
16 | #include <boost/array.hpp> | |
17 | #include "functor.hpp" | |
18 | ||
19 | #include "handle_test_result.hpp" | |
20 | #include "table_type.hpp" | |
21 | ||
22 | #ifndef SC_ | |
23 | #define SC_(x) static_cast<typename table_type<T>::type>(BOOST_JOIN(x, L)) | |
24 | #endif | |
25 | ||
26 | template <class Real, class T> | |
27 | void test_inverses(const T& data) | |
28 | { | |
29 | using namespace std; | |
30 | //typedef typename T::value_type row_type; | |
31 | typedef Real value_type; | |
32 | ||
33 | value_type precision = static_cast<value_type>(ldexp(1.0, 1-boost::math::policies::digits<value_type, boost::math::policies::policy<> >()/2)) * 100; | |
34 | if(boost::math::policies::digits<value_type, boost::math::policies::policy<> >() < 50) | |
35 | precision = 1; // 1% or two decimal digits, all we can hope for when the input is truncated | |
36 | ||
37 | for(unsigned i = 0; i < data.size(); ++i) | |
38 | { | |
39 | // | |
40 | // These inverse tests are thrown off if the output of the | |
41 | // incomplete beta is too close to 1: basically there is insuffient | |
42 | // information left in the value we're using as input to the inverse | |
43 | // to be able to get back to the original value. | |
44 | // | |
45 | if(Real(data[i][5]) == 0) | |
46 | BOOST_CHECK_EQUAL(boost::math::ibeta_inv(Real(data[i][0]), Real(data[i][1]), Real(data[i][5])), value_type(0)); | |
47 | else if((1 - Real(data[i][5]) > 0.001) | |
48 | && (fabs(Real(data[i][5])) > 2 * boost::math::tools::min_value<value_type>()) | |
49 | && (fabs(Real(data[i][5])) > 2 * boost::math::tools::min_value<double>())) | |
50 | { | |
51 | value_type inv = boost::math::ibeta_inv(Real(data[i][0]), Real(data[i][1]), Real(data[i][5])); | |
52 | BOOST_CHECK_CLOSE(Real(data[i][2]), inv, precision); | |
53 | } | |
54 | else if(1 == Real(data[i][5])) | |
55 | BOOST_CHECK_EQUAL(boost::math::ibeta_inv(Real(data[i][0]), Real(data[i][1]), Real(data[i][5])), value_type(1)); | |
56 | ||
57 | if(Real(data[i][6]) == 0) | |
58 | BOOST_CHECK_EQUAL(boost::math::ibetac_inv(Real(data[i][0]), Real(data[i][1]), Real(data[i][6])), value_type(1)); | |
59 | else if((1 - Real(data[i][6]) > 0.001) | |
60 | && (fabs(Real(data[i][6])) > 2 * boost::math::tools::min_value<value_type>()) | |
61 | && (fabs(Real(data[i][6])) > 2 * boost::math::tools::min_value<double>())) | |
62 | { | |
63 | value_type inv = boost::math::ibetac_inv(Real(data[i][0]), Real(data[i][1]), Real(data[i][6])); | |
64 | BOOST_CHECK_CLOSE(Real(data[i][2]), inv, precision); | |
65 | } | |
66 | else if(Real(data[i][6]) == 1) | |
67 | BOOST_CHECK_EQUAL(boost::math::ibetac_inv(Real(data[i][0]), Real(data[i][1]), Real(data[i][6])), value_type(0)); | |
68 | } | |
69 | } | |
70 | ||
71 | template <class Real, class T> | |
72 | void test_inverses2(const T& data, const char* type_name, const char* test_name) | |
73 | { | |
74 | #if !(defined(ERROR_REPORTING_MODE) && !defined(IBETA_INV_FUNCTION_TO_TEST)) | |
75 | //typedef typename T::value_type row_type; | |
76 | typedef Real value_type; | |
77 | ||
78 | typedef value_type (*pg)(value_type, value_type, value_type); | |
79 | #ifdef IBETA_INV_FUNCTION_TO_TEST | |
80 | pg funcp = IBETA_INV_FUNCTION_TO_TEST; | |
81 | #elif defined(BOOST_MATH_NO_DEDUCED_FUNCTION_POINTERS) | |
82 | pg funcp = boost::math::ibeta_inv<value_type, value_type, value_type>; | |
83 | #else | |
84 | pg funcp = boost::math::ibeta_inv; | |
85 | #endif | |
86 | ||
87 | boost::math::tools::test_result<value_type> result; | |
88 | ||
89 | std::cout << "Testing " << test_name << " with type " << type_name | |
90 | << "\n~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~\n"; | |
91 | ||
92 | // | |
93 | // test ibeta_inv(T, T, T) against data: | |
94 | // | |
95 | result = boost::math::tools::test_hetero<Real>( | |
96 | data, | |
97 | bind_func<Real>(funcp, 0, 1, 2), | |
98 | extract_result<Real>(3)); | |
99 | handle_test_result(result, data[result.worst()], result.worst(), type_name, "ibeta_inv", test_name); | |
100 | // | |
101 | // test ibetac_inv(T, T, T) against data: | |
102 | // | |
103 | #ifdef IBETAC_INV_FUNCTION_TO_TEST | |
104 | funcp = IBETAC_INV_FUNCTION_TO_TEST; | |
105 | #elif defined(BOOST_MATH_NO_DEDUCED_FUNCTION_POINTERS) | |
106 | funcp = boost::math::ibetac_inv<value_type, value_type, value_type>; | |
107 | #else | |
108 | funcp = boost::math::ibetac_inv; | |
109 | #endif | |
110 | result = boost::math::tools::test_hetero<Real>( | |
111 | data, | |
112 | bind_func<Real>(funcp, 0, 1, 2), | |
113 | extract_result<Real>(4)); | |
114 | handle_test_result(result, data[result.worst()], result.worst(), type_name, "ibetac_inv", test_name); | |
115 | #endif | |
116 | } | |
117 | ||
118 | ||
119 | template <class T> | |
120 | void test_beta(T, const char* name) | |
121 | { | |
122 | #if !defined(ERROR_REPORTING_MODE) | |
123 | (void)name; | |
124 | // | |
125 | // The actual test data is rather verbose, so it's in a separate file | |
126 | // | |
127 | // The contents are as follows, each row of data contains | |
128 | // five items, input value a, input value b, integration limits x, beta(a, b, x) and ibeta(a, b, x): | |
129 | // | |
130 | #if !defined(TEST_DATA) || (TEST_DATA == 1) | |
131 | # include "ibeta_small_data.ipp" | |
132 | ||
133 | test_inverses<T>(ibeta_small_data); | |
134 | #endif | |
135 | ||
136 | #if !defined(TEST_DATA) || (TEST_DATA == 2) | |
137 | # include "ibeta_data.ipp" | |
138 | ||
139 | test_inverses<T>(ibeta_data); | |
140 | #endif | |
141 | ||
142 | #if !defined(TEST_DATA) || (TEST_DATA == 3) | |
143 | # include "ibeta_large_data.ipp" | |
144 | ||
145 | test_inverses<T>(ibeta_large_data); | |
146 | #endif | |
147 | ||
148 | #endif | |
149 | ||
150 | #if !defined(TEST_DATA) || (TEST_DATA == 4) | |
151 | # include "ibeta_inv_data.ipp" | |
152 | ||
153 | test_inverses2<T>(ibeta_inv_data, name, "Inverse incomplete beta"); | |
154 | #endif | |
155 | } | |
156 | ||
157 | template <class T> | |
158 | void test_spots(T) | |
159 | { | |
160 | BOOST_MATH_STD_USING | |
161 | // | |
162 | // basic sanity checks, tolerance is 100 epsilon expressed as a percentage: | |
163 | // | |
164 | T tolerance = boost::math::tools::epsilon<T>() * 10000; | |
165 | BOOST_CHECK_CLOSE( | |
166 | ::boost::math::ibeta_inv( | |
167 | static_cast<T>(1), | |
168 | static_cast<T>(2), | |
169 | static_cast<T>(0.5)), | |
170 | static_cast<T>(0.29289321881345247559915563789515096071516406231153L), tolerance); | |
171 | BOOST_CHECK_CLOSE( | |
172 | ::boost::math::ibeta_inv( | |
173 | static_cast<T>(3), | |
174 | static_cast<T>(0.5), | |
175 | static_cast<T>(0.5)), | |
176 | static_cast<T>(0.92096723292382700385142816696980724853063433975470L), tolerance); | |
177 | BOOST_CHECK_CLOSE( | |
178 | ::boost::math::ibeta_inv( | |
179 | static_cast<T>(20.125), | |
180 | static_cast<T>(0.5), | |
181 | static_cast<T>(0.5)), | |
182 | static_cast<T>(0.98862133312917003480022776106012775747685870929920L), tolerance); | |
183 | BOOST_CHECK_CLOSE( | |
184 | ::boost::math::ibeta_inv( | |
185 | static_cast<T>(40), | |
186 | static_cast<T>(80), | |
187 | static_cast<T>(0.5)), | |
188 | static_cast<T>(0.33240456430025026300937492802591128972548660643778L), tolerance); | |
189 | BOOST_CHECK_CLOSE( | |
190 | ::boost::math::ibeta_inv( | |
191 | static_cast<T>(40), | |
192 | static_cast<T>(0.5), | |
193 | ldexp(T(1), -30)), | |
194 | static_cast<T>(0.624305407878048788716096298053941618358257550305573588792717L), tolerance); | |
195 | BOOST_CHECK_CLOSE( | |
196 | ::boost::math::ibeta_inv( | |
197 | static_cast<T>(40), | |
198 | static_cast<T>(0.5), | |
199 | static_cast<T>(1 - ldexp(T(1), -30))), | |
200 | static_cast<T>(0.99999999999999999998286262026583217516676792408012252456039L), tolerance); | |
201 | BOOST_CHECK_CLOSE( | |
202 | ::boost::math::ibeta_inv( | |
203 | static_cast<T>(0.5), | |
204 | static_cast<T>(40), | |
205 | static_cast<T>(ldexp(T(1), -30))), | |
206 | static_cast<T>(1.713737973416782483323207591987747543960774485649459249e-20L), tolerance); | |
207 | BOOST_CHECK_CLOSE( | |
208 | ::boost::math::ibeta_inv( | |
209 | static_cast<T>(0.5), | |
210 | static_cast<T>(0.75), | |
211 | static_cast<T>(ldexp(T(1), -30))), | |
212 | static_cast<T>(1.245132488513853853809715434621955746959615015005382639e-18L), tolerance); | |
213 | BOOST_CHECK_CLOSE( | |
214 | ::boost::math::ibeta_inv( | |
215 | static_cast<T>(0.5), | |
216 | static_cast<T>(0.5), | |
217 | static_cast<T>(0.25)), | |
218 | static_cast<T>(0.1464466094067262377995778189475754803575820311557629L), tolerance); | |
219 | BOOST_CHECK_CLOSE( | |
220 | ::boost::math::ibeta_inv( | |
221 | static_cast<T>(0.5), | |
222 | static_cast<T>(0.5), | |
223 | static_cast<T>(0.75)), | |
224 | static_cast<T>(0.853553390593273762200422181052424519642417968844237018294169L), tolerance); | |
225 | BOOST_CHECK_CLOSE( | |
226 | ::boost::math::ibeta_inv( | |
227 | static_cast<T>(1), | |
228 | static_cast<T>(5), | |
229 | static_cast<T>(0.125)), | |
230 | static_cast<T>(0.026352819384831863473794894078665766580641189002729204514544L), tolerance); | |
231 | BOOST_CHECK_CLOSE( | |
232 | ::boost::math::ibeta_inv( | |
233 | static_cast<T>(5), | |
234 | static_cast<T>(1), | |
235 | static_cast<T>(0.125)), | |
236 | static_cast<T>(0.659753955386447129687000985614820066516734506596709340752903L), tolerance); | |
237 | BOOST_CHECK_CLOSE( | |
238 | ::boost::math::ibeta_inv( | |
239 | static_cast<T>(1), | |
240 | static_cast<T>(0.125), | |
241 | static_cast<T>(0.125)), | |
242 | static_cast<T>(0.656391084194183349609374999999999999999999999999999999999999L), tolerance); | |
243 | BOOST_CHECK_CLOSE( | |
244 | ::boost::math::ibeta_inv( | |
245 | static_cast<T>(0.125), | |
246 | static_cast<T>(1), | |
247 | static_cast<T>(0.125)), | |
248 | static_cast<T>(5.960464477539062500000e-8), tolerance); | |
249 | BOOST_CHECK_CLOSE( | |
250 | ::boost::math::ibetac_inv( | |
251 | static_cast<T>(5), | |
252 | static_cast<T>(1), | |
253 | static_cast<T>(0.125)), | |
254 | static_cast<T>(0.973647180615168136526205105921334233419358810997270795485455L), tolerance); | |
255 | BOOST_CHECK_CLOSE( | |
256 | ::boost::math::ibetac_inv( | |
257 | static_cast<T>(1), | |
258 | static_cast<T>(5), | |
259 | static_cast<T>(0.125)), | |
260 | static_cast<T>(0.340246044613552870312999014385179933483265493403290659247096L), tolerance); | |
261 | BOOST_CHECK_CLOSE( | |
262 | ::boost::math::ibetac_inv( | |
263 | static_cast<T>(0.125), | |
264 | static_cast<T>(1), | |
265 | static_cast<T>(0.125)), | |
266 | static_cast<T>(0.343608915805816650390625000000000000000000000000000000000000L), tolerance); | |
267 | BOOST_CHECK_CLOSE( | |
268 | ::boost::math::ibetac_inv( | |
269 | static_cast<T>(1), | |
270 | static_cast<T>(0.125), | |
271 | static_cast<T>(0.125)), | |
272 | static_cast<T>(0.99999994039535522460937500000000000000000000000L), tolerance); | |
273 | } | |
274 |