All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
test_mgr.h
Go to the documentation of this file.
1 /*
2  -------------------------------------------------------------------
3 
4  Copyright (C) 2006-2014, Andrew W. Steiner
5 
6  This file is part of O2scl.
7 
8  O2scl is free software; you can redistribute it and/or modify
9  it under the terms of the GNU General Public License as published by
10  the Free Software Foundation; either version 3 of the License, or
11  (at your option) any later version.
12 
13  O2scl is distributed in the hope that it will be useful,
14  but WITHOUT ANY WARRANTY; without even the implied warranty of
15  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16  GNU General Public License for more details.
17 
18  You should have received a copy of the GNU General Public License
19  along with O2scl. If not, see <http://www.gnu.org/licenses/>.
20 
21  -------------------------------------------------------------------
22 */
23 #ifndef O2SCL_TEST_MGR_H
24 #define O2SCL_TEST_MGR_H
25 
26 /** \file test_mgr.h
27  \brief File defining \ref o2scl::test_mgr
28 */
29 
30 #include <string>
31 
32 #include <o2scl/string_conv.h>
33 #include <o2scl/misc.h>
34 
35 #include <gsl/gsl_vector.h>
36 #include <gsl/gsl_sys.h>
37 #include <gsl/gsl_matrix.h>
38 
39 #ifndef DOXYGEN_NO_O2NS
40 namespace o2scl {
41 #endif
42 
43  /** \brief A class to manage testing and record success and failure
44 
45  \future test_mgr::success and test_mgr::last_fail should be protected,
46  but that breaks the operator+() function. Can this be fixed?
47  */
48  class test_mgr {
49 
50  protected:
51 
52 #ifndef DOXYGEN_INTERNAL
53 
54  /// The number of tests performed
55  int ntests;
56 
57  /// The output level
59 
60  /// A helper function for processing tests
61  void process_test(bool ret, std::string d2, std::string description);
62 
63 #endif
64 
65  public:
66 
67  test_mgr() {
68  success=true;
69  ntests=0;
70  output_level=1;
71  }
72 
73  /** \brief Provide a report of all tests so far.
74 
75  Returns true if all tests have passed and false if at least
76  one test failed.
77  */
78  bool report();
79 
80  /// Returns the description of the last test that failed.
81  std::string get_last_fail() {return last_fail;};
82 
83  /** \brief Set the output level
84 
85  Possible values:
86  - 0 = No output
87  - 1 = Output only tests that fail
88  - 2 = Output all tests
89  */
90  void set_output_level(int l) { output_level=l; };
91 
92  /// Return the number of tests performed so far
93  int get_ntests() {return ntests;};
94 
95  /// \name The testing methods
96  ///@{
97 
98  /** \brief Test for \f$|\mathrm{result}-\mathrm{expected}|/
99  \mathrm{expected}<\mathrm{rel\_error}\f$
100  */
101  bool test_rel(double result, double expected, double rel_error,
102  std::string description);
103 
104  /** \brief Test for \f$|\mathrm{result}-\mathrm{expected}|/
105  <\mathrm{abs\_error}\f$
106  */
107  bool test_abs(double result, double expected, double abs_error,
108  std::string description);
109 
110  /** \brief Test for \f$1/\mathrm{factor} < \mathrm{result/expected}
111  < \mathrm{factor}\f$
112  */
113  bool test_fact(double result, double expected, double factor,
114  std::string description);
115 
116  /// Test for \f$\mathrm{result}=\mathrm{expected}\f$
117  bool test_str(std::string result, std::string expected,
118  std::string description);
119 
120  /// Test for \f$\mathrm{result}=\mathrm{expected}\f$
121  bool test_gen(bool value, std::string description);
122 
123  /** \brief Test for \f$|\mathrm{result}-\mathrm{expected}|/
124  \mathrm{expected}<\mathrm{rel\_error}\f$ over each element
125  of an array
126  */
127  template<class vec_t, class vec2_t>
128  bool test_rel_arr(int nv, const vec_t &result, const vec2_t &expected,
129  double rel_error, std::string description) {
130  bool ret=true;
131  double max=0.0;
132  int i;
133 
134  for(i=0;i<nv;i++) {
135  if (o2scl::is_nan(expected[i])) {
136  ret=(ret && (o2scl::is_nan(expected[i])==o2scl::is_nan(result[i])));
137  } else if (o2scl::is_inf(expected[i])) {
138  ret=(ret && (o2scl::is_inf(expected[i])==o2scl::is_inf(result[i])));
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]);
143  }
144  } else {
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]);
149  }
150  }
151  }
152 
153  description=((std::string)"max=")+o2scl::dtos(max)+
154  "\n "+description;
155  process_test(ret,"relative array",description);
156 
157  return ret;
158  }
159 
160  /** \brief Test for \f$|\mathrm{result}-\mathrm{expected}|/
161  \mathrm{expected}<\mathrm{rel\_error}\f$ over each element
162  of an array
163  */
164  bool test_rel_arrgslgsl(int nv, gsl_vector *result, gsl_vector *expected,
165  double rel_error, std::string description) {
166  bool ret=true;
167  double max=0.0;
168  int i;
169 
170  for(i=0;i<nv;i++) {
171  if (o2scl::is_nan(gsl_vector_get(expected,i))) {
172  ret=(ret && (o2scl::is_nan(gsl_vector_get(expected,i))==
173  o2scl::is_nan(gsl_vector_get(result,i))));
174  } else if (o2scl::is_inf(gsl_vector_get(expected,i))) {
175  ret=(ret && (o2scl::is_inf(gsl_vector_get(expected,i))==
176  o2scl::is_inf(gsl_vector_get(result,i))));
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));
184  }
185  } else {
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));
195  }
196  }
197  }
198 
199  description=((std::string)"max=")+o2scl::dtos(max)+
200  "\n "+description;
201  process_test(ret,"relative array",description);
202 
203  return ret;
204 
205  }
206 
207  /** \brief Test for \f$|\mathrm{result}-\mathrm{expected}|/
208  \mathrm{expected}<\mathrm{rel\_error}\f$ over each element
209  of an array
210  */
211  template<class vec_t>
212  bool test_rel_arrgsl(int nv, const vec_t &result, gsl_vector *expected,
213  double rel_error, std::string description) {
214  bool ret=true;
215  double max=0.0;
216  int i;
217 
218  for(i=0;i<nv;i++) {
219  if (o2scl::is_nan(gsl_vector_get(expected,i))) {
220  ret=(ret && (o2scl::is_nan(gsl_vector_get(expected,i))==
221  o2scl::is_nan(result[i])));
222  } else if (o2scl::is_inf(gsl_vector_get(expected,i))) {
223  ret=(ret && (o2scl::is_inf(gsl_vector_get(expected,i))==
224  o2scl::is_inf(result[i])));
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));
230  }
231  } else {
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));
238  }
239  }
240  }
241 
242  description=((std::string)"max=")+o2scl::dtos(max)+
243  "\n "+description;
244  process_test(ret,"relative array",description);
245 
246  return ret;
247 
248  }
249 
250  /** \brief Test for \f$|\mathrm{result}-\mathrm{expected}|/
251  \mathrm{expected}<\mathrm{rel\_error}\f$ over each element
252  of an array
253  */
254  bool test_rel_matgslgsl(int nr, int nc, gsl_matrix *result,
255  gsl_matrix *expected,
256  double rel_error, std::string description) {
257  bool ret=true;
258  double max=0.0;
259  int i, j;
260 
261  for(i=0;i<nr;i++) {
262  for(j=0;j<nc;j++) {
263  if (o2scl::is_nan(gsl_matrix_get(expected,i,j))) {
264  ret=(ret && (o2scl::is_nan(gsl_matrix_get(expected,i,j))==
265  o2scl::is_nan(gsl_matrix_get(result,i,j))));
266  } else if (o2scl::is_inf(gsl_matrix_get(expected,i,j))) {
267  ret=(ret && (o2scl::is_inf(gsl_matrix_get(expected,i,j))==
268  o2scl::is_inf(gsl_matrix_get(result,i,j))));
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,
272  description));
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));
277  }
278  } else {
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));
288  }
289  }
290  }
291  }
292 
293  description=((std::string)"max=")+o2scl::dtos(max)+
294  "\n "+description;
295  process_test(ret,"relative matrix",description);
296 
297  return ret;
298 
299  }
300 
301  /** \brief Test for \f$|\mathrm{result}-\mathrm{expected}|/
302  \mathrm{expected}<\mathrm{rel\_error}\f$ over each element
303  of an array
304  */
305  template<class mat_t>
306  bool test_rel_matgsl(int nr, int nc, const mat_t &result,
307  gsl_matrix *expected,
308  double rel_error, std::string description) {
309  bool ret=true;
310  double max=0.0;
311  int i, j;
312 
313  for(i=0;i<nr;i++) {
314  for(j=0;j<nc;j++) {
315  if (o2scl::is_nan(gsl_matrix_get(expected,i,j))) {
316  ret=(ret && (o2scl::is_nan(gsl_matrix_get(expected,i,j))==
317  o2scl::is_nan(result(i,j))));
318  } else if (o2scl::is_inf(gsl_matrix_get(expected,i,j))) {
319  ret=(ret && (o2scl::is_inf(gsl_matrix_get(expected,i,j))==
320  o2scl::is_inf(result(i,j))));
321  } else if (gsl_matrix_get(expected,i,j)==0.0) {
322  ret=(ret && test_abs(result(i,j),
323  gsl_matrix_get(expected,i,j),rel_error,
324  description));
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));
327  }
328  } else {
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));
335  }
336  }
337  }
338  }
339 
340  description=((std::string)"max=")+o2scl::dtos(max)+
341  "\n "+description;
342  process_test(ret,"relative matrix",description);
343 
344  return ret;
345 
346  }
347 
348  /** \brief Test for \f$|\mathrm{result}-\mathrm{expected}|
349  <\mathrm{rel\_error}\f$ over each element
350  of an array
351  */
352  template<class mat_t>
353  bool test_abs_matgsl(int nr, int nc, const mat_t &result,
354  gsl_matrix *expected,
355  double rel_error, std::string description) {
356  bool ret=true;
357  double max=0.0;
358  int i, j;
359 
360  for(i=0;i<nr;i++) {
361  for(j=0;j<nc;j++) {
362  if (o2scl::is_nan(gsl_matrix_get(expected,i,j))) {
363  ret=(ret && (o2scl::is_nan(gsl_matrix_get(expected,i,j))==
364  o2scl::is_nan(result(i,j))));
365  } else if (o2scl::is_inf(gsl_matrix_get(expected,i,j))) {
366  ret=(ret && (o2scl::is_inf(gsl_matrix_get(expected,i,j))==
367  o2scl::is_inf(result(i,j))));
368  } else {
369  ret=(ret && ((fabs(gsl_matrix_get(expected,i,j)-result(i,j)))
370  <rel_error));
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));
373  }
374  }
375  }
376  }
377 
378  description=((std::string)"max=")+o2scl::dtos(max)+
379  "\n "+description;
380  process_test(ret,"relative matrix",description);
381 
382  return ret;
383 
384  }
385 
386  /** \brief Test for \f$|\mathrm{result}-\mathrm{expected}|/
387  \mathrm{expected}<\mathrm{rel\_error}\f$ over each element
388  of an array
389  */
390  template<class mat_t, class mat2_t>
391  bool test_rel_mat(int nr, int nc, const mat_t &result,
392  const mat2_t &expected,
393  double rel_error, std::string description) {
394  bool ret=true;
395  double max=0.0;
396  int i, j;
397 
398  for(i=0;i<nr;i++) {
399  for(j=0;j<nc;j++) {
400  if (o2scl::is_nan(expected(i,j))) {
401  ret=(ret && (o2scl::is_nan(expected(i,j))==
402  o2scl::is_nan(result(i,j))));
403  } else if (o2scl::is_inf(expected(i,j))) {
404  ret=(ret && (o2scl::is_inf(expected(i,j))==
405  o2scl::is_inf(result(i,j))));
406  } else if (expected(i,j)==0.0) {
407  ret=(ret && test_abs(result(i,j),expected(i,j),rel_error,
408  description));
409  if (fabs(result(i,j)-expected(i,j))>max) {
410  max=fabs(result(i,j)-expected(i,j));
411  }
412  } else {
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));
417  }
418  }
419  }
420  }
421 
422  description=((std::string)"max=")+o2scl::dtos(max)+
423  "\n "+description;
424  process_test(ret,"relative matrix",description);
425 
426  return ret;
427 
428  }
429 
430  /** \brief Test for \f$|\mathrm{result}-\mathrm{expected}|/
431  <\mathrm{abs\_error}\f$ over each element
432  of an array
433  */
434  template<class vec_t, class vec2_t>
435  bool test_abs_arr(int nv, const vec_t &result, const vec2_t &expected,
436  double rel_error, std::string description) {
437  bool ret=true;
438  int i;
439 
440  for(i=0;i<nv;i++) {
441  if (o2scl::is_nan(expected[i])) {
442  ret=(ret && (o2scl::is_nan(expected[i])==o2scl::is_nan(result[i])));
443  } else if (o2scl::is_inf(expected[i])) {
444  ret=(ret && (o2scl::is_inf(expected[i])==o2scl::is_inf(result[i])));
445  } else {
446  ret=(ret && (fabs(expected[i]-result[i])<rel_error));
447  }
448  }
449 
450  description="\n "+description;
451  process_test(ret,"absolute array",description);
452 
453  return ret;
454  }
455 
456  /** \brief Test for \f$ 1/factor < result/expected < factor \f$
457  over each element of an array
458  */
459  template<class vec_t, class vec2_t>
460  bool test_fact_arr(int nv, const vec_t &result, const vec2_t &expected,
461  double factor, std::string description) {
462  bool ret=true;
463  int i;
464  double ratio;
465 
466  for(i=0;i<nv;i++) {
467  if (o2scl::is_nan(expected[i])) {
468  ret=(ret && (o2scl::is_nan(expected[i])==o2scl::is_nan(result[i])));
469  } else if (o2scl::is_inf(expected[i])) {
470  ret=(ret && (o2scl::is_inf(expected[i])==o2scl::is_inf(result[i])));
471  } else {
472  ratio=expected[i]/result[i];
473  ret=(ret && (ratio<factor && ratio>1.0/factor));
474  }
475  }
476 
477  description="\n "+description;
478  process_test(ret,"factor array",description);
479 
480  return ret;
481  }
482 
483  /// Test for equality of a generic array
484  template<class vec_t>
485  bool test_gen_arr(int nv, const vec_t &result, const vec_t &expected,
486  std::string description) {
487  bool ret=true;
488  int i;
489 
490  for(i=0;i<nv;i++) {
491  ret=(ret && (result[i]==expected[i]));
492  }
493 
494  description="\n "+description;
495  process_test(ret,"generic array",description);
496 
497  return ret;
498  }
499  ///@}
500 
501  /// Add two test_mgr objects (if either failed, the sum fails)
502  friend const test_mgr operator+(const test_mgr& left,
503  const test_mgr& right);
504 
505  /// True if all tests have passed
506  bool success;
507 
508  /// The description of the last failed test
509  std::string last_fail;
510 
511  };
512 
513 #ifndef DOXYGEN_NO_O2NS
514 }
515 #endif
516 
517 #endif
bool success
True if all tests have passed.
Definition: test_mgr.h:506
int get_ntests()
Return the number of tests performed so far.
Definition: test_mgr.h:93
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.
Definition: test_mgr.h:254
A class to manage testing and record success and failure.
Definition: test_mgr.h:48
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.
Definition: test_mgr.h:306
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.
Definition: test_mgr.h:391
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.
Definition: test_mgr.h:460
bool is_nan(double x)
Return true if x is not a number.
int ntests
The number of tests performed.
Definition: test_mgr.h:55
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.
Definition: test_mgr.h:212
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.
Definition: test_mgr.h:509
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.
Definition: test_mgr.h:435
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.
Definition: test_mgr.h:353
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.
Definition: test_mgr.h:485
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.
Definition: test_mgr.h:164
std::string get_last_fail()
Returns the description of the last test that failed.
Definition: test_mgr.h:81
void set_output_level(int l)
Set the output level.
Definition: test_mgr.h:90
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.
Definition: test_mgr.h:128
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.
Definition: test_mgr.h:58

Documentation generated with Doxygen. Provided under the GNU Free Documentation License (see License Information).
Hosted at Get Object-oriented Scientific Computing
Lib at SourceForge.net. Fast, secure and Free Open Source software
downloads..