23 #ifndef O2SCL_TEST_MGR_H
24 #define O2SCL_TEST_MGR_H
33 #include <boost/multiprecision/number.hpp>
36 #include <o2scl/string_conv.h>
37 #include <o2scl/misc.h>
38 #include <o2scl/table_units.h>
40 #include <gsl/gsl_vector.h>
41 #include <gsl/gsl_sys.h>
42 #include <gsl/gsl_matrix.h>
44 #ifndef DOXYGEN_NO_O2NS
54 #ifndef DOXYGEN_INTERNAL
63 void process_test(
bool ret, std::string d2, std::string description);
76 test_mgr(
bool success_l=
true, std::string last_fail_l=
"",
77 int ntests_l=0,
int output_level_l=1) {
130 template<
class data_t=
double>
131 bool test_abs(data_t result, data_t expected, data_t abs_error,
132 std::string description) {
134 if (std::isnan(expected)) {
135 ret=(std::isnan(expected)==std::isnan(result));
136 description=dtos<data_t>(result)+
" vs. "+ dtos<data_t>(expected)+
138 }
else if (std::isinf(expected)) {
139 ret=(std::isinf(expected)==std::isinf(result));
140 description=dtos<data_t>(result)+
" vs. "+ dtos<data_t>(expected)+
143 ret=(std::abs(expected-result)<abs_error);
145 description=dtos<data_t>(result)+
" vs. "+
146 dtos<data_t>(expected)+
" : "
147 +dtos<data_t>(std::abs(expected-result))+
" < "+
148 dtos<data_t>(abs_error)+
151 description=dtos<data_t>(result)+
" vs. "+
152 dtos<data_t>(expected)+
" : "
153 +dtos<data_t>(std::abs(expected-result))+
" > "+
154 dtos<data_t>(abs_error)+
167 template<
class data_t=
double>
168 bool test_rel(data_t result, data_t expected, data_t rel_error,
169 std::string description) {
171 if (std::isnan(expected)) {
172 ret=(std::isnan(expected)==std::isnan(result));
173 description=dtos<data_t>(result)+
" vs. "+ dtos<data_t>(expected)+
175 }
else if (std::isinf(expected)) {
176 ret=(std::isinf(expected)==std::isinf(result));
177 description=dtos<data_t>(result)+
" vs. "+ dtos<data_t>(expected)+
179 }
else if (expected==0.0) {
180 ret=test_abs<data_t>(result,expected,rel_error,description);
183 double obs_err=std::abs(expected-result)/std::abs(expected);
184 ret=(obs_err<rel_error);
186 description=dtos<data_t>(result)+
" vs. "+dtos<data_t>(expected)+
187 " : "+dtos<data_t>(obs_err)+
188 " < "+dtos<data_t>(rel_error)+
"\n "+description;
190 description=dtos<data_t>(result)+
" vs. "+dtos<data_t>(expected)+
191 " : "+dtos<data_t>(obs_err)+
192 " > "+dtos<data_t>(rel_error)+
"\n "+description;
200 #if defined(O2SCL_LD_TYPES) || defined(DOXYGEN)
208 template<
class data_t=
double>
210 std::string description) {
212 if (boost::multiprecision::isnan(expected)) {
213 ret=(boost::multiprecision::isnan(expected)==
214 boost::multiprecision::isnan(result));
215 description=dtos<data_t>(result)+
" vs. "+ dtos<data_t>(expected)+
217 }
else if (boost::multiprecision::isinf(expected)) {
218 ret=(boost::multiprecision::isinf(expected)==
219 boost::multiprecision::isinf(result));
220 description=dtos<data_t>(result)+
" vs. "+ dtos<data_t>(expected)+
223 ret=(boost::multiprecision::abs(expected-result)<abs_error);
225 description=dtos<data_t>(result)+
" vs. "+
226 dtos<data_t>(expected)+
" : "
227 +dtos<data_t>(boost::multiprecision::abs(expected-result))+
" < "+
228 dtos<data_t>(abs_error)+
231 description=dtos<data_t>(result)+
" vs. "+
232 dtos<data_t>(expected)+
" : "
233 +dtos<data_t>(boost::multiprecision::abs(expected-result))+
" > "+
234 dtos<data_t>(abs_error)+
250 template<
class data_t=
double>
252 std::string description) {
254 if (boost::multiprecision::isnan(expected)) {
255 ret=(boost::multiprecision::isnan(expected)==
256 boost::multiprecision::isnan(result));
257 description=dtos<data_t>(result)+
" vs. "+ dtos<data_t>(expected)+
259 }
else if (boost::multiprecision::isinf(expected)) {
260 ret=(boost::multiprecision::isinf(expected)==
261 boost::multiprecision::isinf(result));
262 description=dtos<data_t>(result)+
" vs. "+ dtos<data_t>(expected)+
264 }
else if (expected==0.0) {
265 ret=test_abs_boost<data_t>(result,expected,rel_error,description);
268 double obs_err=boost::multiprecision::abs(expected-result)/boost::multiprecision::abs(expected);
269 ret=(obs_err<rel_error);
271 description=dtos<data_t>(result)+
" vs. "+dtos<data_t>(expected)+
272 " : "+dtos<data_t>(obs_err)+
273 " < "+dtos<data_t>(rel_error)+
"\n "+description;
275 description=dtos<data_t>(result)+
" vs. "+dtos<data_t>(expected)+
276 " : "+dtos<data_t>(obs_err)+
277 " > "+dtos<data_t>(rel_error)+
"\n "+description;
290 template<
class data_t>
291 bool test_fact(data_t result, data_t expected, data_t factor,
292 std::string description) {
295 if (std::isnan(expected)) {
296 ret=(std::isnan(expected)==std::isnan(result));
297 }
else if (std::isinf(expected)) {
298 ret=(std::isinf(expected)==std::isinf(result));
300 ratio=expected/result;
301 ret=(ratio<factor && ratio>1.0/factor);
304 description=
dtos(result)+
" vs. "+
dtos(expected)+
"\n "+
312 bool test_str(std::string result, std::string expected,
313 std::string description);
316 bool test_gen(
bool value, std::string description);
325 template<
class vec_t,
class vec2_t,
class data_t>
327 data_t rel_error, std::string description) {
333 if (std::isnan(expected[i])) {
334 ret=(ret && (std::isnan(expected[i])==std::isnan(result[i])));
335 }
else if (std::isinf(expected[i])) {
336 ret=(ret && (std::isinf(expected[i])==std::isinf(result[i])));
337 }
else if (expected[i]==0.0) {
338 ret=(ret &&
test_abs(result[i],expected[i],rel_error,description));
339 if (std::abs(result[i]-expected[i])>max) {
340 max=std::abs(result[i]-expected[i]);
343 ret=(ret && ((std::abs(expected[i]-result[i]))/
344 std::abs(expected[i])<rel_error));
345 if (std::abs(expected[i]-result[i])/std::abs(expected[i])>max) {
346 max=std::abs(expected[i]-result[i])/std::abs(expected[i]);
351 description=((std::string)
"max=")+
o2scl::dtos(max)+
362 template<
class vec_t,
class vec2_t,
class data_t>
364 data_t abs_error, std::string description) {
369 if (std::isnan(expected[i])) {
370 ret=(ret && (std::isnan(expected[i])==std::isnan(result[i])));
371 }
else if (std::isinf(expected[i])) {
372 ret=(ret && (std::isinf(expected[i])==std::isinf(result[i])));
374 ret=(ret && (std::abs(expected[i]-result[i])<abs_error));
378 description=
"\n "+description;
387 template<
class vec_t,
class vec2_t,
class data_t>
389 data_t factor, std::string description) {
395 if (std::isnan(expected[i])) {
396 ret=(ret && (std::isnan(expected[i])==std::isnan(result[i])));
397 }
else if (std::isinf(expected[i])) {
398 ret=(ret && (std::isinf(expected[i])==std::isinf(result[i])));
400 ratio=expected[i]/result[i];
401 ret=(ret && (ratio<factor && ratio>1.0/factor));
405 description=
"\n "+description;
412 template<
class vec_t>
414 std::string description) {
419 ret=(ret && (result[i]==expected[i]));
422 description=
"\n "+description;
435 template<
class mat_t,
class mat2_t,
class data_t>
437 const mat2_t &expected,
438 data_t rel_error, std::string description) {
445 if (std::isnan(expected(i,j))) {
446 ret=(ret && (std::isnan(expected(i,j))==
447 std::isnan(result(i,j))));
448 }
else if (std::isinf(expected(i,j))) {
449 ret=(ret && (std::isinf(expected(i,j))==
450 std::isinf(result(i,j))));
451 }
else if (expected(i,j)==0.0) {
452 ret=(ret &&
test_abs(result(i,j),expected(i,j),rel_error,
454 if (std::abs(result(i,j)-expected(i,j))>max) {
455 max=std::abs(result(i,j)-expected(i,j));
458 ret=(ret && ((std::abs(expected(i,j)-result(i,j)))/
459 std::abs(expected(i,j))<rel_error));
460 if (std::abs(expected(i,j)-result(i,j))/
461 std::abs(expected(i,j))>max) {
462 max=std::abs(expected(i,j)-result(i,j))/
463 std::abs(expected(i,j));
469 description=((std::string)
"max=")+
o2scl::dtos(max)+
481 template<
class mat_t,
class mat2_t,
class data_t>
483 const mat2_t &expected,
484 data_t error, data_t zero_tol,
485 std::string description) {
492 if (std::isnan(expected(i,j))) {
493 ret=(ret && (std::isnan(expected(i,j))==
494 std::isnan(result(i,j))));
495 }
else if (std::isinf(expected(i,j))) {
496 ret=(ret && (std::isinf(expected(i,j))==
497 std::isinf(result(i,j))));
498 }
else if (expected(i,j)<zero_tol) {
499 ret=(ret &&
test_abs(result(i,j),expected(i,j),error,
501 if (std::abs(result(i,j)-expected(i,j))>max) {
502 max=std::abs(result(i,j)-expected(i,j));
505 ret=(ret && ((std::abs(expected(i,j)-result(i,j)))/
506 std::abs(expected(i,j))<error));
507 if (std::abs(expected(i,j)-result(i,j))/
508 std::abs(expected(i,j))>max) {
509 max=std::abs(expected(i,j)-result(i,j))/
510 std::abs(expected(i,j));
516 description=((std::string)
"max=")+
o2scl::dtos(max)+
527 template<
class mat_t,
class mat2_t,
class data_t>
529 const mat2_t &expected, data_t abs_error,
530 std::string description) {
538 if (std::isnan(expected(i,j))) {
539 ret=(ret && (std::isnan(expected(i,j))==
540 std::isnan(result(i,j))));
541 }
else if (std::isinf(expected(i,j))) {
542 ret=(ret && (std::isinf(expected(i,j))==
543 std::isinf(result(i,j))));
544 }
else if (expected(i,j)==0.0) {
545 ret=(ret &&
test_abs(result(i,j),expected(i,j),abs_error,
547 if (std::abs(result(i,j)-expected(i,j))>max) {
548 max=std::abs(result(i,j)-expected(i,j));
551 ret=(ret && ((std::abs(expected(i,j)-result(i,j)))<abs_error));
552 if (std::abs(expected(i,j)-result(i,j))>max) {
553 max=std::abs(expected(i,j)-result(i,j));
559 description=((std::string)
"max=")+
o2scl::dtos(max)+
577 template<
class vec_t,
class data_t>
580 data_t error, data_t zero_tol,
581 std::string description) {
584 int nr=result.get_nlines();
585 int nc=result.get_ncolumns();
586 std::vector<double> max(nc);
589 std::string col_name=expected.get_column_name(i);
591 std::string desc1=description+
" col: "+col_name+
" row: "+
593 if (std::isnan(expected.get(i,j))) {
594 bool ret1=
test_gen(std::isnan(expected.get(i,j))==
595 std::isnan(result.get(i,j)),desc1);
597 }
else if (std::isinf(expected.get(i,j))) {
598 bool ret1=
test_gen(std::isinf(expected.get(i,j))==
599 std::isinf(result.get(i,j)),desc1);
601 }
else if (expected.get(i,j)<zero_tol) {
602 bool ret1=
test_abs(result.get(i,j),expected.get(i,j),error,
604 if (std::abs(result.get(i,j)-expected.get(i,j))>max[i]) {
605 max[i]=std::abs(result.get(i,j)-expected.get(i,j));
608 bool ret1=
test_rel(result.get(i,j),expected.get(i,j),error,
611 if (std::abs(expected.get(i,j)-result.get(i,j))/
612 std::abs(expected.get(i,j))>max[i]) {
613 max[i]=std::abs(expected.get(i,j)-result.get(i,j))/
614 std::abs(expected.get(i,j));
622 std::cout <<
"Max for " << expected.get_column_name(i) <<
" is "
623 << max[i] << std::endl;
645 #ifndef DOXYGEN_NO_O2NS