23 #ifndef O2SCL_TEST_MGR_H
24 #define O2SCL_TEST_MGR_H
32 #include <o2scl/string_conv.h>
33 #include <o2scl/misc.h>
35 #include <gsl/gsl_vector.h>
36 #include <gsl/gsl_sys.h>
37 #include <gsl/gsl_matrix.h>
39 #ifndef DOXYGEN_NO_O2NS
52 #ifndef DOXYGEN_INTERNAL
61 void process_test(
bool ret, std::string d2, std::string description);
101 bool test_rel(
double result,
double expected,
double rel_error,
102 std::string description);
107 bool test_abs(
double result,
double expected,
double abs_error,
108 std::string description);
113 bool test_fact(
double result,
double expected,
double factor,
114 std::string description);
117 bool test_str(std::string result, std::string expected,
118 std::string description);
121 bool test_gen(
bool value, std::string description);
127 template<
class vec_t,
class vec2_t>
129 double rel_error, std::string description) {
139 }
else if (expected[i]==0.0) {
140 ret=(ret &&
test_abs(result[i],expected[i],rel_error,description));
141 if (fabs(result[i]-expected[i])>max) {
142 max=fabs(result[i]-expected[i]);
145 ret=(ret && ((fabs(expected[i]-result[i]))/
146 fabs(expected[i])<rel_error));
147 if (fabs(expected[i]-result[i])/fabs(expected[i])>max) {
148 max=fabs(expected[i]-result[i])/fabs(expected[i]);
153 description=((std::string)
"max=")+
o2scl::dtos(max)+
165 double rel_error, std::string description) {
177 }
else if (gsl_vector_get(expected,i)==0.0) {
178 ret=(ret &&
test_abs(gsl_vector_get(result,i),
179 gsl_vector_get(expected,i),
180 rel_error,description));
181 if (fabs(gsl_vector_get(result,i)-
182 gsl_vector_get(expected,i))>max) {
183 max=fabs(gsl_vector_get(result,i)-gsl_vector_get(expected,i));
186 ret=(ret && ((fabs(gsl_vector_get(expected,i)-
187 gsl_vector_get(result,i)))/
188 fabs(gsl_vector_get(expected,i))<rel_error));
189 if (fabs(gsl_vector_get(expected,i)-
190 gsl_vector_get(result,i))/
191 fabs(gsl_vector_get(expected,i))>max) {
192 max=fabs(gsl_vector_get(expected,i)-
193 gsl_vector_get(result,i))/
194 fabs(gsl_vector_get(expected,i));
199 description=((std::string)
"max=")+
o2scl::dtos(max)+
211 template<
class vec_t>
213 double rel_error, std::string description) {
225 }
else if (gsl_vector_get(expected,i)==0.0) {
226 ret=(ret &&
test_abs(result[i],gsl_vector_get(expected,i),
227 rel_error,description));
228 if (fabs(result[i]-gsl_vector_get(expected,i))>max) {
229 max=fabs(result[i]-gsl_vector_get(expected,i));
232 ret=(ret && ((fabs(gsl_vector_get(expected,i)-result[i]))/
233 fabs(gsl_vector_get(expected,i))<rel_error));
234 if (fabs(gsl_vector_get(expected,i)-result[i])/
235 fabs(gsl_vector_get(expected,i))>max) {
236 max=fabs(gsl_vector_get(expected,i)-result[i])/
237 fabs(gsl_vector_get(expected,i));
242 description=((std::string)
"max=")+
o2scl::dtos(max)+
255 gsl_matrix *expected,
256 double rel_error, std::string description) {
269 }
else if (gsl_matrix_get(expected,i,j)==0.0) {
270 ret=(ret &&
test_abs(gsl_matrix_get(result,i,j),
271 gsl_matrix_get(expected,i,j),rel_error,
273 if (fabs(gsl_matrix_get(result,i,j)-
274 gsl_matrix_get(expected,i,j))>max) {
275 max=fabs(gsl_matrix_get(result,i,j)-
276 gsl_matrix_get(expected,i,j));
279 ret=(ret && ((fabs(gsl_matrix_get(expected,i,j)-
280 gsl_matrix_get(result,i,j)))/
281 fabs(gsl_matrix_get(expected,i,j))<rel_error));
282 if (fabs(gsl_matrix_get(expected,i,j)-
283 gsl_matrix_get(result,i,j))/
284 fabs(gsl_matrix_get(expected,i,j))>max) {
285 max=fabs(gsl_matrix_get(expected,i,j)-
286 gsl_matrix_get(result,i,j))/
287 fabs(gsl_matrix_get(expected,i,j));
293 description=((std::string)
"max=")+
o2scl::dtos(max)+
305 template<
class mat_t>
307 gsl_matrix *expected,
308 double rel_error, std::string description) {
321 }
else if (gsl_matrix_get(expected,i,j)==0.0) {
323 gsl_matrix_get(expected,i,j),rel_error,
325 if (fabs(result(i,j)-gsl_matrix_get(expected,i,j))>max) {
326 max=fabs(result(i,j)-gsl_matrix_get(expected,i,j));
329 ret=(ret && ((fabs(gsl_matrix_get(expected,i,j)-result(i,j)))/
330 fabs(gsl_matrix_get(expected,i,j))<rel_error));
331 if (fabs(gsl_matrix_get(expected,i,j)-
332 result(i,j))/fabs(gsl_matrix_get(expected,i,j))>max) {
333 max=fabs(gsl_matrix_get(expected,i,j)-
334 result(i,j))/fabs(gsl_matrix_get(expected,i,j));
340 description=((std::string)
"max=")+
o2scl::dtos(max)+
352 template<
class mat_t>
354 gsl_matrix *expected,
355 double rel_error, std::string description) {
369 ret=(ret && ((fabs(gsl_matrix_get(expected,i,j)-result(i,j)))
371 if (fabs(gsl_matrix_get(expected,i,j)-result(i,j))>max) {
372 max=fabs(gsl_matrix_get(expected,i,j)-result(i,j));
378 description=((std::string)
"max=")+
o2scl::dtos(max)+
390 template<
class mat_t,
class mat2_t>
392 const mat2_t &expected,
393 double rel_error, std::string description) {
406 }
else if (expected(i,j)==0.0) {
407 ret=(ret &&
test_abs(result(i,j),expected(i,j),rel_error,
409 if (fabs(result(i,j)-expected(i,j))>max) {
410 max=fabs(result(i,j)-expected(i,j));
413 ret=(ret && ((fabs(expected(i,j)-result(i,j)))/
414 fabs(expected(i,j))<rel_error));
415 if (fabs(expected(i,j)-result(i,j))/fabs(expected(i,j))>max) {
416 max=fabs(expected(i,j)-result(i,j))/fabs(expected(i,j));
422 description=((std::string)
"max=")+
o2scl::dtos(max)+
434 template<
class vec_t,
class vec2_t>
436 double rel_error, std::string description) {
446 ret=(ret && (fabs(expected[i]-result[i])<rel_error));
450 description=
"\n "+description;
459 template<
class vec_t,
class vec2_t>
461 double factor, std::string description) {
472 ratio=expected[i]/result[i];
473 ret=(ret && (ratio<factor && ratio>1.0/factor));
477 description=
"\n "+description;
484 template<
class vec_t>
486 std::string description) {
491 ret=(ret && (result[i]==expected[i]));
494 description=
"\n "+description;
513 #ifndef DOXYGEN_NO_O2NS
bool success
True if all tests have passed.
int get_ntests()
Return the number of tests performed so far.
bool test_rel_matgslgsl(int nr, int nc, gsl_matrix *result, gsl_matrix *expected, double rel_error, std::string description)
Test for over each element of an array.
A class to manage testing and record success and failure.
bool test_fact(double result, double expected, double factor, std::string description)
Test for .
bool test_rel_matgsl(int nr, int nc, const mat_t &result, gsl_matrix *expected, double rel_error, std::string description)
Test for over each element of an array.
bool test_rel_mat(int nr, int nc, const mat_t &result, const mat2_t &expected, double rel_error, std::string description)
Test for over each element of an array.
bool test_fact_arr(int nv, const vec_t &result, const vec2_t &expected, double factor, std::string description)
Test for over each element of an array.
bool is_nan(double x)
Return true if x is not a number.
int ntests
The number of tests performed.
bool test_rel_arrgsl(int nv, const vec_t &result, gsl_vector *expected, double rel_error, std::string description)
Test for over each element of an array.
void process_test(bool ret, std::string d2, std::string description)
A helper function for processing tests.
std::string last_fail
The description of the last failed test.
std::string dtos(double x, int prec=6, bool auto_prec=false)
Convert a double to a string.
bool test_abs_arr(int nv, const vec_t &result, const vec2_t &expected, double rel_error, std::string description)
Test for over each element of an array.
bool test_str(std::string result, std::string expected, std::string description)
Test for .
bool test_abs_matgsl(int nr, int nc, const mat_t &result, gsl_matrix *expected, double rel_error, std::string description)
Test for over each element of an array.
bool test_abs(double result, double expected, double abs_error, std::string description)
Test for .
bool report()
Provide a report of all tests so far.
bool test_gen_arr(int nv, const vec_t &result, const vec_t &expected, std::string description)
Test for equality of a generic array.
bool test_rel_arrgslgsl(int nv, gsl_vector *result, gsl_vector *expected, double rel_error, std::string description)
Test for over each element of an array.
std::string get_last_fail()
Returns the description of the last test that failed.
void set_output_level(int l)
Set the output level.
bool test_rel_arr(int nv, const vec_t &result, const vec2_t &expected, double rel_error, std::string description)
Test for over each element of an array.
bool is_inf(double x)
Return true if x is infinite.
friend const test_mgr operator+(const test_mgr &left, const test_mgr &right)
Add two test_mgr objects (if either failed, the sum fails)
bool test_rel(double result, double expected, double rel_error, std::string description)
Test for .
bool test_gen(bool value, std::string description)
Test for .
int output_level
The output level.