[426] | 1 | #include "../bdm/math/square_mat.h" |
---|
| 2 | #include "../bdm/math/chmat.h" |
---|
| 3 | #include "mat_checks.h" |
---|
| 4 | #include "UnitTest++.h" |
---|
| 5 | #include "TestReporterStdout.h" |
---|
| 6 | #include <iostream> |
---|
| 7 | #include <iomanip> |
---|
| 8 | #include <stdlib.h> |
---|
| 9 | #include <string.h> |
---|
| 10 | |
---|
| 11 | using std::cout; |
---|
| 12 | using std::cerr; |
---|
| 13 | using std::endl; |
---|
| 14 | |
---|
| 15 | double epsilon = 0.00001; |
---|
| 16 | |
---|
| 17 | bool fast = false; |
---|
| 18 | |
---|
| 19 | namespace UnitTest |
---|
| 20 | { |
---|
| 21 | |
---|
| 22 | void CheckClose(TestResults &results, const itpp::mat &expected, |
---|
| 23 | const itpp::mat &actual, double tolerance, |
---|
| 24 | TestDetails const& details) { |
---|
| 25 | if (!AreClose(expected, actual, tolerance)) { |
---|
| 26 | MemoryOutStream stream; |
---|
| 27 | stream << "failed at " << expected.rows() |
---|
| 28 | << " x " << expected.cols(); |
---|
| 29 | |
---|
| 30 | results.OnTestFailure(details, stream.GetText()); |
---|
| 31 | } |
---|
| 32 | } |
---|
| 33 | |
---|
| 34 | } |
---|
| 35 | |
---|
| 36 | template<typename TMatrix> |
---|
| 37 | void test_until_overflow() { |
---|
| 38 | Real_Timer tt; |
---|
| 39 | int sz = 7; |
---|
| 40 | while (true) { |
---|
| 41 | mat A0 = randu(sz, sz); |
---|
| 42 | mat A = A0 * A0.T(); |
---|
| 43 | |
---|
| 44 | tt.tic(); |
---|
| 45 | TMatrix sqmat(A); |
---|
| 46 | double elapsed = tt.toc(); |
---|
| 47 | cout << "ctor(" << sz << " x " << sz << "): " << elapsed << " s" << endl; |
---|
| 48 | |
---|
| 49 | tt.tic(); |
---|
| 50 | mat res = sqmat.to_mat(); |
---|
| 51 | elapsed = tt.toc(); |
---|
| 52 | |
---|
| 53 | if (!fast) { |
---|
| 54 | CHECK_CLOSE(A, res, epsilon); |
---|
| 55 | } |
---|
| 56 | |
---|
| 57 | cout << "to_mat: " << elapsed << " s" << endl; |
---|
| 58 | |
---|
| 59 | vec v = randu(sz); |
---|
| 60 | double w = randu(); |
---|
| 61 | TMatrix sqmat2 = sqmat; |
---|
| 62 | |
---|
| 63 | tt.tic(); |
---|
| 64 | sqmat2.opupdt(v, w); |
---|
| 65 | elapsed = tt.toc(); |
---|
| 66 | |
---|
| 67 | if (!fast) { |
---|
| 68 | mat expA = A + w * outer_product(v, v); |
---|
| 69 | CHECK_CLOSE(expA, sqmat2.to_mat(), epsilon); |
---|
| 70 | } |
---|
| 71 | |
---|
| 72 | cout << "opupdt: " << elapsed << " s" << endl; |
---|
| 73 | |
---|
[438] | 74 | TMatrix invmat(sz); |
---|
| 75 | |
---|
| 76 | tt.tic(); |
---|
| 77 | sqmat.inv(invmat); |
---|
| 78 | elapsed = tt.toc(); |
---|
| 79 | |
---|
| 80 | if (!fast) { |
---|
| 81 | mat invA = inv(A); |
---|
| 82 | CHECK_CLOSE(invA, invmat.to_mat(), epsilon); |
---|
| 83 | } |
---|
| 84 | |
---|
| 85 | cout << "inv: " << elapsed << " s" << endl; |
---|
| 86 | |
---|
[426] | 87 | sz *= 7; |
---|
| 88 | } |
---|
| 89 | } |
---|
| 90 | |
---|
| 91 | SUITE(ldmat) { |
---|
| 92 | TEST(cycle) { |
---|
| 93 | test_until_overflow<ldmat>(); |
---|
| 94 | } |
---|
| 95 | } |
---|
| 96 | |
---|
| 97 | SUITE(fsqmat) { |
---|
| 98 | TEST(cycle) { |
---|
| 99 | test_until_overflow<fsqmat>(); |
---|
| 100 | } |
---|
| 101 | } |
---|
| 102 | |
---|
| 103 | SUITE(chmat) { |
---|
| 104 | TEST(cycle) { |
---|
| 105 | test_until_overflow<chmat>(); |
---|
| 106 | } |
---|
| 107 | } |
---|
| 108 | |
---|
| 109 | int main(int argc, char const *argv[]) { |
---|
| 110 | bool unknown = false; |
---|
| 111 | int update_next = 0; // 1 suite, 2 epsilon |
---|
| 112 | const char *suite = "ldmat"; |
---|
| 113 | const char **param = argv + 1; |
---|
| 114 | while (*param && !unknown) { |
---|
| 115 | if (update_next) { |
---|
| 116 | if (update_next == 1) { |
---|
| 117 | suite = *param; |
---|
| 118 | } else { |
---|
| 119 | double eps = atof(*param); |
---|
| 120 | if (eps > 0) { |
---|
| 121 | epsilon = eps; |
---|
| 122 | } else { |
---|
| 123 | cerr << "invalid epsilon value ignored" << endl; |
---|
| 124 | } |
---|
| 125 | } |
---|
| 126 | |
---|
| 127 | update_next = 0; |
---|
| 128 | } else { |
---|
| 129 | if (!strcmp(*param, "-c")) { |
---|
| 130 | update_next = 1; |
---|
| 131 | } else if (!strcmp(*param, "-e")) { |
---|
| 132 | update_next = 2; |
---|
| 133 | } else if (!strcmp(*param, "-f")) { |
---|
| 134 | fast = true; |
---|
| 135 | } else { |
---|
| 136 | unknown = true; |
---|
| 137 | } |
---|
| 138 | } |
---|
| 139 | |
---|
| 140 | ++param; |
---|
| 141 | } |
---|
| 142 | |
---|
| 143 | if (unknown || update_next) { |
---|
| 144 | cerr << "usage: " << argv[0] << " [ -f ] [ -e epsilon ] [ -c class ]" << endl; |
---|
| 145 | } else { |
---|
| 146 | UnitTest::TestReporterStdout reporter; |
---|
| 147 | UnitTest::TestRunner runner(reporter); |
---|
| 148 | return runner.RunTestsIf(UnitTest::Test::GetTestList(), |
---|
| 149 | suite, |
---|
| 150 | UnitTest::True(), |
---|
| 151 | 0); |
---|
| 152 | } |
---|
| 153 | } |
---|