Robotsystem of the Robot Group Leipzig: matrix.tests.hpp Source File

00001 /*************************************************************************** 00002 matrix.tests.cpp - description 00003 ------------------- 00004 email : georg.martius@web.de 00005 ***************************************************************************/ 00006 // 00007 // $Log: matrix.tests.hpp,v $ 00008 // Revision 1.3 2006/07/20 17:14:36 martius 00009 // removed std namespace from matrix.h 00010 // storable interface 00011 // abstract model and invertablemodel as superclasses for networks 00012 // 00013 // Revision 1.2 2006/07/14 12:24:01 martius 00014 // selforg becomes HEAD 00015 // 00016 // Revision 1.1.2.1 2006/07/10 12:01:02 martius 00017 // Matrixlib moved to selforg 00018 // 00019 // Revision 1.15 2005/10/21 11:58:25 martius 00020 // map2 (similar to map but for 2 matrices) 00021 // changed naming of functions to be more consistent. 00022 // Functions with "to" in front indicate the change of this. (Still not consistent with add, mult ...) 00023 // 00024 // Revision 1.14 2005/10/06 17:10:06 martius 00025 // convertToList 00026 // above and toAbove 00027 // 00028 // Revision 1.13 2005/08/06 20:47:36 martius 00029 // Commented 00030 // 00031 // Revision 1.12 2005/06/17 15:18:09 martius 00032 // 2x2 Inversion tested 00033 // 00034 // Revision 1.11 2005/06/10 08:21:50 martius 00035 // mult???wise are copy operations now! 00036 // toMult???wise are inplace instead 00037 // 00038 // Revision 1.10 2005/06/09 11:52:03 martius 00039 // multMT (M * M^T) and multTM (M^T * M) 00040 // 00041 // Revision 1.9 2005/06/02 22:49:33 martius 00042 // speed tests extended 00043 // 00044 // Revision 1.8 2005/06/02 08:49:26 martius 00045 // mult_row/column_wise 00046 // 00047 // Revision 1.7 2005/05/30 22:40:56 martius 00048 // map becomes toMap and the new map returns a new matrix 00049 // exp becomes toExp 00050 // 00051 // Revision 1.6 2005/05/30 21:46:54 martius 00052 // *** empty log message *** 00053 // 00054 // Revision 1.5 2005/05/30 16:43:02 martius 00055 // map function included (component-wise function application) 00056 // 00057 // Revision 1.4 2005/05/30 10:15:47 martius 00058 // proper log entry in header 00059 // 00060 /***************************************************************************/ 00061 00062 // ////////////////////////////////////////////////////////////////////////////// 00063 // ////////////////// UNIT TESTS for matrix lib ///////////////////////////////////////////////// 00064 00065 #ifdef UNITTEST 00066 #include "unit_test.hpp" 00067 using namespace matrix; 00068 using namespace std; 00069 00070 const D EPS=1e-9; 00071 bool comparetoidentity(const Matrix& m) { 00072 int worstdiagonal = 0; 00073 D maxunitydeviation = 0.0; 00074 D currentunitydeviation; 00075 for ( int i = 0; i < m.getM(); i++ ) { 00076 currentunitydeviation = m.val(i,i) - 1.; 00077 if ( currentunitydeviation < 0.0) currentunitydeviation *= -1.; 00078 if ( currentunitydeviation > maxunitydeviation ) { 00079 maxunitydeviation = currentunitydeviation; 00080 worstdiagonal = i; 00081 } 00082 } 00083 int worstoffdiagonalrow = 0; 00084 int worstoffdiagonalcolumn = 0; 00085 D maxzerodeviation = 0.0; 00086 D currentzerodeviation ; 00087 for ( int i = 0; i < m.getM(); i++ ) { 00088 for ( int j = 0; j < m.getN(); j++ ) { 00089 if ( i == j ) continue; // we look only at non-diagonal terms 00090 currentzerodeviation = m.val(i,j); 00091 if ( currentzerodeviation < 0.0) currentzerodeviation *= -1.0; 00092 if ( currentzerodeviation > maxzerodeviation ) { 00093 maxzerodeviation = currentzerodeviation; 00094 worstoffdiagonalrow = i; 00095 worstoffdiagonalcolumn = j; 00096 } 00097 00098 } 00099 } 00100 // cout << "\tWorst diagonal value deviation from unity: " 00101 // << maxunitydeviation << " at row/column " << worstdiagonal << endl; 00102 // cout << "\tWorst off-diagonal value deviation from zero: " 00103 // << maxzerodeviation << " at row = " << worstoffdiagonalrow 00104 // << ", column = " << worstoffdiagonalcolumn << endl; 00105 return (maxunitydeviation < EPS && maxzerodeviation < EPS); 00106 } 00107 00108 00109 UNIT_TEST_DEFINES 00110 00111 DEFINE_TEST( check_creation ) { 00112 cout << "\n -[ Creation and Access ]-\n"; 00113 Matrix M1(3,3); 00114 D testdata[9]={1,0,0, 0,1,0, 0,0,1}; 00115 Matrix M2(3,3); 00116 M2.set(testdata); 00117 M1.toId(); 00118 unit_assert( "id=identity_set", M1 == M2 ); 00119 D testdata2[12]={1,2,3, 4,5,6, 7,8,9, 10,11,12}; 00120 Matrix M3(4,3, testdata2); 00121 unit_assert( "dimension of matrix", M3.getM() == 4 && M3.getN() == 3 ); 00122 unit_assert( "field check", 00123 M3.val(0,2) == 3 && M3.val(2,1) == 8 && M3.val(3,0) == 10 ); 00124 Matrix M4(M3); 00125 unit_assert( "copy constructor", M3 == M4 ); 00126 Matrix M5(1,3,testdata2+3); 00127 unit_assert( "row", M3.row(1) == M5 ); 00128 00129 unit_pass(); 00130 } 00131 00132 DEFINE_TEST( check_vector_operation ) { 00133 cout << "\n -[ Vector Operations ]-\n"; 00134 D testdata[3]={-1,3,2}; 00135 const Matrix V1(1,3, testdata); 00136 const Matrix V2(3,1, testdata); 00137 Matrix V3(V1); 00138 00139 V3.toTranspose(); 00140 unit_assert( "transpose", V3 == V2 ); 00141 V3.toTranspose(); 00142 unit_assert( "double transpose", V3 == V1 ); 00143 D testdata2[3]={-2,6,4}; 00144 Matrix V4(1,3,testdata2); 00145 V3.add(V1,V1); 00146 unit_assert( "add", V3 == V4 ); 00147 D testdata3[3]={1,-3,-2}; 00148 Matrix V5(1,3,testdata3); 00149 V4.sub(V1,V3); 00150 unit_assert( "sub", V4 == V5 ); 00151 00152 D testdata4[3]={3,-9,-6}; 00153 V5.set(testdata4); 00154 V4.copy(V1); 00155 V4.toMult(-3.0); 00156 unit_assert( "mult with scalar I", V4 == V5 ); 00157 V4.mult(V1,-3.0); 00158 unit_assert( "mult with scalar II", V4 == V5 ); 00159 00160 double f = 14; 00161 Matrix V6(1,1,&f); 00162 V3.copy(V1); 00163 V3.toTranspose(); 00164 V4.mult(V1,V3); 00165 unit_assert( "scalarproduct", V4 == V6 ); 00166 V4.copy(V1); 00167 unit_assert( "scalarproduct with exp(2)", V4.multMT() == V6 ); 00168 V4.mult(V3,V1); 00169 unit_assert( "vector^T*vector=matrix", V4.getM() == 3); 00170 unit_assert( "vector^T*vector=exp(2)", V3.multMT().getM() == 3); 00171 00172 unit_pass(); 00173 } 00174 00175 DEFINE_TEST( check_matrix_operation ) { 00176 cout << "\n -[ Matrix Operations ]-\n"; 00177 D testdata[6]={1,2,3, 4,5,6 }; 00178 const Matrix M1(2,3, testdata); 00179 D testdata2[6]={1,4, 2,5, 3,6 }; 00180 const Matrix M2(3,2, testdata2); 00181 Matrix M3(M1); 00182 00183 M3.toTranspose(); 00184 unit_assert( "transpose", M3 == M2 ); 00185 D testdata3[6]={2,4,6, 8,10,12}; 00186 Matrix M4(2,3,testdata3); 00187 M3.add(M1,M1); 00188 unit_assert( "add", M3 == M4 ); 00189 00190 D testdata4[6]={-0.5, -1, -1.5, -2, -2.5, -3}; 00191 M3.set(testdata4); 00192 M4.copy(M1); 00193 M4.toMult(-0.5); 00194 unit_assert( "mult with scalar", M4 == M3 ); 00195 00196 D testdata5[6] = {2, -1, 0, 3, 2, -2}; 00197 D testdata6[4] = {8, -1, 20, -1}; 00198 Matrix M5 (3,2, testdata5); 00199 Matrix M6(2,2,testdata6); 00200 M4.mult(M1,M5); 00201 unit_assert( "mult(matrix, matrix)", M4 == M6 ); 00202 00203 M3.copy(M1); 00204 M4.copy(M1); 00205 M4.toExp(1); 00206 unit_assert( "exp(1)", M3 == M4 ); 00207 M3.toTranspose(); 00208 M4.toExp(T); 00209 unit_assert( "exp(T)=transpose", M3 == M4 ); 00210 M3.toId(); 00211 M4.toExp(0); 00212 unit_assert( "exp(0)=id", M3 == M4 ); 00213 00214 D testdata7[16] = {1,2,3,4, -4,2,1,3, 0.3,-0.9, 4, -3, 1,0.5,0.3,5.0}; 00215 Matrix M7(4,4,testdata7); 00216 Matrix M8(M7); 00217 M7.toExp(-1); 00218 M4.mult(M7,M8); 00219 unit_assert( "exp(-1)*exp(1)=id", comparetoidentity(M4) ); 00220 00221 D testdata9[6] = {sin(1.0),sin(2.0),sin(3.0), sin(4.0),sin(5.0),sin(6.0) }; 00222 Matrix M9(2,3,testdata9); 00223 M4.copy(M1); 00224 M4.toMap(sin); 00225 unit_assert( "map(sin)", M4 == M9 ); 00226 00227 D testdata10[6] = {2,4,6, -0.4,-0.5,-0.6 }; 00228 D testdata11[2] = {2,-0.1}; 00229 Matrix M10(2,3,testdata10); 00230 Matrix M11(2,1,testdata11); 00231 M4.copy(M1); 00232 M4.toMultrowwise(M11); 00233 unit_assert( "multrowwise()", M4 == M10 ); 00234 D testdata12[6] = {2,1,0, 8, 2.5, 0 }; 00235 D testdata13[3] = {2, 0.5, 0}; 00236 M10.set(2,3,testdata12); 00237 Matrix M12(3,1,testdata13); 00238 M4.copy(M1); 00239 M4.toMultcolwise(M12); 00240 unit_assert( "multrowwise()", M4 == M10 ); 00241 00242 M3.copy(M1); 00243 M4.copy(M1); 00244 M4.toTranspose(); 00245 M5 = M3.multMT(); 00246 M6.mult(M1,M4); 00247 unit_assert( "multMT() ", M5 == M6 ); 00248 M5 = M3.multTM(); 00249 M6.mult(M4,M1); 00250 unit_assert( "multTM() ", M5 == M6 ); 00251 00252 D testdata20[12]={1,2,3, 4,5,6, 1,2,3, 4,5,6 }; 00253 const Matrix M20(4,3, testdata20); 00254 const Matrix M21 = M1.above(M1); 00255 unit_assert( "above() ", M20 == M21 ); 00256 00257 00258 00259 unit_pass(); 00260 } 00261 00262 DEFINE_TEST( check_matrix_operators ) { 00263 cout << "\n -[ Matrix Operators (+ - * ^)]-\n"; 00264 D testdata[6]={1,2,3, 4,5,6 }; 00265 const Matrix M1(2,3, testdata); 00266 D testdata2[6]={1,4, 2,5, 3,6 }; 00267 const Matrix M2(3,2, testdata2); 00268 unit_assert( "^T ", (M1^T) == M2 ); 00269 D testdata3[6]={2,4,6, 8,10,12}; 00270 Matrix M4(2,3,testdata3); 00271 unit_assert( "+ ", M1+M1 == M4 ); 00272 unit_assert( "- ", M1+M1-M1 == M1 ); 00273 00274 D testdata4[6]={-0.5, -1, -1.5, -2, -2.5, -3}; 00275 Matrix M3(2,3, testdata4); 00276 unit_assert( "* scalar", M1*(-0.5) == M3 ); 00277 00278 D testdata5[6] = {2, -1, 0, 3, 2, -2}; 00279 D testdata6[4] = {8, -1, 20, -1}; 00280 Matrix M5 (3,2, testdata5); 00281 Matrix M6(2,2,testdata6); 00282 unit_assert( "* ", M1*M5 == M6 ); 00283 00284 unit_assert( "^1 ", (M1^1) == M1 ); 00285 M3.toId(); 00286 unit_assert( "^0=id ", (M1^0) == M3 ); 00287 00288 D testdata7[16] = {1,2,3,4, -4,2,1,3, 0.3,-0.9, 4, -3, 1,0.5,0.3,5.0}; 00289 Matrix M7(4,4,testdata7); 00290 unit_assert( "^1 * ^-1=id ", comparetoidentity(M7*(M7^-1)) ); 00291 unit_pass(); 00292 } 00293 00294 DEFINE_TEST( speed ) { 00295 cout << "\n -[ Speed: Inverion]-\n"; 00296 #ifndef NDEBUG 00297 cout << " DEBUG MODE! use -DNDEBUG -O3 (not -g) to get full performance\n"; 00298 #endif 00299 Matrix M1; 00300 srand(time(0)); 00301 D testdata0[9] = {1,2, -4,2}; 00302 Matrix M2(2,2,testdata0); 00303 UNIT_MEASURE_START("2x2 Matrix inversion", 100000) 00304 M1 = (M2^-1); 00305 UNIT_MEASURE_STOP(""); 00306 unit_assert( "validation", comparetoidentity(M1*M2)); 00307 /* LU version: 555648/s */ 00308 /* Explicit : 1428775/s */ 00309 D testdata1[9] = {1,2,3, -4,2,1, 0.3,-0.9}; 00310 Matrix M3(3,3,testdata1); 00311 UNIT_MEASURE_START("3x3 Matrix inversion", 100000) 00312 M1 = (M3^-1); 00313 UNIT_MEASURE_STOP(""); 00314 unit_assert( "validation", comparetoidentity(M1*M3)); 00315 00316 D testdata2[16] = {1,2,3,4, -4,2,1,3, 0.3,-0.9, 4, -3, 1,0.5,0.3,5.0}; 00317 Matrix M4(4,4,testdata2); 00318 UNIT_MEASURE_START("4x4 Matrix inversion", 100000) 00319 M1 = (M4^-1); 00320 UNIT_MEASURE_STOP(""); 00321 unit_assert( "validation", comparetoidentity(M1*M4)); 00322 00323 Matrix M20(20,20); 00324 for (int i=0; i < M20.getM(); i++) // define random values for initial matrix 00325 for (int j=0; j < M20.getN(); j++) { 00326 M20.val(i,j) = -22+(100. * rand())/RAND_MAX; 00327 } 00328 UNIT_MEASURE_START("20x20 Matrix inversion",1000) 00329 M1 = (M20^-1); 00330 UNIT_MEASURE_STOP(""); 00331 unit_assert( "validation", comparetoidentity(M1*M20)); 00332 00333 Matrix M200(200,200); 00334 rand(); // eliminates the first (= zero) call 00335 for (int i=0; i < M200.getM(); i++) // define random values for initial matrix 00336 for (int j=0; j < M200.getN(); j++) { 00337 M200.val(i,j) = -22+(100. * rand())/RAND_MAX; 00338 } 00339 UNIT_MEASURE_START("200x200 Matrix inversion",2) 00340 M1 = (M200^-1); 00341 UNIT_MEASURE_STOP(""); 00342 unit_assert( "validation", comparetoidentity(M1*M200)); 00343 00344 cout << "\n -[ Speed: Other Operations]-\n"; 00345 UNIT_MEASURE_START("20x20 Matrix multiplication with assignment",5000) 00346 M1 = M20*M20; 00347 UNIT_MEASURE_STOP(""); 00348 UNIT_MEASURE_START("20x20 Matrix addition with assignment",100000) 00349 M1= (M1 + M20); 00350 UNIT_MEASURE_STOP(""); 00351 UNIT_MEASURE_START("20x20 Matrix inplace addition",100000) 00352 M1 += M20; 00353 UNIT_MEASURE_STOP(""); 00354 UNIT_MEASURE_START("20x20 Matrix transposition",100000) 00355 M1 += M20; 00356 UNIT_MEASURE_STOP(""); 00357 00358 unit_pass(); 00359 } 00360 00361 UNIT_TEST_RUN( "Matrix Tests" ) 00362 ADD_TEST( check_creation ) 00363 ADD_TEST( check_vector_operation ) 00364 ADD_TEST( check_matrix_operation ) 00365 ADD_TEST( check_matrix_operators ) 00366 ADD_TEST( speed ) 00367 00368 UNIT_TEST_END 00369 00370 #endif // UNITTEST