Line data Source code
1 : // Copyright (c) 2018-2022 Made to Order Software Corp. All Rights Reserved
2 : //
3 : // https://snapwebsites.org/
4 : // contact@m2osw.com
5 : //
6 : // This program is free software; you can redistribute it and/or modify
7 : // it under the terms of the GNU General Public License as published by
8 : // the Free Software Foundation; either version 2 of the License, or
9 : // (at your option) any later version.
10 : //
11 : // This program is distributed in the hope that it will be useful,
12 : // but WITHOUT ANY WARRANTY; without even the implied warranty of
13 : // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 : // GNU General Public License for more details.
15 : //
16 : // You should have received a copy of the GNU General Public License along
17 : // with this program; if not, write to the Free Software Foundation, Inc.,
18 : // 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
19 :
20 :
21 : /** \file
22 : * \brief Verify the various matrix implementation.
23 : *
24 : * This file implements tests to verify that the matrix templates
25 : * work as expected.
26 : */
27 :
28 : // self
29 : //
30 : #include "catch_main.h"
31 :
32 :
33 :
34 : // ignore the == and != against float warnings
35 : //
36 : #pragma GCC diagnostic ignored "-Wfloat-equal"
37 :
38 :
39 : // snapdev lib
40 : //
41 : #include <snapdev/matrix.h>
42 :
43 :
44 : // last include
45 : //
46 : #include <snapdev/poison.h>
47 :
48 :
49 :
50 :
51 :
52 : namespace
53 : {
54 700 : auto frand = []()
55 : {
56 700 : return static_cast<double>(rand()) / static_cast<double>(rand());
57 700 : };
58 : }
59 : // no name namespace
60 :
61 :
62 :
63 6 : CATCH_TEST_CASE("matrix_init", "[matrix]")
64 : {
65 : // constructor/copy
66 : // and
67 : // zero/identity
68 : //
69 8 : CATCH_GIVEN("constructor")
70 : {
71 8 : CATCH_START_SECTION("matrix: empty")
72 : {
73 2 : snapdev::matrix<double> empty;
74 :
75 1 : CATCH_REQUIRE(empty.rows() == 0);
76 1 : CATCH_REQUIRE(empty.columns() == 0);
77 :
78 2 : snapdev::matrix<double> copy(empty);
79 :
80 1 : CATCH_REQUIRE(empty.rows() == 0);
81 1 : CATCH_REQUIRE(empty.columns() == 0);
82 : }
83 : CATCH_END_SECTION()
84 :
85 8 : CATCH_START_SECTION("matrix: 2x2")
86 : {
87 2 : snapdev::matrix<double> m(2, 2);
88 :
89 1 : CATCH_REQUIRE(m.rows() == 2);
90 1 : CATCH_REQUIRE(m.columns() == 2);
91 :
92 : // by default we get an identity
93 : //
94 1 : CATCH_REQUIRE(m[0][0] == 1.0);
95 1 : CATCH_REQUIRE(m[0][1] == 0.0);
96 1 : CATCH_REQUIRE(m[1][0] == 0.0);
97 1 : CATCH_REQUIRE(m[1][1] == 1.0);
98 :
99 1 : double const r00(frand());
100 1 : double const r01(frand());
101 1 : double const r10(frand());
102 1 : double const r11(frand());
103 :
104 1 : m[0][0] = r00;
105 1 : m[0][1] = r01;
106 1 : m[1][0] = r10;
107 1 : m[1][1] = r11;
108 :
109 1 : CATCH_REQUIRE(m[0][0] == r00);
110 1 : CATCH_REQUIRE(m[0][1] == r01);
111 1 : CATCH_REQUIRE(m[1][0] == r10);
112 1 : CATCH_REQUIRE(m[1][1] == r11);
113 :
114 2 : snapdev::matrix<double> copy(m);
115 2 : snapdev::matrix<double> c2;
116 :
117 1 : CATCH_REQUIRE(copy[0][0] == r00);
118 1 : CATCH_REQUIRE(copy[0][1] == r01);
119 1 : CATCH_REQUIRE(copy[1][0] == r10);
120 1 : CATCH_REQUIRE(copy[1][1] == r11);
121 :
122 1 : m.clear();
123 :
124 1 : CATCH_REQUIRE(m[0][0] == 0.0);
125 1 : CATCH_REQUIRE(m[0][1] == 0.0);
126 1 : CATCH_REQUIRE(m[1][0] == 0.0);
127 1 : CATCH_REQUIRE(m[1][1] == 0.0);
128 :
129 : // copy is still intact
130 : //
131 1 : CATCH_REQUIRE(copy[0][0] == r00);
132 1 : CATCH_REQUIRE(copy[0][1] == r01);
133 1 : CATCH_REQUIRE(copy[1][0] == r10);
134 1 : CATCH_REQUIRE(copy[1][1] == r11);
135 :
136 1 : CATCH_REQUIRE(c2.rows() == 0);
137 1 : CATCH_REQUIRE(c2.columns() == 0);
138 :
139 1 : c2 = copy;
140 :
141 1 : CATCH_REQUIRE(c2[0][0] == r00);
142 1 : CATCH_REQUIRE(c2[0][1] == r01);
143 1 : CATCH_REQUIRE(c2[1][0] == r10);
144 1 : CATCH_REQUIRE(c2[1][1] == r11);
145 :
146 1 : c2.swap(m);
147 :
148 1 : CATCH_REQUIRE(m[0][0] == r00);
149 1 : CATCH_REQUIRE(m[0][1] == r01);
150 1 : CATCH_REQUIRE(m[1][0] == r10);
151 1 : CATCH_REQUIRE(m[1][1] == r11);
152 :
153 1 : CATCH_REQUIRE(c2[0][0] == 0.0);
154 1 : CATCH_REQUIRE(c2[0][1] == 0.0);
155 1 : CATCH_REQUIRE(c2[1][0] == 0.0);
156 1 : CATCH_REQUIRE(c2[1][1] == 0.0);
157 : }
158 : CATCH_END_SECTION()
159 :
160 8 : CATCH_START_SECTION("matrix: 3x3")
161 : {
162 2 : snapdev::matrix<double> m(3, 3);
163 :
164 1 : CATCH_REQUIRE(m.rows() == 3);
165 1 : CATCH_REQUIRE(m.columns() == 3);
166 :
167 : // by default we get an identity
168 : //
169 1 : CATCH_REQUIRE(m[0][0] == 1.0);
170 1 : CATCH_REQUIRE(m[0][1] == 0.0);
171 1 : CATCH_REQUIRE(m[0][2] == 0.0);
172 1 : CATCH_REQUIRE(m[1][0] == 0.0);
173 1 : CATCH_REQUIRE(m[1][1] == 1.0);
174 1 : CATCH_REQUIRE(m[1][2] == 0.0);
175 1 : CATCH_REQUIRE(m[2][0] == 0.0);
176 1 : CATCH_REQUIRE(m[2][1] == 0.0);
177 1 : CATCH_REQUIRE(m[2][2] == 1.0);
178 :
179 1 : double const r00(frand());
180 1 : double const r01(frand());
181 1 : double const r02(frand());
182 1 : double const r10(frand());
183 1 : double const r11(frand());
184 1 : double const r12(frand());
185 1 : double const r20(frand());
186 1 : double const r21(frand());
187 1 : double const r22(frand());
188 :
189 1 : m[0][0] = r00;
190 1 : m[0][1] = r01;
191 1 : m[0][2] = r02;
192 1 : m[1][0] = r10;
193 1 : m[1][1] = r11;
194 1 : m[1][2] = r12;
195 1 : m[2][0] = r20;
196 1 : m[2][1] = r21;
197 1 : m[2][2] = r22;
198 :
199 1 : CATCH_REQUIRE(m[0][0] == r00);
200 1 : CATCH_REQUIRE(m[0][1] == r01);
201 1 : CATCH_REQUIRE(m[0][2] == r02);
202 1 : CATCH_REQUIRE(m[1][0] == r10);
203 1 : CATCH_REQUIRE(m[1][1] == r11);
204 1 : CATCH_REQUIRE(m[1][2] == r12);
205 1 : CATCH_REQUIRE(m[2][0] == r20);
206 1 : CATCH_REQUIRE(m[2][1] == r21);
207 1 : CATCH_REQUIRE(m[2][2] == r22);
208 :
209 2 : snapdev::matrix<double> copy(m);
210 2 : snapdev::matrix<double> c2;
211 :
212 1 : CATCH_REQUIRE(copy[0][0] == r00);
213 1 : CATCH_REQUIRE(copy[0][1] == r01);
214 1 : CATCH_REQUIRE(copy[0][2] == r02);
215 1 : CATCH_REQUIRE(copy[1][0] == r10);
216 1 : CATCH_REQUIRE(copy[1][1] == r11);
217 1 : CATCH_REQUIRE(copy[1][2] == r12);
218 1 : CATCH_REQUIRE(copy[2][0] == r20);
219 1 : CATCH_REQUIRE(copy[2][1] == r21);
220 1 : CATCH_REQUIRE(copy[2][2] == r22);
221 :
222 1 : m.clear();
223 :
224 1 : CATCH_REQUIRE(m[0][0] == 0.0);
225 1 : CATCH_REQUIRE(m[0][1] == 0.0);
226 1 : CATCH_REQUIRE(m[0][2] == 0.0);
227 1 : CATCH_REQUIRE(m[1][0] == 0.0);
228 1 : CATCH_REQUIRE(m[1][1] == 0.0);
229 1 : CATCH_REQUIRE(m[1][2] == 0.0);
230 1 : CATCH_REQUIRE(m[2][0] == 0.0);
231 1 : CATCH_REQUIRE(m[2][1] == 0.0);
232 1 : CATCH_REQUIRE(m[2][2] == 0.0);
233 :
234 : // copy is still intact
235 : //
236 1 : CATCH_REQUIRE(copy[0][0] == r00);
237 1 : CATCH_REQUIRE(copy[0][1] == r01);
238 1 : CATCH_REQUIRE(copy[0][2] == r02);
239 1 : CATCH_REQUIRE(copy[1][0] == r10);
240 1 : CATCH_REQUIRE(copy[1][1] == r11);
241 1 : CATCH_REQUIRE(copy[1][2] == r12);
242 1 : CATCH_REQUIRE(copy[2][0] == r20);
243 1 : CATCH_REQUIRE(copy[2][1] == r21);
244 1 : CATCH_REQUIRE(copy[2][2] == r22);
245 :
246 1 : CATCH_REQUIRE(c2.rows() == 0);
247 1 : CATCH_REQUIRE(c2.columns() == 0);
248 :
249 1 : c2 = copy;
250 :
251 1 : CATCH_REQUIRE(c2[0][0] == r00);
252 1 : CATCH_REQUIRE(c2[0][1] == r01);
253 1 : CATCH_REQUIRE(c2[0][2] == r02);
254 1 : CATCH_REQUIRE(c2[1][0] == r10);
255 1 : CATCH_REQUIRE(c2[1][1] == r11);
256 1 : CATCH_REQUIRE(c2[1][2] == r12);
257 1 : CATCH_REQUIRE(c2[2][0] == r20);
258 1 : CATCH_REQUIRE(c2[2][1] == r21);
259 1 : CATCH_REQUIRE(c2[2][2] == r22);
260 :
261 1 : c2.swap(m);
262 :
263 1 : CATCH_REQUIRE(m[0][0] == r00);
264 1 : CATCH_REQUIRE(m[0][1] == r01);
265 1 : CATCH_REQUIRE(m[0][2] == r02);
266 1 : CATCH_REQUIRE(m[1][0] == r10);
267 1 : CATCH_REQUIRE(m[1][1] == r11);
268 1 : CATCH_REQUIRE(m[1][2] == r12);
269 1 : CATCH_REQUIRE(m[2][0] == r20);
270 1 : CATCH_REQUIRE(m[2][1] == r21);
271 1 : CATCH_REQUIRE(m[2][2] == r22);
272 :
273 1 : CATCH_REQUIRE(c2[0][0] == 0.0);
274 1 : CATCH_REQUIRE(c2[0][1] == 0.0);
275 1 : CATCH_REQUIRE(c2[0][2] == 0.0);
276 1 : CATCH_REQUIRE(c2[1][0] == 0.0);
277 1 : CATCH_REQUIRE(c2[1][1] == 0.0);
278 1 : CATCH_REQUIRE(c2[1][2] == 0.0);
279 1 : CATCH_REQUIRE(c2[2][0] == 0.0);
280 1 : CATCH_REQUIRE(c2[2][1] == 0.0);
281 1 : CATCH_REQUIRE(c2[2][2] == 0.0);
282 : }
283 : CATCH_END_SECTION()
284 :
285 8 : CATCH_START_SECTION("matrix: 4x4")
286 : {
287 2 : snapdev::matrix<double> m(4, 4);
288 :
289 1 : CATCH_REQUIRE(m.rows() == 4);
290 1 : CATCH_REQUIRE(m.columns() == 4);
291 :
292 : // by default we get an identity
293 : //
294 1 : CATCH_REQUIRE(m[0][0] == 1.0);
295 1 : CATCH_REQUIRE(m[0][1] == 0.0);
296 1 : CATCH_REQUIRE(m[0][2] == 0.0);
297 1 : CATCH_REQUIRE(m[0][3] == 0.0);
298 1 : CATCH_REQUIRE(m[1][0] == 0.0);
299 1 : CATCH_REQUIRE(m[1][1] == 1.0);
300 1 : CATCH_REQUIRE(m[1][2] == 0.0);
301 1 : CATCH_REQUIRE(m[1][3] == 0.0);
302 1 : CATCH_REQUIRE(m[2][0] == 0.0);
303 1 : CATCH_REQUIRE(m[2][1] == 0.0);
304 1 : CATCH_REQUIRE(m[2][2] == 1.0);
305 1 : CATCH_REQUIRE(m[2][3] == 0.0);
306 1 : CATCH_REQUIRE(m[3][0] == 0.0);
307 1 : CATCH_REQUIRE(m[3][1] == 0.0);
308 1 : CATCH_REQUIRE(m[3][2] == 0.0);
309 1 : CATCH_REQUIRE(m[3][3] == 1.0);
310 :
311 1 : double const r00(frand());
312 1 : double const r01(frand());
313 1 : double const r02(frand());
314 1 : double const r03(frand());
315 1 : double const r10(frand());
316 1 : double const r11(frand());
317 1 : double const r12(frand());
318 1 : double const r13(frand());
319 1 : double const r20(frand());
320 1 : double const r21(frand());
321 1 : double const r22(frand());
322 1 : double const r23(frand());
323 1 : double const r30(frand());
324 1 : double const r31(frand());
325 1 : double const r32(frand());
326 1 : double const r33(frand());
327 :
328 1 : m[0][0] = r00;
329 1 : m[0][1] = r01;
330 1 : m[0][2] = r02;
331 1 : m[0][3] = r03;
332 1 : m[1][0] = r10;
333 1 : m[1][1] = r11;
334 1 : m[1][2] = r12;
335 1 : m[1][3] = r13;
336 1 : m[2][0] = r20;
337 1 : m[2][1] = r21;
338 1 : m[2][2] = r22;
339 1 : m[2][3] = r23;
340 1 : m[3][0] = r30;
341 1 : m[3][1] = r31;
342 1 : m[3][2] = r32;
343 1 : m[3][3] = r33;
344 :
345 1 : CATCH_REQUIRE(m[0][0] == r00);
346 1 : CATCH_REQUIRE(m[0][1] == r01);
347 1 : CATCH_REQUIRE(m[0][2] == r02);
348 1 : CATCH_REQUIRE(m[0][3] == r03);
349 1 : CATCH_REQUIRE(m[1][0] == r10);
350 1 : CATCH_REQUIRE(m[1][1] == r11);
351 1 : CATCH_REQUIRE(m[1][2] == r12);
352 1 : CATCH_REQUIRE(m[1][3] == r13);
353 1 : CATCH_REQUIRE(m[2][0] == r20);
354 1 : CATCH_REQUIRE(m[2][1] == r21);
355 1 : CATCH_REQUIRE(m[2][2] == r22);
356 1 : CATCH_REQUIRE(m[2][3] == r23);
357 1 : CATCH_REQUIRE(m[3][0] == r30);
358 1 : CATCH_REQUIRE(m[3][1] == r31);
359 1 : CATCH_REQUIRE(m[3][2] == r32);
360 1 : CATCH_REQUIRE(m[3][3] == r33);
361 :
362 2 : snapdev::matrix<double> copy(m);
363 2 : snapdev::matrix<double> c2;
364 :
365 1 : CATCH_REQUIRE(copy[0][0] == r00);
366 1 : CATCH_REQUIRE(copy[0][1] == r01);
367 1 : CATCH_REQUIRE(copy[0][2] == r02);
368 1 : CATCH_REQUIRE(copy[0][3] == r03);
369 1 : CATCH_REQUIRE(copy[1][0] == r10);
370 1 : CATCH_REQUIRE(copy[1][1] == r11);
371 1 : CATCH_REQUIRE(copy[1][2] == r12);
372 1 : CATCH_REQUIRE(copy[1][3] == r13);
373 1 : CATCH_REQUIRE(copy[2][0] == r20);
374 1 : CATCH_REQUIRE(copy[2][1] == r21);
375 1 : CATCH_REQUIRE(copy[2][2] == r22);
376 1 : CATCH_REQUIRE(copy[2][3] == r23);
377 1 : CATCH_REQUIRE(copy[3][0] == r30);
378 1 : CATCH_REQUIRE(copy[3][1] == r31);
379 1 : CATCH_REQUIRE(copy[3][2] == r32);
380 1 : CATCH_REQUIRE(copy[3][3] == r33);
381 :
382 1 : m.clear();
383 :
384 1 : CATCH_REQUIRE(m[0][0] == 0.0);
385 1 : CATCH_REQUIRE(m[0][1] == 0.0);
386 1 : CATCH_REQUIRE(m[0][2] == 0.0);
387 1 : CATCH_REQUIRE(m[0][3] == 0.0);
388 1 : CATCH_REQUIRE(m[1][0] == 0.0);
389 1 : CATCH_REQUIRE(m[1][1] == 0.0);
390 1 : CATCH_REQUIRE(m[1][2] == 0.0);
391 1 : CATCH_REQUIRE(m[1][3] == 0.0);
392 1 : CATCH_REQUIRE(m[2][0] == 0.0);
393 1 : CATCH_REQUIRE(m[2][1] == 0.0);
394 1 : CATCH_REQUIRE(m[2][2] == 0.0);
395 1 : CATCH_REQUIRE(m[2][3] == 0.0);
396 :
397 : // copy is still intact
398 : //
399 1 : CATCH_REQUIRE(copy[0][0] == r00);
400 1 : CATCH_REQUIRE(copy[0][1] == r01);
401 1 : CATCH_REQUIRE(copy[0][2] == r02);
402 1 : CATCH_REQUIRE(copy[0][3] == r03);
403 1 : CATCH_REQUIRE(copy[1][0] == r10);
404 1 : CATCH_REQUIRE(copy[1][1] == r11);
405 1 : CATCH_REQUIRE(copy[1][2] == r12);
406 1 : CATCH_REQUIRE(copy[1][3] == r13);
407 1 : CATCH_REQUIRE(copy[2][0] == r20);
408 1 : CATCH_REQUIRE(copy[2][1] == r21);
409 1 : CATCH_REQUIRE(copy[2][2] == r22);
410 1 : CATCH_REQUIRE(copy[2][3] == r23);
411 1 : CATCH_REQUIRE(copy[3][0] == r30);
412 1 : CATCH_REQUIRE(copy[3][1] == r31);
413 1 : CATCH_REQUIRE(copy[3][2] == r32);
414 1 : CATCH_REQUIRE(copy[3][3] == r33);
415 :
416 1 : CATCH_REQUIRE(c2.rows() == 0);
417 1 : CATCH_REQUIRE(c2.columns() == 0);
418 :
419 1 : c2 = copy;
420 :
421 1 : CATCH_REQUIRE(c2[0][0] == r00);
422 1 : CATCH_REQUIRE(c2[0][1] == r01);
423 1 : CATCH_REQUIRE(c2[0][2] == r02);
424 1 : CATCH_REQUIRE(c2[0][3] == r03);
425 1 : CATCH_REQUIRE(c2[1][0] == r10);
426 1 : CATCH_REQUIRE(c2[1][1] == r11);
427 1 : CATCH_REQUIRE(c2[1][2] == r12);
428 1 : CATCH_REQUIRE(c2[1][3] == r13);
429 1 : CATCH_REQUIRE(c2[2][0] == r20);
430 1 : CATCH_REQUIRE(c2[2][1] == r21);
431 1 : CATCH_REQUIRE(c2[2][2] == r22);
432 1 : CATCH_REQUIRE(c2[2][3] == r23);
433 1 : CATCH_REQUIRE(c2[3][0] == r30);
434 1 : CATCH_REQUIRE(c2[3][1] == r31);
435 1 : CATCH_REQUIRE(c2[3][2] == r32);
436 1 : CATCH_REQUIRE(c2[3][3] == r33);
437 :
438 : //std::cout << c2 << std::endl;
439 :
440 1 : c2.swap(m);
441 :
442 1 : CATCH_REQUIRE(m[0][0] == r00);
443 1 : CATCH_REQUIRE(m[0][1] == r01);
444 1 : CATCH_REQUIRE(m[0][2] == r02);
445 1 : CATCH_REQUIRE(m[0][3] == r03);
446 1 : CATCH_REQUIRE(m[1][0] == r10);
447 1 : CATCH_REQUIRE(m[1][1] == r11);
448 1 : CATCH_REQUIRE(m[1][2] == r12);
449 1 : CATCH_REQUIRE(m[1][3] == r13);
450 1 : CATCH_REQUIRE(m[2][0] == r20);
451 1 : CATCH_REQUIRE(m[2][1] == r21);
452 1 : CATCH_REQUIRE(m[2][2] == r22);
453 1 : CATCH_REQUIRE(m[2][3] == r23);
454 1 : CATCH_REQUIRE(m[3][0] == r30);
455 1 : CATCH_REQUIRE(m[3][1] == r31);
456 1 : CATCH_REQUIRE(m[3][2] == r32);
457 1 : CATCH_REQUIRE(m[3][3] == r33);
458 :
459 1 : CATCH_REQUIRE(c2[0][0] == 0.0);
460 1 : CATCH_REQUIRE(c2[0][1] == 0.0);
461 1 : CATCH_REQUIRE(c2[0][2] == 0.0);
462 1 : CATCH_REQUIRE(c2[0][3] == 0.0);
463 1 : CATCH_REQUIRE(c2[1][0] == 0.0);
464 1 : CATCH_REQUIRE(c2[1][1] == 0.0);
465 1 : CATCH_REQUIRE(c2[1][2] == 0.0);
466 1 : CATCH_REQUIRE(c2[1][3] == 0.0);
467 1 : CATCH_REQUIRE(c2[2][0] == 0.0);
468 1 : CATCH_REQUIRE(c2[2][1] == 0.0);
469 1 : CATCH_REQUIRE(c2[2][2] == 0.0);
470 1 : CATCH_REQUIRE(c2[2][3] == 0.0);
471 : }
472 : CATCH_END_SECTION()
473 : }
474 4 : }
475 :
476 :
477 10 : CATCH_TEST_CASE("matrix_additive", "[matrix]")
478 : {
479 : // create two random 4x4 matrices and make sure the add works
480 : //
481 16 : CATCH_GIVEN("add")
482 : {
483 8 : CATCH_START_SECTION("matrix: a+=<scalar>")
484 : {
485 : // setup A
486 : //
487 2 : snapdev::matrix<double> a(4, 4);
488 :
489 1 : CATCH_REQUIRE(a.rows() == 4);
490 1 : CATCH_REQUIRE(a.columns() == 4);
491 :
492 1 : double const a00(frand());
493 1 : double const a01(frand());
494 1 : double const a02(frand());
495 1 : double const a03(frand());
496 1 : double const a10(frand());
497 1 : double const a11(frand());
498 1 : double const a12(frand());
499 1 : double const a13(frand());
500 1 : double const a20(frand());
501 1 : double const a21(frand());
502 1 : double const a22(frand());
503 1 : double const a23(frand());
504 1 : double const a30(frand());
505 1 : double const a31(frand());
506 1 : double const a32(frand());
507 1 : double const a33(frand());
508 :
509 1 : a[0][0] = a00;
510 1 : a[0][1] = a01;
511 1 : a[0][2] = a02;
512 1 : a[0][3] = a03;
513 1 : a[1][0] = a10;
514 1 : a[1][1] = a11;
515 1 : a[1][2] = a12;
516 1 : a[1][3] = a13;
517 1 : a[2][0] = a20;
518 1 : a[2][1] = a21;
519 1 : a[2][2] = a22;
520 1 : a[2][3] = a23;
521 1 : a[3][0] = a30;
522 1 : a[3][1] = a31;
523 1 : a[3][2] = a32;
524 1 : a[3][3] = a33;
525 :
526 1 : CATCH_REQUIRE(a[0][0] == a00);
527 1 : CATCH_REQUIRE(a[0][1] == a01);
528 1 : CATCH_REQUIRE(a[0][2] == a02);
529 1 : CATCH_REQUIRE(a[0][3] == a03);
530 1 : CATCH_REQUIRE(a[1][0] == a10);
531 1 : CATCH_REQUIRE(a[1][1] == a11);
532 1 : CATCH_REQUIRE(a[1][2] == a12);
533 1 : CATCH_REQUIRE(a[1][3] == a13);
534 1 : CATCH_REQUIRE(a[2][0] == a20);
535 1 : CATCH_REQUIRE(a[2][1] == a21);
536 1 : CATCH_REQUIRE(a[2][2] == a22);
537 1 : CATCH_REQUIRE(a[2][3] == a23);
538 1 : CATCH_REQUIRE(a[3][0] == a30);
539 1 : CATCH_REQUIRE(a[3][1] == a31);
540 1 : CATCH_REQUIRE(a[3][2] == a32);
541 1 : CATCH_REQUIRE(a[3][3] == a33);
542 :
543 : // create a scalar for our test
544 : //
545 1 : double scalar(frand());
546 1 : while(scalar == 0.0)
547 : {
548 0 : scalar = frand();
549 : }
550 :
551 : // run operation A += <scalar>
552 : //
553 1 : a += scalar;
554 :
555 : // this can't fail because we ensure scalar != 0
556 : //
557 1 : CATCH_REQUIRE_FALSE(a[0][0] == a00);
558 1 : CATCH_REQUIRE_FALSE(a[0][1] == a01);
559 1 : CATCH_REQUIRE_FALSE(a[0][2] == a02);
560 1 : CATCH_REQUIRE_FALSE(a[0][3] == a03);
561 1 : CATCH_REQUIRE_FALSE(a[1][0] == a10);
562 1 : CATCH_REQUIRE_FALSE(a[1][1] == a11);
563 1 : CATCH_REQUIRE_FALSE(a[1][2] == a12);
564 1 : CATCH_REQUIRE_FALSE(a[1][3] == a13);
565 1 : CATCH_REQUIRE_FALSE(a[2][0] == a20);
566 1 : CATCH_REQUIRE_FALSE(a[2][1] == a21);
567 1 : CATCH_REQUIRE_FALSE(a[2][2] == a22);
568 1 : CATCH_REQUIRE_FALSE(a[2][3] == a23);
569 1 : CATCH_REQUIRE_FALSE(a[3][0] == a30);
570 1 : CATCH_REQUIRE_FALSE(a[3][1] == a31);
571 1 : CATCH_REQUIRE_FALSE(a[3][2] == a32);
572 1 : CATCH_REQUIRE_FALSE(a[3][3] == a33);
573 :
574 1 : CATCH_REQUIRE(a[0][0] == a00 + scalar);
575 1 : CATCH_REQUIRE(a[0][1] == a01 + scalar);
576 1 : CATCH_REQUIRE(a[0][2] == a02 + scalar);
577 1 : CATCH_REQUIRE(a[0][3] == a03 + scalar);
578 1 : CATCH_REQUIRE(a[1][0] == a10 + scalar);
579 1 : CATCH_REQUIRE(a[1][1] == a11 + scalar);
580 1 : CATCH_REQUIRE(a[1][2] == a12 + scalar);
581 1 : CATCH_REQUIRE(a[1][3] == a13 + scalar);
582 1 : CATCH_REQUIRE(a[2][0] == a20 + scalar);
583 1 : CATCH_REQUIRE(a[2][1] == a21 + scalar);
584 1 : CATCH_REQUIRE(a[2][2] == a22 + scalar);
585 1 : CATCH_REQUIRE(a[2][3] == a23 + scalar);
586 1 : CATCH_REQUIRE(a[3][0] == a30 + scalar);
587 1 : CATCH_REQUIRE(a[3][1] == a31 + scalar);
588 1 : CATCH_REQUIRE(a[3][2] == a32 + scalar);
589 1 : CATCH_REQUIRE(a[3][3] == a33 + scalar);
590 : }
591 : CATCH_END_SECTION()
592 :
593 8 : CATCH_START_SECTION("matrix: b=a+<scalar>")
594 : {
595 : // setup A
596 : //
597 2 : snapdev::matrix<double> a(4, 4);
598 :
599 1 : CATCH_REQUIRE(a.rows() == 4);
600 1 : CATCH_REQUIRE(a.columns() == 4);
601 :
602 1 : double const a00(frand());
603 1 : double const a01(frand());
604 1 : double const a02(frand());
605 1 : double const a03(frand());
606 1 : double const a10(frand());
607 1 : double const a11(frand());
608 1 : double const a12(frand());
609 1 : double const a13(frand());
610 1 : double const a20(frand());
611 1 : double const a21(frand());
612 1 : double const a22(frand());
613 1 : double const a23(frand());
614 1 : double const a30(frand());
615 1 : double const a31(frand());
616 1 : double const a32(frand());
617 1 : double const a33(frand());
618 :
619 1 : a[0][0] = a00;
620 1 : a[0][1] = a01;
621 1 : a[0][2] = a02;
622 1 : a[0][3] = a03;
623 1 : a[1][0] = a10;
624 1 : a[1][1] = a11;
625 1 : a[1][2] = a12;
626 1 : a[1][3] = a13;
627 1 : a[2][0] = a20;
628 1 : a[2][1] = a21;
629 1 : a[2][2] = a22;
630 1 : a[2][3] = a23;
631 1 : a[3][0] = a30;
632 1 : a[3][1] = a31;
633 1 : a[3][2] = a32;
634 1 : a[3][3] = a33;
635 :
636 1 : CATCH_REQUIRE(a[0][0] == a00);
637 1 : CATCH_REQUIRE(a[0][1] == a01);
638 1 : CATCH_REQUIRE(a[0][2] == a02);
639 1 : CATCH_REQUIRE(a[0][3] == a03);
640 1 : CATCH_REQUIRE(a[1][0] == a10);
641 1 : CATCH_REQUIRE(a[1][1] == a11);
642 1 : CATCH_REQUIRE(a[1][2] == a12);
643 1 : CATCH_REQUIRE(a[1][3] == a13);
644 1 : CATCH_REQUIRE(a[2][0] == a20);
645 1 : CATCH_REQUIRE(a[2][1] == a21);
646 1 : CATCH_REQUIRE(a[2][2] == a22);
647 1 : CATCH_REQUIRE(a[2][3] == a23);
648 1 : CATCH_REQUIRE(a[3][0] == a30);
649 1 : CATCH_REQUIRE(a[3][1] == a31);
650 1 : CATCH_REQUIRE(a[3][2] == a32);
651 1 : CATCH_REQUIRE(a[3][3] == a33);
652 :
653 : // setup B
654 : //
655 2 : snapdev::matrix<double> b(4, 4);
656 :
657 1 : CATCH_REQUIRE(b.rows() == 4);
658 1 : CATCH_REQUIRE(b.columns() == 4);
659 :
660 1 : double const b00(frand());
661 1 : double const b01(frand());
662 1 : double const b02(frand());
663 1 : double const b03(frand());
664 1 : double const b10(frand());
665 1 : double const b11(frand());
666 1 : double const b12(frand());
667 1 : double const b13(frand());
668 1 : double const b20(frand());
669 1 : double const b21(frand());
670 1 : double const b22(frand());
671 1 : double const b23(frand());
672 1 : double const b30(frand());
673 1 : double const b31(frand());
674 1 : double const b32(frand());
675 1 : double const b33(frand());
676 :
677 1 : b[0][0] = b00;
678 1 : b[0][1] = b01;
679 1 : b[0][2] = b02;
680 1 : b[0][3] = b03;
681 1 : b[1][0] = b10;
682 1 : b[1][1] = b11;
683 1 : b[1][2] = b12;
684 1 : b[1][3] = b13;
685 1 : b[2][0] = b20;
686 1 : b[2][1] = b21;
687 1 : b[2][2] = b22;
688 1 : b[2][3] = b23;
689 1 : b[3][0] = b30;
690 1 : b[3][1] = b31;
691 1 : b[3][2] = b32;
692 1 : b[3][3] = b33;
693 :
694 1 : CATCH_REQUIRE(b[0][0] == b00);
695 1 : CATCH_REQUIRE(b[0][1] == b01);
696 1 : CATCH_REQUIRE(b[0][2] == b02);
697 1 : CATCH_REQUIRE(b[0][3] == b03);
698 1 : CATCH_REQUIRE(b[1][0] == b10);
699 1 : CATCH_REQUIRE(b[1][1] == b11);
700 1 : CATCH_REQUIRE(b[1][2] == b12);
701 1 : CATCH_REQUIRE(b[1][3] == b13);
702 1 : CATCH_REQUIRE(b[2][0] == b20);
703 1 : CATCH_REQUIRE(b[2][1] == b21);
704 1 : CATCH_REQUIRE(b[2][2] == b22);
705 1 : CATCH_REQUIRE(b[2][3] == b23);
706 1 : CATCH_REQUIRE(b[3][0] == b30);
707 1 : CATCH_REQUIRE(b[3][1] == b31);
708 1 : CATCH_REQUIRE(b[3][2] == b32);
709 1 : CATCH_REQUIRE(b[3][3] == b33);
710 :
711 : // create a scalar for our test
712 : //
713 1 : double scalar(frand());
714 1 : while(scalar == 0.0)
715 : {
716 0 : scalar = frand();
717 : }
718 :
719 : // run operation B = A + <scalar>
720 : //
721 1 : b = a + scalar;
722 :
723 1 : CATCH_REQUIRE(a[0][0] == a00);
724 1 : CATCH_REQUIRE(a[0][1] == a01);
725 1 : CATCH_REQUIRE(a[0][2] == a02);
726 1 : CATCH_REQUIRE(a[0][3] == a03);
727 1 : CATCH_REQUIRE(a[1][0] == a10);
728 1 : CATCH_REQUIRE(a[1][1] == a11);
729 1 : CATCH_REQUIRE(a[1][2] == a12);
730 1 : CATCH_REQUIRE(a[1][3] == a13);
731 1 : CATCH_REQUIRE(a[2][0] == a20);
732 1 : CATCH_REQUIRE(a[2][1] == a21);
733 1 : CATCH_REQUIRE(a[2][2] == a22);
734 1 : CATCH_REQUIRE(a[2][3] == a23);
735 1 : CATCH_REQUIRE(a[3][0] == a30);
736 1 : CATCH_REQUIRE(a[3][1] == a31);
737 1 : CATCH_REQUIRE(a[3][2] == a32);
738 1 : CATCH_REQUIRE(a[3][3] == a33);
739 :
740 : // this can't fail because we ensure scalar != 0
741 : //
742 1 : CATCH_REQUIRE_FALSE(b[0][0] == b00);
743 1 : CATCH_REQUIRE_FALSE(b[0][1] == b01);
744 1 : CATCH_REQUIRE_FALSE(b[0][2] == b02);
745 1 : CATCH_REQUIRE_FALSE(b[0][3] == b03);
746 1 : CATCH_REQUIRE_FALSE(b[1][0] == b10);
747 1 : CATCH_REQUIRE_FALSE(b[1][1] == b11);
748 1 : CATCH_REQUIRE_FALSE(b[1][2] == b12);
749 1 : CATCH_REQUIRE_FALSE(b[1][3] == b13);
750 1 : CATCH_REQUIRE_FALSE(b[2][0] == b20);
751 1 : CATCH_REQUIRE_FALSE(b[2][1] == b21);
752 1 : CATCH_REQUIRE_FALSE(b[2][2] == b22);
753 1 : CATCH_REQUIRE_FALSE(b[2][3] == b23);
754 1 : CATCH_REQUIRE_FALSE(b[3][0] == b30);
755 1 : CATCH_REQUIRE_FALSE(b[3][1] == b31);
756 1 : CATCH_REQUIRE_FALSE(b[3][2] == b32);
757 1 : CATCH_REQUIRE_FALSE(b[3][3] == b33);
758 :
759 1 : CATCH_REQUIRE(b[0][0] == a00 + scalar);
760 1 : CATCH_REQUIRE(b[0][1] == a01 + scalar);
761 1 : CATCH_REQUIRE(b[0][2] == a02 + scalar);
762 1 : CATCH_REQUIRE(b[0][3] == a03 + scalar);
763 1 : CATCH_REQUIRE(b[1][0] == a10 + scalar);
764 1 : CATCH_REQUIRE(b[1][1] == a11 + scalar);
765 1 : CATCH_REQUIRE(b[1][2] == a12 + scalar);
766 1 : CATCH_REQUIRE(b[1][3] == a13 + scalar);
767 1 : CATCH_REQUIRE(b[2][0] == a20 + scalar);
768 1 : CATCH_REQUIRE(b[2][1] == a21 + scalar);
769 1 : CATCH_REQUIRE(b[2][2] == a22 + scalar);
770 1 : CATCH_REQUIRE(b[2][3] == a23 + scalar);
771 1 : CATCH_REQUIRE(b[3][0] == a30 + scalar);
772 1 : CATCH_REQUIRE(b[3][1] == a31 + scalar);
773 1 : CATCH_REQUIRE(b[3][2] == a32 + scalar);
774 1 : CATCH_REQUIRE(b[3][3] == a33 + scalar);
775 : }
776 : CATCH_END_SECTION()
777 :
778 8 : CATCH_START_SECTION("matrix: c=a+b")
779 : {
780 : // setup A
781 : //
782 2 : snapdev::matrix<double> a(4, 4);
783 :
784 1 : CATCH_REQUIRE(a.rows() == 4);
785 1 : CATCH_REQUIRE(a.columns() == 4);
786 :
787 1 : double const a00(frand());
788 1 : double const a01(frand());
789 1 : double const a02(frand());
790 1 : double const a03(frand());
791 1 : double const a10(frand());
792 1 : double const a11(frand());
793 1 : double const a12(frand());
794 1 : double const a13(frand());
795 1 : double const a20(frand());
796 1 : double const a21(frand());
797 1 : double const a22(frand());
798 1 : double const a23(frand());
799 1 : double const a30(frand());
800 1 : double const a31(frand());
801 1 : double const a32(frand());
802 1 : double const a33(frand());
803 :
804 1 : a[0][0] = a00;
805 1 : a[0][1] = a01;
806 1 : a[0][2] = a02;
807 1 : a[0][3] = a03;
808 1 : a[1][0] = a10;
809 1 : a[1][1] = a11;
810 1 : a[1][2] = a12;
811 1 : a[1][3] = a13;
812 1 : a[2][0] = a20;
813 1 : a[2][1] = a21;
814 1 : a[2][2] = a22;
815 1 : a[2][3] = a23;
816 1 : a[3][0] = a30;
817 1 : a[3][1] = a31;
818 1 : a[3][2] = a32;
819 1 : a[3][3] = a33;
820 :
821 1 : CATCH_REQUIRE(a[0][0] == a00);
822 1 : CATCH_REQUIRE(a[0][1] == a01);
823 1 : CATCH_REQUIRE(a[0][2] == a02);
824 1 : CATCH_REQUIRE(a[0][3] == a03);
825 1 : CATCH_REQUIRE(a[1][0] == a10);
826 1 : CATCH_REQUIRE(a[1][1] == a11);
827 1 : CATCH_REQUIRE(a[1][2] == a12);
828 1 : CATCH_REQUIRE(a[1][3] == a13);
829 1 : CATCH_REQUIRE(a[2][0] == a20);
830 1 : CATCH_REQUIRE(a[2][1] == a21);
831 1 : CATCH_REQUIRE(a[2][2] == a22);
832 1 : CATCH_REQUIRE(a[2][3] == a23);
833 1 : CATCH_REQUIRE(a[3][0] == a30);
834 1 : CATCH_REQUIRE(a[3][1] == a31);
835 1 : CATCH_REQUIRE(a[3][2] == a32);
836 1 : CATCH_REQUIRE(a[3][3] == a33);
837 :
838 : // setup B
839 : //
840 2 : snapdev::matrix<double> b(4, 4);
841 :
842 1 : CATCH_REQUIRE(b.rows() == 4);
843 1 : CATCH_REQUIRE(b.columns() == 4);
844 :
845 1 : double const b00(frand());
846 1 : double const b01(frand());
847 1 : double const b02(frand());
848 1 : double const b03(frand());
849 1 : double const b10(frand());
850 1 : double const b11(frand());
851 1 : double const b12(frand());
852 1 : double const b13(frand());
853 1 : double const b20(frand());
854 1 : double const b21(frand());
855 1 : double const b22(frand());
856 1 : double const b23(frand());
857 1 : double const b30(frand());
858 1 : double const b31(frand());
859 1 : double const b32(frand());
860 1 : double const b33(frand());
861 :
862 1 : b[0][0] = b00;
863 1 : b[0][1] = b01;
864 1 : b[0][2] = b02;
865 1 : b[0][3] = b03;
866 1 : b[1][0] = b10;
867 1 : b[1][1] = b11;
868 1 : b[1][2] = b12;
869 1 : b[1][3] = b13;
870 1 : b[2][0] = b20;
871 1 : b[2][1] = b21;
872 1 : b[2][2] = b22;
873 1 : b[2][3] = b23;
874 1 : b[3][0] = b30;
875 1 : b[3][1] = b31;
876 1 : b[3][2] = b32;
877 1 : b[3][3] = b33;
878 :
879 1 : CATCH_REQUIRE(b[0][0] == b00);
880 1 : CATCH_REQUIRE(b[0][1] == b01);
881 1 : CATCH_REQUIRE(b[0][2] == b02);
882 1 : CATCH_REQUIRE(b[0][3] == b03);
883 1 : CATCH_REQUIRE(b[1][0] == b10);
884 1 : CATCH_REQUIRE(b[1][1] == b11);
885 1 : CATCH_REQUIRE(b[1][2] == b12);
886 1 : CATCH_REQUIRE(b[1][3] == b13);
887 1 : CATCH_REQUIRE(b[2][0] == b20);
888 1 : CATCH_REQUIRE(b[2][1] == b21);
889 1 : CATCH_REQUIRE(b[2][2] == b22);
890 1 : CATCH_REQUIRE(b[2][3] == b23);
891 1 : CATCH_REQUIRE(b[3][0] == b30);
892 1 : CATCH_REQUIRE(b[3][1] == b31);
893 1 : CATCH_REQUIRE(b[3][2] == b32);
894 1 : CATCH_REQUIRE(b[3][3] == b33);
895 :
896 : // setup C
897 : //
898 2 : snapdev::matrix<double> c(4, 4);
899 :
900 1 : CATCH_REQUIRE(c.rows() == 4);
901 1 : CATCH_REQUIRE(c.columns() == 4);
902 :
903 1 : double const c00(frand());
904 1 : double const c01(frand());
905 1 : double const c02(frand());
906 1 : double const c03(frand());
907 1 : double const c10(frand());
908 1 : double const c11(frand());
909 1 : double const c12(frand());
910 1 : double const c13(frand());
911 1 : double const c20(frand());
912 1 : double const c21(frand());
913 1 : double const c22(frand());
914 1 : double const c23(frand());
915 1 : double const c30(frand());
916 1 : double const c31(frand());
917 1 : double const c32(frand());
918 1 : double const c33(frand());
919 :
920 1 : c[0][0] = c00;
921 1 : c[0][1] = c01;
922 1 : c[0][2] = c02;
923 1 : c[0][3] = c03;
924 1 : c[1][0] = c10;
925 1 : c[1][1] = c11;
926 1 : c[1][2] = c12;
927 1 : c[1][3] = c13;
928 1 : c[2][0] = c20;
929 1 : c[2][1] = c21;
930 1 : c[2][2] = c22;
931 1 : c[2][3] = c23;
932 1 : c[3][0] = c30;
933 1 : c[3][1] = c31;
934 1 : c[3][2] = c32;
935 1 : c[3][3] = c33;
936 :
937 1 : CATCH_REQUIRE(c[0][0] == c00);
938 1 : CATCH_REQUIRE(c[0][1] == c01);
939 1 : CATCH_REQUIRE(c[0][2] == c02);
940 1 : CATCH_REQUIRE(c[0][3] == c03);
941 1 : CATCH_REQUIRE(c[1][0] == c10);
942 1 : CATCH_REQUIRE(c[1][1] == c11);
943 1 : CATCH_REQUIRE(c[1][2] == c12);
944 1 : CATCH_REQUIRE(c[1][3] == c13);
945 1 : CATCH_REQUIRE(c[2][0] == c20);
946 1 : CATCH_REQUIRE(c[2][1] == c21);
947 1 : CATCH_REQUIRE(c[2][2] == c22);
948 1 : CATCH_REQUIRE(c[2][3] == c23);
949 1 : CATCH_REQUIRE(c[3][0] == c30);
950 1 : CATCH_REQUIRE(c[3][1] == c31);
951 1 : CATCH_REQUIRE(c[3][2] == c32);
952 1 : CATCH_REQUIRE(c[3][3] == c33);
953 :
954 : // run operation C = A + B
955 : //
956 1 : c = a + b;
957 :
958 1 : CATCH_REQUIRE(a[0][0] == a00);
959 1 : CATCH_REQUIRE(a[0][1] == a01);
960 1 : CATCH_REQUIRE(a[0][2] == a02);
961 1 : CATCH_REQUIRE(a[0][3] == a03);
962 1 : CATCH_REQUIRE(a[1][0] == a10);
963 1 : CATCH_REQUIRE(a[1][1] == a11);
964 1 : CATCH_REQUIRE(a[1][2] == a12);
965 1 : CATCH_REQUIRE(a[1][3] == a13);
966 1 : CATCH_REQUIRE(a[2][0] == a20);
967 1 : CATCH_REQUIRE(a[2][1] == a21);
968 1 : CATCH_REQUIRE(a[2][2] == a22);
969 1 : CATCH_REQUIRE(a[2][3] == a23);
970 1 : CATCH_REQUIRE(a[3][0] == a30);
971 1 : CATCH_REQUIRE(a[3][1] == a31);
972 1 : CATCH_REQUIRE(a[3][2] == a32);
973 1 : CATCH_REQUIRE(a[3][3] == a33);
974 :
975 1 : CATCH_REQUIRE(b[0][0] == b00);
976 1 : CATCH_REQUIRE(b[0][1] == b01);
977 1 : CATCH_REQUIRE(b[0][2] == b02);
978 1 : CATCH_REQUIRE(b[0][3] == b03);
979 1 : CATCH_REQUIRE(b[1][0] == b10);
980 1 : CATCH_REQUIRE(b[1][1] == b11);
981 1 : CATCH_REQUIRE(b[1][2] == b12);
982 1 : CATCH_REQUIRE(b[1][3] == b13);
983 1 : CATCH_REQUIRE(b[2][0] == b20);
984 1 : CATCH_REQUIRE(b[2][1] == b21);
985 1 : CATCH_REQUIRE(b[2][2] == b22);
986 1 : CATCH_REQUIRE(b[2][3] == b23);
987 1 : CATCH_REQUIRE(b[3][0] == b30);
988 1 : CATCH_REQUIRE(b[3][1] == b31);
989 1 : CATCH_REQUIRE(b[3][2] == b32);
990 1 : CATCH_REQUIRE(b[3][3] == b33);
991 :
992 : // this could fail
993 : //CATCH_REQUIRE_FALSE(c[0][0] == c00);
994 : //CATCH_REQUIRE_FALSE(c[0][1] == c01);
995 : //CATCH_REQUIRE_FALSE(c[0][2] == c02);
996 : //CATCH_REQUIRE_FALSE(c[0][3] == c03);
997 : //CATCH_REQUIRE_FALSE(c[1][0] == c10);
998 : //CATCH_REQUIRE_FALSE(c[1][1] == c11);
999 : //CATCH_REQUIRE_FALSE(c[1][2] == c12);
1000 : //CATCH_REQUIRE_FALSE(c[1][3] == c13);
1001 : //CATCH_REQUIRE_FALSE(c[2][0] == c20);
1002 : //CATCH_REQUIRE_FALSE(c[2][1] == c21);
1003 : //CATCH_REQUIRE_FALSE(c[2][2] == c22);
1004 : //CATCH_REQUIRE_FALSE(c[2][3] == c23);
1005 : //CATCH_REQUIRE_FALSE(c[3][0] == c30);
1006 : //CATCH_REQUIRE_FALSE(c[3][1] == c31);
1007 : //CATCH_REQUIRE_FALSE(c[3][2] == c32);
1008 : //CATCH_REQUIRE_FALSE(c[3][3] == c33);
1009 :
1010 1 : CATCH_REQUIRE(c[0][0] == a00 + b00);
1011 1 : CATCH_REQUIRE(c[0][1] == a01 + b01);
1012 1 : CATCH_REQUIRE(c[0][2] == a02 + b02);
1013 1 : CATCH_REQUIRE(c[0][3] == a03 + b03);
1014 1 : CATCH_REQUIRE(c[1][0] == a10 + b10);
1015 1 : CATCH_REQUIRE(c[1][1] == a11 + b11);
1016 1 : CATCH_REQUIRE(c[1][2] == a12 + b12);
1017 1 : CATCH_REQUIRE(c[1][3] == a13 + b13);
1018 1 : CATCH_REQUIRE(c[2][0] == a20 + b20);
1019 1 : CATCH_REQUIRE(c[2][1] == a21 + b21);
1020 1 : CATCH_REQUIRE(c[2][2] == a22 + b22);
1021 1 : CATCH_REQUIRE(c[2][3] == a23 + b23);
1022 1 : CATCH_REQUIRE(c[3][0] == a30 + b30);
1023 1 : CATCH_REQUIRE(c[3][1] == a31 + b31);
1024 1 : CATCH_REQUIRE(c[3][2] == a32 + b32);
1025 1 : CATCH_REQUIRE(c[3][3] == a33 + b33);
1026 : }
1027 : CATCH_END_SECTION()
1028 :
1029 8 : CATCH_START_SECTION("matrix: a+=b")
1030 : {
1031 : // setup A
1032 : //
1033 2 : snapdev::matrix<double> a(4, 4);
1034 :
1035 1 : CATCH_REQUIRE(a.rows() == 4);
1036 1 : CATCH_REQUIRE(a.columns() == 4);
1037 :
1038 1 : double const a00(frand());
1039 1 : double const a01(frand());
1040 1 : double const a02(frand());
1041 1 : double const a03(frand());
1042 1 : double const a10(frand());
1043 1 : double const a11(frand());
1044 1 : double const a12(frand());
1045 1 : double const a13(frand());
1046 1 : double const a20(frand());
1047 1 : double const a21(frand());
1048 1 : double const a22(frand());
1049 1 : double const a23(frand());
1050 1 : double const a30(frand());
1051 1 : double const a31(frand());
1052 1 : double const a32(frand());
1053 1 : double const a33(frand());
1054 :
1055 1 : a[0][0] = a00;
1056 1 : a[0][1] = a01;
1057 1 : a[0][2] = a02;
1058 1 : a[0][3] = a03;
1059 1 : a[1][0] = a10;
1060 1 : a[1][1] = a11;
1061 1 : a[1][2] = a12;
1062 1 : a[1][3] = a13;
1063 1 : a[2][0] = a20;
1064 1 : a[2][1] = a21;
1065 1 : a[2][2] = a22;
1066 1 : a[2][3] = a23;
1067 1 : a[3][0] = a30;
1068 1 : a[3][1] = a31;
1069 1 : a[3][2] = a32;
1070 1 : a[3][3] = a33;
1071 :
1072 1 : CATCH_REQUIRE(a[0][0] == a00);
1073 1 : CATCH_REQUIRE(a[0][1] == a01);
1074 1 : CATCH_REQUIRE(a[0][2] == a02);
1075 1 : CATCH_REQUIRE(a[0][3] == a03);
1076 1 : CATCH_REQUIRE(a[1][0] == a10);
1077 1 : CATCH_REQUIRE(a[1][1] == a11);
1078 1 : CATCH_REQUIRE(a[1][2] == a12);
1079 1 : CATCH_REQUIRE(a[1][3] == a13);
1080 1 : CATCH_REQUIRE(a[2][0] == a20);
1081 1 : CATCH_REQUIRE(a[2][1] == a21);
1082 1 : CATCH_REQUIRE(a[2][2] == a22);
1083 1 : CATCH_REQUIRE(a[2][3] == a23);
1084 1 : CATCH_REQUIRE(a[3][0] == a30);
1085 1 : CATCH_REQUIRE(a[3][1] == a31);
1086 1 : CATCH_REQUIRE(a[3][2] == a32);
1087 1 : CATCH_REQUIRE(a[3][3] == a33);
1088 :
1089 : // setup B
1090 : //
1091 2 : snapdev::matrix<double> b(4, 4);
1092 :
1093 1 : CATCH_REQUIRE(b.rows() == 4);
1094 1 : CATCH_REQUIRE(b.columns() == 4);
1095 :
1096 1 : double const b00(frand());
1097 1 : double const b01(frand());
1098 1 : double const b02(frand());
1099 1 : double const b03(frand());
1100 1 : double const b10(frand());
1101 1 : double const b11(frand());
1102 1 : double const b12(frand());
1103 1 : double const b13(frand());
1104 1 : double const b20(frand());
1105 1 : double const b21(frand());
1106 1 : double const b22(frand());
1107 1 : double const b23(frand());
1108 1 : double const b30(frand());
1109 1 : double const b31(frand());
1110 1 : double const b32(frand());
1111 1 : double const b33(frand());
1112 :
1113 1 : b[0][0] = b00;
1114 1 : b[0][1] = b01;
1115 1 : b[0][2] = b02;
1116 1 : b[0][3] = b03;
1117 1 : b[1][0] = b10;
1118 1 : b[1][1] = b11;
1119 1 : b[1][2] = b12;
1120 1 : b[1][3] = b13;
1121 1 : b[2][0] = b20;
1122 1 : b[2][1] = b21;
1123 1 : b[2][2] = b22;
1124 1 : b[2][3] = b23;
1125 1 : b[3][0] = b30;
1126 1 : b[3][1] = b31;
1127 1 : b[3][2] = b32;
1128 1 : b[3][3] = b33;
1129 :
1130 1 : CATCH_REQUIRE(b[0][0] == b00);
1131 1 : CATCH_REQUIRE(b[0][1] == b01);
1132 1 : CATCH_REQUIRE(b[0][2] == b02);
1133 1 : CATCH_REQUIRE(b[0][3] == b03);
1134 1 : CATCH_REQUIRE(b[1][0] == b10);
1135 1 : CATCH_REQUIRE(b[1][1] == b11);
1136 1 : CATCH_REQUIRE(b[1][2] == b12);
1137 1 : CATCH_REQUIRE(b[1][3] == b13);
1138 1 : CATCH_REQUIRE(b[2][0] == b20);
1139 1 : CATCH_REQUIRE(b[2][1] == b21);
1140 1 : CATCH_REQUIRE(b[2][2] == b22);
1141 1 : CATCH_REQUIRE(b[2][3] == b23);
1142 1 : CATCH_REQUIRE(b[3][0] == b30);
1143 1 : CATCH_REQUIRE(b[3][1] == b31);
1144 1 : CATCH_REQUIRE(b[3][2] == b32);
1145 1 : CATCH_REQUIRE(b[3][3] == b33);
1146 :
1147 : // run operation A += B
1148 : //
1149 1 : a += b;
1150 :
1151 : // this could fail if any bXX is 0.0
1152 : //
1153 : //CATCH_REQUIRE_FALSE(a[0][0] == a00);
1154 : //CATCH_REQUIRE_FALSE(a[0][1] == a01);
1155 : //CATCH_REQUIRE_FALSE(a[0][2] == a02);
1156 : //CATCH_REQUIRE_FALSE(a[0][3] == a03);
1157 : //CATCH_REQUIRE_FALSE(a[1][0] == a10);
1158 : //CATCH_REQUIRE_FALSE(a[1][1] == a11);
1159 : //CATCH_REQUIRE_FALSE(a[1][2] == a12);
1160 : //CATCH_REQUIRE_FALSE(a[1][3] == a13);
1161 : //CATCH_REQUIRE_FALSE(a[2][0] == a20);
1162 : //CATCH_REQUIRE_FALSE(a[2][1] == a21);
1163 : //CATCH_REQUIRE_FALSE(a[2][2] == a22);
1164 : //CATCH_REQUIRE_FALSE(a[2][3] == a23);
1165 : //CATCH_REQUIRE_FALSE(a[3][0] == a30);
1166 : //CATCH_REQUIRE_FALSE(a[3][1] == a31);
1167 : //CATCH_REQUIRE_FALSE(a[3][2] == a32);
1168 : //CATCH_REQUIRE_FALSE(a[3][3] == a33);
1169 :
1170 1 : CATCH_REQUIRE(b[0][0] == b00);
1171 1 : CATCH_REQUIRE(b[0][1] == b01);
1172 1 : CATCH_REQUIRE(b[0][2] == b02);
1173 1 : CATCH_REQUIRE(b[0][3] == b03);
1174 1 : CATCH_REQUIRE(b[1][0] == b10);
1175 1 : CATCH_REQUIRE(b[1][1] == b11);
1176 1 : CATCH_REQUIRE(b[1][2] == b12);
1177 1 : CATCH_REQUIRE(b[1][3] == b13);
1178 1 : CATCH_REQUIRE(b[2][0] == b20);
1179 1 : CATCH_REQUIRE(b[2][1] == b21);
1180 1 : CATCH_REQUIRE(b[2][2] == b22);
1181 1 : CATCH_REQUIRE(b[2][3] == b23);
1182 1 : CATCH_REQUIRE(b[3][0] == b30);
1183 1 : CATCH_REQUIRE(b[3][1] == b31);
1184 1 : CATCH_REQUIRE(b[3][2] == b32);
1185 1 : CATCH_REQUIRE(b[3][3] == b33);
1186 :
1187 1 : CATCH_REQUIRE(a[0][0] == a00 + b00);
1188 1 : CATCH_REQUIRE(a[0][1] == a01 + b01);
1189 1 : CATCH_REQUIRE(a[0][2] == a02 + b02);
1190 1 : CATCH_REQUIRE(a[0][3] == a03 + b03);
1191 1 : CATCH_REQUIRE(a[1][0] == a10 + b10);
1192 1 : CATCH_REQUIRE(a[1][1] == a11 + b11);
1193 1 : CATCH_REQUIRE(a[1][2] == a12 + b12);
1194 1 : CATCH_REQUIRE(a[1][3] == a13 + b13);
1195 1 : CATCH_REQUIRE(a[2][0] == a20 + b20);
1196 1 : CATCH_REQUIRE(a[2][1] == a21 + b21);
1197 1 : CATCH_REQUIRE(a[2][2] == a22 + b22);
1198 1 : CATCH_REQUIRE(a[2][3] == a23 + b23);
1199 1 : CATCH_REQUIRE(a[3][0] == a30 + b30);
1200 1 : CATCH_REQUIRE(a[3][1] == a31 + b31);
1201 1 : CATCH_REQUIRE(a[3][2] == a32 + b32);
1202 1 : CATCH_REQUIRE(a[3][3] == a33 + b33);
1203 : }
1204 : CATCH_END_SECTION()
1205 : }
1206 :
1207 : // create two random 4x4 matrices and make sure the add works
1208 : //
1209 16 : CATCH_GIVEN("sub")
1210 : {
1211 8 : CATCH_START_SECTION("matrix: b=a-<scalar>")
1212 : {
1213 : // setup A
1214 : //
1215 2 : snapdev::matrix<double> a(4, 4);
1216 :
1217 1 : CATCH_REQUIRE(a.rows() == 4);
1218 1 : CATCH_REQUIRE(a.columns() == 4);
1219 :
1220 1 : double const a00(frand());
1221 1 : double const a01(frand());
1222 1 : double const a02(frand());
1223 1 : double const a03(frand());
1224 1 : double const a10(frand());
1225 1 : double const a11(frand());
1226 1 : double const a12(frand());
1227 1 : double const a13(frand());
1228 1 : double const a20(frand());
1229 1 : double const a21(frand());
1230 1 : double const a22(frand());
1231 1 : double const a23(frand());
1232 1 : double const a30(frand());
1233 1 : double const a31(frand());
1234 1 : double const a32(frand());
1235 1 : double const a33(frand());
1236 :
1237 1 : a[0][0] = a00;
1238 1 : a[0][1] = a01;
1239 1 : a[0][2] = a02;
1240 1 : a[0][3] = a03;
1241 1 : a[1][0] = a10;
1242 1 : a[1][1] = a11;
1243 1 : a[1][2] = a12;
1244 1 : a[1][3] = a13;
1245 1 : a[2][0] = a20;
1246 1 : a[2][1] = a21;
1247 1 : a[2][2] = a22;
1248 1 : a[2][3] = a23;
1249 1 : a[3][0] = a30;
1250 1 : a[3][1] = a31;
1251 1 : a[3][2] = a32;
1252 1 : a[3][3] = a33;
1253 :
1254 1 : CATCH_REQUIRE(a[0][0] == a00);
1255 1 : CATCH_REQUIRE(a[0][1] == a01);
1256 1 : CATCH_REQUIRE(a[0][2] == a02);
1257 1 : CATCH_REQUIRE(a[0][3] == a03);
1258 1 : CATCH_REQUIRE(a[1][0] == a10);
1259 1 : CATCH_REQUIRE(a[1][1] == a11);
1260 1 : CATCH_REQUIRE(a[1][2] == a12);
1261 1 : CATCH_REQUIRE(a[1][3] == a13);
1262 1 : CATCH_REQUIRE(a[2][0] == a20);
1263 1 : CATCH_REQUIRE(a[2][1] == a21);
1264 1 : CATCH_REQUIRE(a[2][2] == a22);
1265 1 : CATCH_REQUIRE(a[2][3] == a23);
1266 1 : CATCH_REQUIRE(a[3][0] == a30);
1267 1 : CATCH_REQUIRE(a[3][1] == a31);
1268 1 : CATCH_REQUIRE(a[3][2] == a32);
1269 1 : CATCH_REQUIRE(a[3][3] == a33);
1270 :
1271 : // setup B
1272 : //
1273 2 : snapdev::matrix<double> b(4, 4);
1274 :
1275 1 : CATCH_REQUIRE(b.rows() == 4);
1276 1 : CATCH_REQUIRE(b.columns() == 4);
1277 :
1278 1 : double const b00(frand());
1279 1 : double const b01(frand());
1280 1 : double const b02(frand());
1281 1 : double const b03(frand());
1282 1 : double const b10(frand());
1283 1 : double const b11(frand());
1284 1 : double const b12(frand());
1285 1 : double const b13(frand());
1286 1 : double const b20(frand());
1287 1 : double const b21(frand());
1288 1 : double const b22(frand());
1289 1 : double const b23(frand());
1290 1 : double const b30(frand());
1291 1 : double const b31(frand());
1292 1 : double const b32(frand());
1293 1 : double const b33(frand());
1294 :
1295 1 : b[0][0] = b00;
1296 1 : b[0][1] = b01;
1297 1 : b[0][2] = b02;
1298 1 : b[0][3] = b03;
1299 1 : b[1][0] = b10;
1300 1 : b[1][1] = b11;
1301 1 : b[1][2] = b12;
1302 1 : b[1][3] = b13;
1303 1 : b[2][0] = b20;
1304 1 : b[2][1] = b21;
1305 1 : b[2][2] = b22;
1306 1 : b[2][3] = b23;
1307 1 : b[3][0] = b30;
1308 1 : b[3][1] = b31;
1309 1 : b[3][2] = b32;
1310 1 : b[3][3] = b33;
1311 :
1312 1 : CATCH_REQUIRE(b[0][0] == b00);
1313 1 : CATCH_REQUIRE(b[0][1] == b01);
1314 1 : CATCH_REQUIRE(b[0][2] == b02);
1315 1 : CATCH_REQUIRE(b[0][3] == b03);
1316 1 : CATCH_REQUIRE(b[1][0] == b10);
1317 1 : CATCH_REQUIRE(b[1][1] == b11);
1318 1 : CATCH_REQUIRE(b[1][2] == b12);
1319 1 : CATCH_REQUIRE(b[1][3] == b13);
1320 1 : CATCH_REQUIRE(b[2][0] == b20);
1321 1 : CATCH_REQUIRE(b[2][1] == b21);
1322 1 : CATCH_REQUIRE(b[2][2] == b22);
1323 1 : CATCH_REQUIRE(b[2][3] == b23);
1324 1 : CATCH_REQUIRE(b[3][0] == b30);
1325 1 : CATCH_REQUIRE(b[3][1] == b31);
1326 1 : CATCH_REQUIRE(b[3][2] == b32);
1327 1 : CATCH_REQUIRE(b[3][3] == b33);
1328 :
1329 : // create a scalar for our test
1330 : //
1331 1 : double scalar(frand());
1332 1 : while(scalar == 0.0)
1333 : {
1334 0 : scalar = frand();
1335 : }
1336 :
1337 : // run operation B = A - <scalar>
1338 : //
1339 1 : b = a - scalar;
1340 :
1341 1 : CATCH_REQUIRE(a[0][0] == a00);
1342 1 : CATCH_REQUIRE(a[0][1] == a01);
1343 1 : CATCH_REQUIRE(a[0][2] == a02);
1344 1 : CATCH_REQUIRE(a[0][3] == a03);
1345 1 : CATCH_REQUIRE(a[1][0] == a10);
1346 1 : CATCH_REQUIRE(a[1][1] == a11);
1347 1 : CATCH_REQUIRE(a[1][2] == a12);
1348 1 : CATCH_REQUIRE(a[1][3] == a13);
1349 1 : CATCH_REQUIRE(a[2][0] == a20);
1350 1 : CATCH_REQUIRE(a[2][1] == a21);
1351 1 : CATCH_REQUIRE(a[2][2] == a22);
1352 1 : CATCH_REQUIRE(a[2][3] == a23);
1353 1 : CATCH_REQUIRE(a[3][0] == a30);
1354 1 : CATCH_REQUIRE(a[3][1] == a31);
1355 1 : CATCH_REQUIRE(a[3][2] == a32);
1356 1 : CATCH_REQUIRE(a[3][3] == a33);
1357 :
1358 : // this can't fail because we ensure scalar != 0
1359 : //
1360 1 : CATCH_REQUIRE_FALSE(b[0][0] == b00);
1361 1 : CATCH_REQUIRE_FALSE(b[0][1] == b01);
1362 1 : CATCH_REQUIRE_FALSE(b[0][2] == b02);
1363 1 : CATCH_REQUIRE_FALSE(b[0][3] == b03);
1364 1 : CATCH_REQUIRE_FALSE(b[1][0] == b10);
1365 1 : CATCH_REQUIRE_FALSE(b[1][1] == b11);
1366 1 : CATCH_REQUIRE_FALSE(b[1][2] == b12);
1367 1 : CATCH_REQUIRE_FALSE(b[1][3] == b13);
1368 1 : CATCH_REQUIRE_FALSE(b[2][0] == b20);
1369 1 : CATCH_REQUIRE_FALSE(b[2][1] == b21);
1370 1 : CATCH_REQUIRE_FALSE(b[2][2] == b22);
1371 1 : CATCH_REQUIRE_FALSE(b[2][3] == b23);
1372 1 : CATCH_REQUIRE_FALSE(b[3][0] == b30);
1373 1 : CATCH_REQUIRE_FALSE(b[3][1] == b31);
1374 1 : CATCH_REQUIRE_FALSE(b[3][2] == b32);
1375 1 : CATCH_REQUIRE_FALSE(b[3][3] == b33);
1376 :
1377 1 : CATCH_REQUIRE(b[0][0] == a00 - scalar);
1378 1 : CATCH_REQUIRE(b[0][1] == a01 - scalar);
1379 1 : CATCH_REQUIRE(b[0][2] == a02 - scalar);
1380 1 : CATCH_REQUIRE(b[0][3] == a03 - scalar);
1381 1 : CATCH_REQUIRE(b[1][0] == a10 - scalar);
1382 1 : CATCH_REQUIRE(b[1][1] == a11 - scalar);
1383 1 : CATCH_REQUIRE(b[1][2] == a12 - scalar);
1384 1 : CATCH_REQUIRE(b[1][3] == a13 - scalar);
1385 1 : CATCH_REQUIRE(b[2][0] == a20 - scalar);
1386 1 : CATCH_REQUIRE(b[2][1] == a21 - scalar);
1387 1 : CATCH_REQUIRE(b[2][2] == a22 - scalar);
1388 1 : CATCH_REQUIRE(b[2][3] == a23 - scalar);
1389 1 : CATCH_REQUIRE(b[3][0] == a30 - scalar);
1390 1 : CATCH_REQUIRE(b[3][1] == a31 - scalar);
1391 1 : CATCH_REQUIRE(b[3][2] == a32 - scalar);
1392 1 : CATCH_REQUIRE(b[3][3] == a33 - scalar);
1393 : }
1394 : CATCH_END_SECTION()
1395 :
1396 8 : CATCH_START_SECTION("matrix: a-=<scalar>")
1397 : {
1398 : // setup A
1399 : //
1400 2 : snapdev::matrix<double> a(4, 4);
1401 :
1402 1 : CATCH_REQUIRE(a.rows() == 4);
1403 1 : CATCH_REQUIRE(a.columns() == 4);
1404 :
1405 1 : double const a00(frand());
1406 1 : double const a01(frand());
1407 1 : double const a02(frand());
1408 1 : double const a03(frand());
1409 1 : double const a10(frand());
1410 1 : double const a11(frand());
1411 1 : double const a12(frand());
1412 1 : double const a13(frand());
1413 1 : double const a20(frand());
1414 1 : double const a21(frand());
1415 1 : double const a22(frand());
1416 1 : double const a23(frand());
1417 1 : double const a30(frand());
1418 1 : double const a31(frand());
1419 1 : double const a32(frand());
1420 1 : double const a33(frand());
1421 :
1422 1 : a[0][0] = a00;
1423 1 : a[0][1] = a01;
1424 1 : a[0][2] = a02;
1425 1 : a[0][3] = a03;
1426 1 : a[1][0] = a10;
1427 1 : a[1][1] = a11;
1428 1 : a[1][2] = a12;
1429 1 : a[1][3] = a13;
1430 1 : a[2][0] = a20;
1431 1 : a[2][1] = a21;
1432 1 : a[2][2] = a22;
1433 1 : a[2][3] = a23;
1434 1 : a[3][0] = a30;
1435 1 : a[3][1] = a31;
1436 1 : a[3][2] = a32;
1437 1 : a[3][3] = a33;
1438 :
1439 1 : CATCH_REQUIRE(a[0][0] == a00);
1440 1 : CATCH_REQUIRE(a[0][1] == a01);
1441 1 : CATCH_REQUIRE(a[0][2] == a02);
1442 1 : CATCH_REQUIRE(a[0][3] == a03);
1443 1 : CATCH_REQUIRE(a[1][0] == a10);
1444 1 : CATCH_REQUIRE(a[1][1] == a11);
1445 1 : CATCH_REQUIRE(a[1][2] == a12);
1446 1 : CATCH_REQUIRE(a[1][3] == a13);
1447 1 : CATCH_REQUIRE(a[2][0] == a20);
1448 1 : CATCH_REQUIRE(a[2][1] == a21);
1449 1 : CATCH_REQUIRE(a[2][2] == a22);
1450 1 : CATCH_REQUIRE(a[2][3] == a23);
1451 1 : CATCH_REQUIRE(a[3][0] == a30);
1452 1 : CATCH_REQUIRE(a[3][1] == a31);
1453 1 : CATCH_REQUIRE(a[3][2] == a32);
1454 1 : CATCH_REQUIRE(a[3][3] == a33);
1455 :
1456 : // create a scalar for our test
1457 : //
1458 1 : double scalar(frand());
1459 1 : while(scalar == 0.0)
1460 : {
1461 0 : scalar = frand();
1462 : }
1463 :
1464 : // run operation A -= <scalar>
1465 : //
1466 1 : a -= scalar;
1467 :
1468 : // this can't fail because we ensure scalar != 0
1469 : //
1470 1 : CATCH_REQUIRE_FALSE(a[0][0] == a00);
1471 1 : CATCH_REQUIRE_FALSE(a[0][1] == a01);
1472 1 : CATCH_REQUIRE_FALSE(a[0][2] == a02);
1473 1 : CATCH_REQUIRE_FALSE(a[0][3] == a03);
1474 1 : CATCH_REQUIRE_FALSE(a[1][0] == a10);
1475 1 : CATCH_REQUIRE_FALSE(a[1][1] == a11);
1476 1 : CATCH_REQUIRE_FALSE(a[1][2] == a12);
1477 1 : CATCH_REQUIRE_FALSE(a[1][3] == a13);
1478 1 : CATCH_REQUIRE_FALSE(a[2][0] == a20);
1479 1 : CATCH_REQUIRE_FALSE(a[2][1] == a21);
1480 1 : CATCH_REQUIRE_FALSE(a[2][2] == a22);
1481 1 : CATCH_REQUIRE_FALSE(a[2][3] == a23);
1482 1 : CATCH_REQUIRE_FALSE(a[3][0] == a30);
1483 1 : CATCH_REQUIRE_FALSE(a[3][1] == a31);
1484 1 : CATCH_REQUIRE_FALSE(a[3][2] == a32);
1485 1 : CATCH_REQUIRE_FALSE(a[3][3] == a33);
1486 :
1487 1 : CATCH_REQUIRE(a[0][0] == a00 - scalar);
1488 1 : CATCH_REQUIRE(a[0][1] == a01 - scalar);
1489 1 : CATCH_REQUIRE(a[0][2] == a02 - scalar);
1490 1 : CATCH_REQUIRE(a[0][3] == a03 - scalar);
1491 1 : CATCH_REQUIRE(a[1][0] == a10 - scalar);
1492 1 : CATCH_REQUIRE(a[1][1] == a11 - scalar);
1493 1 : CATCH_REQUIRE(a[1][2] == a12 - scalar);
1494 1 : CATCH_REQUIRE(a[1][3] == a13 - scalar);
1495 1 : CATCH_REQUIRE(a[2][0] == a20 - scalar);
1496 1 : CATCH_REQUIRE(a[2][1] == a21 - scalar);
1497 1 : CATCH_REQUIRE(a[2][2] == a22 - scalar);
1498 1 : CATCH_REQUIRE(a[2][3] == a23 - scalar);
1499 1 : CATCH_REQUIRE(a[3][0] == a30 - scalar);
1500 1 : CATCH_REQUIRE(a[3][1] == a31 - scalar);
1501 1 : CATCH_REQUIRE(a[3][2] == a32 - scalar);
1502 1 : CATCH_REQUIRE(a[3][3] == a33 - scalar);
1503 : }
1504 : CATCH_END_SECTION()
1505 :
1506 8 : CATCH_START_SECTION("matrix: c=a-b")
1507 : {
1508 : // setup A
1509 : //
1510 2 : snapdev::matrix<double> a(4, 4);
1511 :
1512 1 : CATCH_REQUIRE(a.rows() == 4);
1513 1 : CATCH_REQUIRE(a.columns() == 4);
1514 :
1515 1 : double const a00(frand());
1516 1 : double const a01(frand());
1517 1 : double const a02(frand());
1518 1 : double const a03(frand());
1519 1 : double const a10(frand());
1520 1 : double const a11(frand());
1521 1 : double const a12(frand());
1522 1 : double const a13(frand());
1523 1 : double const a20(frand());
1524 1 : double const a21(frand());
1525 1 : double const a22(frand());
1526 1 : double const a23(frand());
1527 1 : double const a30(frand());
1528 1 : double const a31(frand());
1529 1 : double const a32(frand());
1530 1 : double const a33(frand());
1531 :
1532 1 : a[0][0] = a00;
1533 1 : a[0][1] = a01;
1534 1 : a[0][2] = a02;
1535 1 : a[0][3] = a03;
1536 1 : a[1][0] = a10;
1537 1 : a[1][1] = a11;
1538 1 : a[1][2] = a12;
1539 1 : a[1][3] = a13;
1540 1 : a[2][0] = a20;
1541 1 : a[2][1] = a21;
1542 1 : a[2][2] = a22;
1543 1 : a[2][3] = a23;
1544 1 : a[3][0] = a30;
1545 1 : a[3][1] = a31;
1546 1 : a[3][2] = a32;
1547 1 : a[3][3] = a33;
1548 :
1549 1 : CATCH_REQUIRE(a[0][0] == a00);
1550 1 : CATCH_REQUIRE(a[0][1] == a01);
1551 1 : CATCH_REQUIRE(a[0][2] == a02);
1552 1 : CATCH_REQUIRE(a[0][3] == a03);
1553 1 : CATCH_REQUIRE(a[1][0] == a10);
1554 1 : CATCH_REQUIRE(a[1][1] == a11);
1555 1 : CATCH_REQUIRE(a[1][2] == a12);
1556 1 : CATCH_REQUIRE(a[1][3] == a13);
1557 1 : CATCH_REQUIRE(a[2][0] == a20);
1558 1 : CATCH_REQUIRE(a[2][1] == a21);
1559 1 : CATCH_REQUIRE(a[2][2] == a22);
1560 1 : CATCH_REQUIRE(a[2][3] == a23);
1561 1 : CATCH_REQUIRE(a[3][0] == a30);
1562 1 : CATCH_REQUIRE(a[3][1] == a31);
1563 1 : CATCH_REQUIRE(a[3][2] == a32);
1564 1 : CATCH_REQUIRE(a[3][3] == a33);
1565 :
1566 : // setup B
1567 : //
1568 2 : snapdev::matrix<double> b(4, 4);
1569 :
1570 1 : CATCH_REQUIRE(b.rows() == 4);
1571 1 : CATCH_REQUIRE(b.columns() == 4);
1572 :
1573 1 : double const b00(frand());
1574 1 : double const b01(frand());
1575 1 : double const b02(frand());
1576 1 : double const b03(frand());
1577 1 : double const b10(frand());
1578 1 : double const b11(frand());
1579 1 : double const b12(frand());
1580 1 : double const b13(frand());
1581 1 : double const b20(frand());
1582 1 : double const b21(frand());
1583 1 : double const b22(frand());
1584 1 : double const b23(frand());
1585 1 : double const b30(frand());
1586 1 : double const b31(frand());
1587 1 : double const b32(frand());
1588 1 : double const b33(frand());
1589 :
1590 1 : b[0][0] = b00;
1591 1 : b[0][1] = b01;
1592 1 : b[0][2] = b02;
1593 1 : b[0][3] = b03;
1594 1 : b[1][0] = b10;
1595 1 : b[1][1] = b11;
1596 1 : b[1][2] = b12;
1597 1 : b[1][3] = b13;
1598 1 : b[2][0] = b20;
1599 1 : b[2][1] = b21;
1600 1 : b[2][2] = b22;
1601 1 : b[2][3] = b23;
1602 1 : b[3][0] = b30;
1603 1 : b[3][1] = b31;
1604 1 : b[3][2] = b32;
1605 1 : b[3][3] = b33;
1606 :
1607 1 : CATCH_REQUIRE(b[0][0] == b00);
1608 1 : CATCH_REQUIRE(b[0][1] == b01);
1609 1 : CATCH_REQUIRE(b[0][2] == b02);
1610 1 : CATCH_REQUIRE(b[0][3] == b03);
1611 1 : CATCH_REQUIRE(b[1][0] == b10);
1612 1 : CATCH_REQUIRE(b[1][1] == b11);
1613 1 : CATCH_REQUIRE(b[1][2] == b12);
1614 1 : CATCH_REQUIRE(b[1][3] == b13);
1615 1 : CATCH_REQUIRE(b[2][0] == b20);
1616 1 : CATCH_REQUIRE(b[2][1] == b21);
1617 1 : CATCH_REQUIRE(b[2][2] == b22);
1618 1 : CATCH_REQUIRE(b[2][3] == b23);
1619 1 : CATCH_REQUIRE(b[3][0] == b30);
1620 1 : CATCH_REQUIRE(b[3][1] == b31);
1621 1 : CATCH_REQUIRE(b[3][2] == b32);
1622 1 : CATCH_REQUIRE(b[3][3] == b33);
1623 :
1624 : // setup C
1625 : //
1626 2 : snapdev::matrix<double> c(4, 4);
1627 :
1628 1 : CATCH_REQUIRE(c.rows() == 4);
1629 1 : CATCH_REQUIRE(c.columns() == 4);
1630 :
1631 1 : double const c00(frand());
1632 1 : double const c01(frand());
1633 1 : double const c02(frand());
1634 1 : double const c03(frand());
1635 1 : double const c10(frand());
1636 1 : double const c11(frand());
1637 1 : double const c12(frand());
1638 1 : double const c13(frand());
1639 1 : double const c20(frand());
1640 1 : double const c21(frand());
1641 1 : double const c22(frand());
1642 1 : double const c23(frand());
1643 1 : double const c30(frand());
1644 1 : double const c31(frand());
1645 1 : double const c32(frand());
1646 1 : double const c33(frand());
1647 :
1648 1 : c[0][0] = c00;
1649 1 : c[0][1] = c01;
1650 1 : c[0][2] = c02;
1651 1 : c[0][3] = c03;
1652 1 : c[1][0] = c10;
1653 1 : c[1][1] = c11;
1654 1 : c[1][2] = c12;
1655 1 : c[1][3] = c13;
1656 1 : c[2][0] = c20;
1657 1 : c[2][1] = c21;
1658 1 : c[2][2] = c22;
1659 1 : c[2][3] = c23;
1660 1 : c[3][0] = c30;
1661 1 : c[3][1] = c31;
1662 1 : c[3][2] = c32;
1663 1 : c[3][3] = c33;
1664 :
1665 1 : CATCH_REQUIRE(c[0][0] == c00);
1666 1 : CATCH_REQUIRE(c[0][1] == c01);
1667 1 : CATCH_REQUIRE(c[0][2] == c02);
1668 1 : CATCH_REQUIRE(c[0][3] == c03);
1669 1 : CATCH_REQUIRE(c[1][0] == c10);
1670 1 : CATCH_REQUIRE(c[1][1] == c11);
1671 1 : CATCH_REQUIRE(c[1][2] == c12);
1672 1 : CATCH_REQUIRE(c[1][3] == c13);
1673 1 : CATCH_REQUIRE(c[2][0] == c20);
1674 1 : CATCH_REQUIRE(c[2][1] == c21);
1675 1 : CATCH_REQUIRE(c[2][2] == c22);
1676 1 : CATCH_REQUIRE(c[2][3] == c23);
1677 1 : CATCH_REQUIRE(c[3][0] == c30);
1678 1 : CATCH_REQUIRE(c[3][1] == c31);
1679 1 : CATCH_REQUIRE(c[3][2] == c32);
1680 1 : CATCH_REQUIRE(c[3][3] == c33);
1681 :
1682 : // run operation C = A - B
1683 : //
1684 1 : c = a - b;
1685 :
1686 1 : CATCH_REQUIRE(a[0][0] == a00);
1687 1 : CATCH_REQUIRE(a[0][1] == a01);
1688 1 : CATCH_REQUIRE(a[0][2] == a02);
1689 1 : CATCH_REQUIRE(a[0][3] == a03);
1690 1 : CATCH_REQUIRE(a[1][0] == a10);
1691 1 : CATCH_REQUIRE(a[1][1] == a11);
1692 1 : CATCH_REQUIRE(a[1][2] == a12);
1693 1 : CATCH_REQUIRE(a[1][3] == a13);
1694 1 : CATCH_REQUIRE(a[2][0] == a20);
1695 1 : CATCH_REQUIRE(a[2][1] == a21);
1696 1 : CATCH_REQUIRE(a[2][2] == a22);
1697 1 : CATCH_REQUIRE(a[2][3] == a23);
1698 1 : CATCH_REQUIRE(a[3][0] == a30);
1699 1 : CATCH_REQUIRE(a[3][1] == a31);
1700 1 : CATCH_REQUIRE(a[3][2] == a32);
1701 1 : CATCH_REQUIRE(a[3][3] == a33);
1702 :
1703 1 : CATCH_REQUIRE(b[0][0] == b00);
1704 1 : CATCH_REQUIRE(b[0][1] == b01);
1705 1 : CATCH_REQUIRE(b[0][2] == b02);
1706 1 : CATCH_REQUIRE(b[0][3] == b03);
1707 1 : CATCH_REQUIRE(b[1][0] == b10);
1708 1 : CATCH_REQUIRE(b[1][1] == b11);
1709 1 : CATCH_REQUIRE(b[1][2] == b12);
1710 1 : CATCH_REQUIRE(b[1][3] == b13);
1711 1 : CATCH_REQUIRE(b[2][0] == b20);
1712 1 : CATCH_REQUIRE(b[2][1] == b21);
1713 1 : CATCH_REQUIRE(b[2][2] == b22);
1714 1 : CATCH_REQUIRE(b[2][3] == b23);
1715 1 : CATCH_REQUIRE(b[3][0] == b30);
1716 1 : CATCH_REQUIRE(b[3][1] == b31);
1717 1 : CATCH_REQUIRE(b[3][2] == b32);
1718 1 : CATCH_REQUIRE(b[3][3] == b33);
1719 :
1720 : // this could fail
1721 : //CATCH_REQUIRE_FALSE(c[0][0] == c00);
1722 : //CATCH_REQUIRE_FALSE(c[0][1] == c01);
1723 : //CATCH_REQUIRE_FALSE(c[0][2] == c02);
1724 : //CATCH_REQUIRE_FALSE(c[0][3] == c03);
1725 : //CATCH_REQUIRE_FALSE(c[1][0] == c10);
1726 : //CATCH_REQUIRE_FALSE(c[1][1] == c11);
1727 : //CATCH_REQUIRE_FALSE(c[1][2] == c12);
1728 : //CATCH_REQUIRE_FALSE(c[1][3] == c13);
1729 : //CATCH_REQUIRE_FALSE(c[2][0] == c20);
1730 : //CATCH_REQUIRE_FALSE(c[2][1] == c21);
1731 : //CATCH_REQUIRE_FALSE(c[2][2] == c22);
1732 : //CATCH_REQUIRE_FALSE(c[2][3] == c23);
1733 : //CATCH_REQUIRE_FALSE(c[3][0] == c30);
1734 : //CATCH_REQUIRE_FALSE(c[3][1] == c31);
1735 : //CATCH_REQUIRE_FALSE(c[3][2] == c32);
1736 : //CATCH_REQUIRE_FALSE(c[3][3] == c33);
1737 :
1738 1 : CATCH_REQUIRE(c[0][0] == a00 - b00);
1739 1 : CATCH_REQUIRE(c[0][1] == a01 - b01);
1740 1 : CATCH_REQUIRE(c[0][2] == a02 - b02);
1741 1 : CATCH_REQUIRE(c[0][3] == a03 - b03);
1742 1 : CATCH_REQUIRE(c[1][0] == a10 - b10);
1743 1 : CATCH_REQUIRE(c[1][1] == a11 - b11);
1744 1 : CATCH_REQUIRE(c[1][2] == a12 - b12);
1745 1 : CATCH_REQUIRE(c[1][3] == a13 - b13);
1746 1 : CATCH_REQUIRE(c[2][0] == a20 - b20);
1747 1 : CATCH_REQUIRE(c[2][1] == a21 - b21);
1748 1 : CATCH_REQUIRE(c[2][2] == a22 - b22);
1749 1 : CATCH_REQUIRE(c[2][3] == a23 - b23);
1750 1 : CATCH_REQUIRE(c[3][0] == a30 - b30);
1751 1 : CATCH_REQUIRE(c[3][1] == a31 - b31);
1752 1 : CATCH_REQUIRE(c[3][2] == a32 - b32);
1753 1 : CATCH_REQUIRE(c[3][3] == a33 - b33);
1754 : }
1755 : CATCH_END_SECTION()
1756 :
1757 8 : CATCH_START_SECTION("matrix: a-=b")
1758 : {
1759 : // setup A
1760 : //
1761 2 : snapdev::matrix<double> a(4, 4);
1762 :
1763 1 : CATCH_REQUIRE(a.rows() == 4);
1764 1 : CATCH_REQUIRE(a.columns() == 4);
1765 :
1766 1 : double const a00(frand());
1767 1 : double const a01(frand());
1768 1 : double const a02(frand());
1769 1 : double const a03(frand());
1770 1 : double const a10(frand());
1771 1 : double const a11(frand());
1772 1 : double const a12(frand());
1773 1 : double const a13(frand());
1774 1 : double const a20(frand());
1775 1 : double const a21(frand());
1776 1 : double const a22(frand());
1777 1 : double const a23(frand());
1778 1 : double const a30(frand());
1779 1 : double const a31(frand());
1780 1 : double const a32(frand());
1781 1 : double const a33(frand());
1782 :
1783 1 : a[0][0] = a00;
1784 1 : a[0][1] = a01;
1785 1 : a[0][2] = a02;
1786 1 : a[0][3] = a03;
1787 1 : a[1][0] = a10;
1788 1 : a[1][1] = a11;
1789 1 : a[1][2] = a12;
1790 1 : a[1][3] = a13;
1791 1 : a[2][0] = a20;
1792 1 : a[2][1] = a21;
1793 1 : a[2][2] = a22;
1794 1 : a[2][3] = a23;
1795 1 : a[3][0] = a30;
1796 1 : a[3][1] = a31;
1797 1 : a[3][2] = a32;
1798 1 : a[3][3] = a33;
1799 :
1800 1 : CATCH_REQUIRE(a[0][0] == a00);
1801 1 : CATCH_REQUIRE(a[0][1] == a01);
1802 1 : CATCH_REQUIRE(a[0][2] == a02);
1803 1 : CATCH_REQUIRE(a[0][3] == a03);
1804 1 : CATCH_REQUIRE(a[1][0] == a10);
1805 1 : CATCH_REQUIRE(a[1][1] == a11);
1806 1 : CATCH_REQUIRE(a[1][2] == a12);
1807 1 : CATCH_REQUIRE(a[1][3] == a13);
1808 1 : CATCH_REQUIRE(a[2][0] == a20);
1809 1 : CATCH_REQUIRE(a[2][1] == a21);
1810 1 : CATCH_REQUIRE(a[2][2] == a22);
1811 1 : CATCH_REQUIRE(a[2][3] == a23);
1812 1 : CATCH_REQUIRE(a[3][0] == a30);
1813 1 : CATCH_REQUIRE(a[3][1] == a31);
1814 1 : CATCH_REQUIRE(a[3][2] == a32);
1815 1 : CATCH_REQUIRE(a[3][3] == a33);
1816 :
1817 : // setup B
1818 : //
1819 2 : snapdev::matrix<double> b(4, 4);
1820 :
1821 1 : CATCH_REQUIRE(b.rows() == 4);
1822 1 : CATCH_REQUIRE(b.columns() == 4);
1823 :
1824 1 : double const b00(frand());
1825 1 : double const b01(frand());
1826 1 : double const b02(frand());
1827 1 : double const b03(frand());
1828 1 : double const b10(frand());
1829 1 : double const b11(frand());
1830 1 : double const b12(frand());
1831 1 : double const b13(frand());
1832 1 : double const b20(frand());
1833 1 : double const b21(frand());
1834 1 : double const b22(frand());
1835 1 : double const b23(frand());
1836 1 : double const b30(frand());
1837 1 : double const b31(frand());
1838 1 : double const b32(frand());
1839 1 : double const b33(frand());
1840 :
1841 1 : b[0][0] = b00;
1842 1 : b[0][1] = b01;
1843 1 : b[0][2] = b02;
1844 1 : b[0][3] = b03;
1845 1 : b[1][0] = b10;
1846 1 : b[1][1] = b11;
1847 1 : b[1][2] = b12;
1848 1 : b[1][3] = b13;
1849 1 : b[2][0] = b20;
1850 1 : b[2][1] = b21;
1851 1 : b[2][2] = b22;
1852 1 : b[2][3] = b23;
1853 1 : b[3][0] = b30;
1854 1 : b[3][1] = b31;
1855 1 : b[3][2] = b32;
1856 1 : b[3][3] = b33;
1857 :
1858 1 : CATCH_REQUIRE(b[0][0] == b00);
1859 1 : CATCH_REQUIRE(b[0][1] == b01);
1860 1 : CATCH_REQUIRE(b[0][2] == b02);
1861 1 : CATCH_REQUIRE(b[0][3] == b03);
1862 1 : CATCH_REQUIRE(b[1][0] == b10);
1863 1 : CATCH_REQUIRE(b[1][1] == b11);
1864 1 : CATCH_REQUIRE(b[1][2] == b12);
1865 1 : CATCH_REQUIRE(b[1][3] == b13);
1866 1 : CATCH_REQUIRE(b[2][0] == b20);
1867 1 : CATCH_REQUIRE(b[2][1] == b21);
1868 1 : CATCH_REQUIRE(b[2][2] == b22);
1869 1 : CATCH_REQUIRE(b[2][3] == b23);
1870 1 : CATCH_REQUIRE(b[3][0] == b30);
1871 1 : CATCH_REQUIRE(b[3][1] == b31);
1872 1 : CATCH_REQUIRE(b[3][2] == b32);
1873 1 : CATCH_REQUIRE(b[3][3] == b33);
1874 :
1875 : // run operation A -= B
1876 : //
1877 1 : a -= b;
1878 :
1879 : // this could fail if one of bXX is 0.0
1880 : //
1881 : //CATCH_REQUIRE(a[0][0] == a00);
1882 : //CATCH_REQUIRE(a[0][1] == a01);
1883 : //CATCH_REQUIRE(a[0][2] == a02);
1884 : //CATCH_REQUIRE(a[0][3] == a03);
1885 : //CATCH_REQUIRE(a[1][0] == a10);
1886 : //CATCH_REQUIRE(a[1][1] == a11);
1887 : //CATCH_REQUIRE(a[1][2] == a12);
1888 : //CATCH_REQUIRE(a[1][3] == a13);
1889 : //CATCH_REQUIRE(a[2][0] == a20);
1890 : //CATCH_REQUIRE(a[2][1] == a21);
1891 : //CATCH_REQUIRE(a[2][2] == a22);
1892 : //CATCH_REQUIRE(a[2][3] == a23);
1893 : //CATCH_REQUIRE(a[3][0] == a30);
1894 : //CATCH_REQUIRE(a[3][1] == a31);
1895 : //CATCH_REQUIRE(a[3][2] == a32);
1896 : //CATCH_REQUIRE(a[3][3] == a33);
1897 :
1898 1 : CATCH_REQUIRE(b[0][0] == b00);
1899 1 : CATCH_REQUIRE(b[0][1] == b01);
1900 1 : CATCH_REQUIRE(b[0][2] == b02);
1901 1 : CATCH_REQUIRE(b[0][3] == b03);
1902 1 : CATCH_REQUIRE(b[1][0] == b10);
1903 1 : CATCH_REQUIRE(b[1][1] == b11);
1904 1 : CATCH_REQUIRE(b[1][2] == b12);
1905 1 : CATCH_REQUIRE(b[1][3] == b13);
1906 1 : CATCH_REQUIRE(b[2][0] == b20);
1907 1 : CATCH_REQUIRE(b[2][1] == b21);
1908 1 : CATCH_REQUIRE(b[2][2] == b22);
1909 1 : CATCH_REQUIRE(b[2][3] == b23);
1910 1 : CATCH_REQUIRE(b[3][0] == b30);
1911 1 : CATCH_REQUIRE(b[3][1] == b31);
1912 1 : CATCH_REQUIRE(b[3][2] == b32);
1913 1 : CATCH_REQUIRE(b[3][3] == b33);
1914 :
1915 1 : CATCH_REQUIRE(a[0][0] == a00 - b00);
1916 1 : CATCH_REQUIRE(a[0][1] == a01 - b01);
1917 1 : CATCH_REQUIRE(a[0][2] == a02 - b02);
1918 1 : CATCH_REQUIRE(a[0][3] == a03 - b03);
1919 1 : CATCH_REQUIRE(a[1][0] == a10 - b10);
1920 1 : CATCH_REQUIRE(a[1][1] == a11 - b11);
1921 1 : CATCH_REQUIRE(a[1][2] == a12 - b12);
1922 1 : CATCH_REQUIRE(a[1][3] == a13 - b13);
1923 1 : CATCH_REQUIRE(a[2][0] == a20 - b20);
1924 1 : CATCH_REQUIRE(a[2][1] == a21 - b21);
1925 1 : CATCH_REQUIRE(a[2][2] == a22 - b22);
1926 1 : CATCH_REQUIRE(a[2][3] == a23 - b23);
1927 1 : CATCH_REQUIRE(a[3][0] == a30 - b30);
1928 1 : CATCH_REQUIRE(a[3][1] == a31 - b31);
1929 1 : CATCH_REQUIRE(a[3][2] == a32 - b32);
1930 1 : CATCH_REQUIRE(a[3][3] == a33 - b33);
1931 : }
1932 : CATCH_END_SECTION()
1933 : }
1934 8 : }
1935 :
1936 :
1937 10 : CATCH_TEST_CASE("matrix_util", "[matrix]")
1938 : {
1939 : // various ways to change the data order, minor, enlarge
1940 : //
1941 16 : CATCH_GIVEN("move data")
1942 : {
1943 10 : CATCH_START_SECTION("matrix: minor")
1944 : {
1945 : // setup A
1946 : //
1947 2 : snapdev::matrix<double> a(4, 4);
1948 :
1949 1 : CATCH_REQUIRE(a.rows() == 4);
1950 1 : CATCH_REQUIRE(a.columns() == 4);
1951 :
1952 1 : double const a00(frand());
1953 1 : double const a01(frand());
1954 1 : double const a02(frand());
1955 1 : double const a03(frand());
1956 1 : double const a10(frand());
1957 1 : double const a11(frand());
1958 1 : double const a12(frand());
1959 1 : double const a13(frand());
1960 1 : double const a20(frand());
1961 1 : double const a21(frand());
1962 1 : double const a22(frand());
1963 1 : double const a23(frand());
1964 1 : double const a30(frand());
1965 1 : double const a31(frand());
1966 1 : double const a32(frand());
1967 1 : double const a33(frand());
1968 :
1969 1 : a[0][0] = a00;
1970 1 : a[0][1] = a01;
1971 1 : a[0][2] = a02;
1972 1 : a[0][3] = a03;
1973 1 : a[1][0] = a10;
1974 1 : a[1][1] = a11;
1975 1 : a[1][2] = a12;
1976 1 : a[1][3] = a13;
1977 1 : a[2][0] = a20;
1978 1 : a[2][1] = a21;
1979 1 : a[2][2] = a22;
1980 1 : a[2][3] = a23;
1981 1 : a[3][0] = a30;
1982 1 : a[3][1] = a31;
1983 1 : a[3][2] = a32;
1984 1 : a[3][3] = a33;
1985 :
1986 1 : CATCH_REQUIRE(a[0][0] == a00);
1987 1 : CATCH_REQUIRE(a[0][1] == a01);
1988 1 : CATCH_REQUIRE(a[0][2] == a02);
1989 1 : CATCH_REQUIRE(a[0][3] == a03);
1990 1 : CATCH_REQUIRE(a[1][0] == a10);
1991 1 : CATCH_REQUIRE(a[1][1] == a11);
1992 1 : CATCH_REQUIRE(a[1][2] == a12);
1993 1 : CATCH_REQUIRE(a[1][3] == a13);
1994 1 : CATCH_REQUIRE(a[2][0] == a20);
1995 1 : CATCH_REQUIRE(a[2][1] == a21);
1996 1 : CATCH_REQUIRE(a[2][2] == a22);
1997 1 : CATCH_REQUIRE(a[2][3] == a23);
1998 1 : CATCH_REQUIRE(a[3][0] == a30);
1999 1 : CATCH_REQUIRE(a[3][1] == a31);
2000 1 : CATCH_REQUIRE(a[3][2] == a32);
2001 1 : CATCH_REQUIRE(a[3][3] == a33);
2002 :
2003 : {
2004 2 : snapdev::matrix<double> p = a.minor_matrix(0, 0);
2005 :
2006 1 : CATCH_REQUIRE(p.rows() == 3);
2007 1 : CATCH_REQUIRE(p.columns() == 3);
2008 :
2009 1 : CATCH_REQUIRE(p[0][0] == a11);
2010 1 : CATCH_REQUIRE(p[0][1] == a12);
2011 1 : CATCH_REQUIRE(p[0][2] == a13);
2012 1 : CATCH_REQUIRE(p[1][0] == a21);
2013 1 : CATCH_REQUIRE(p[1][1] == a22);
2014 1 : CATCH_REQUIRE(p[1][2] == a23);
2015 1 : CATCH_REQUIRE(p[2][0] == a31);
2016 1 : CATCH_REQUIRE(p[2][1] == a32);
2017 1 : CATCH_REQUIRE(p[2][2] == a33);
2018 : }
2019 :
2020 : {
2021 2 : snapdev::matrix<double> p = a.minor_matrix(0, 1);
2022 :
2023 1 : CATCH_REQUIRE(p.rows() == 3);
2024 1 : CATCH_REQUIRE(p.columns() == 3);
2025 :
2026 1 : CATCH_REQUIRE(p[0][0] == a10);
2027 1 : CATCH_REQUIRE(p[0][1] == a12);
2028 1 : CATCH_REQUIRE(p[0][2] == a13);
2029 1 : CATCH_REQUIRE(p[1][0] == a20);
2030 1 : CATCH_REQUIRE(p[1][1] == a22);
2031 1 : CATCH_REQUIRE(p[1][2] == a23);
2032 1 : CATCH_REQUIRE(p[2][0] == a30);
2033 1 : CATCH_REQUIRE(p[2][1] == a32);
2034 1 : CATCH_REQUIRE(p[2][2] == a33);
2035 : }
2036 :
2037 : {
2038 2 : snapdev::matrix<double> p = a.minor_matrix(0, 2);
2039 :
2040 1 : CATCH_REQUIRE(p.rows() == 3);
2041 1 : CATCH_REQUIRE(p.columns() == 3);
2042 :
2043 1 : CATCH_REQUIRE(p[0][0] == a10);
2044 1 : CATCH_REQUIRE(p[0][1] == a11);
2045 1 : CATCH_REQUIRE(p[0][2] == a13);
2046 1 : CATCH_REQUIRE(p[1][0] == a20);
2047 1 : CATCH_REQUIRE(p[1][1] == a21);
2048 1 : CATCH_REQUIRE(p[1][2] == a23);
2049 1 : CATCH_REQUIRE(p[2][0] == a30);
2050 1 : CATCH_REQUIRE(p[2][1] == a31);
2051 1 : CATCH_REQUIRE(p[2][2] == a33);
2052 : }
2053 :
2054 : {
2055 2 : snapdev::matrix<double> p = a.minor_matrix(0, 3);
2056 :
2057 1 : CATCH_REQUIRE(p.rows() == 3);
2058 1 : CATCH_REQUIRE(p.columns() == 3);
2059 :
2060 1 : CATCH_REQUIRE(p[0][0] == a10);
2061 1 : CATCH_REQUIRE(p[0][1] == a11);
2062 1 : CATCH_REQUIRE(p[0][2] == a12);
2063 1 : CATCH_REQUIRE(p[1][0] == a20);
2064 1 : CATCH_REQUIRE(p[1][1] == a21);
2065 1 : CATCH_REQUIRE(p[1][2] == a22);
2066 1 : CATCH_REQUIRE(p[2][0] == a30);
2067 1 : CATCH_REQUIRE(p[2][1] == a31);
2068 1 : CATCH_REQUIRE(p[2][2] == a32);
2069 : }
2070 :
2071 : {
2072 2 : snapdev::matrix<double> p = a.minor_matrix(2, 1);
2073 :
2074 1 : CATCH_REQUIRE(p.rows() == 3);
2075 1 : CATCH_REQUIRE(p.columns() == 3);
2076 :
2077 1 : CATCH_REQUIRE(p[0][0] == a00);
2078 1 : CATCH_REQUIRE(p[0][1] == a02);
2079 1 : CATCH_REQUIRE(p[0][2] == a03);
2080 1 : CATCH_REQUIRE(p[1][0] == a10);
2081 1 : CATCH_REQUIRE(p[1][1] == a12);
2082 1 : CATCH_REQUIRE(p[1][2] == a13);
2083 1 : CATCH_REQUIRE(p[2][0] == a30);
2084 1 : CATCH_REQUIRE(p[2][1] == a32);
2085 1 : CATCH_REQUIRE(p[2][2] == a33);
2086 : }
2087 : }
2088 : CATCH_END_SECTION()
2089 :
2090 10 : CATCH_START_SECTION("matrix: transpose 4x4")
2091 : {
2092 : // setup A
2093 : //
2094 2 : snapdev::matrix<double> a(4, 4);
2095 :
2096 1 : CATCH_REQUIRE(a.rows() == 4);
2097 1 : CATCH_REQUIRE(a.columns() == 4);
2098 :
2099 1 : double const a00(frand());
2100 1 : double const a01(frand());
2101 1 : double const a02(frand());
2102 1 : double const a03(frand());
2103 1 : double const a10(frand());
2104 1 : double const a11(frand());
2105 1 : double const a12(frand());
2106 1 : double const a13(frand());
2107 1 : double const a20(frand());
2108 1 : double const a21(frand());
2109 1 : double const a22(frand());
2110 1 : double const a23(frand());
2111 1 : double const a30(frand());
2112 1 : double const a31(frand());
2113 1 : double const a32(frand());
2114 1 : double const a33(frand());
2115 :
2116 1 : a[0][0] = a00;
2117 1 : a[0][1] = a01;
2118 1 : a[0][2] = a02;
2119 1 : a[0][3] = a03;
2120 1 : a[1][0] = a10;
2121 1 : a[1][1] = a11;
2122 1 : a[1][2] = a12;
2123 1 : a[1][3] = a13;
2124 1 : a[2][0] = a20;
2125 1 : a[2][1] = a21;
2126 1 : a[2][2] = a22;
2127 1 : a[2][3] = a23;
2128 1 : a[3][0] = a30;
2129 1 : a[3][1] = a31;
2130 1 : a[3][2] = a32;
2131 1 : a[3][3] = a33;
2132 :
2133 1 : CATCH_REQUIRE(a[0][0] == a00);
2134 1 : CATCH_REQUIRE(a[0][1] == a01);
2135 1 : CATCH_REQUIRE(a[0][2] == a02);
2136 1 : CATCH_REQUIRE(a[0][3] == a03);
2137 1 : CATCH_REQUIRE(a[1][0] == a10);
2138 1 : CATCH_REQUIRE(a[1][1] == a11);
2139 1 : CATCH_REQUIRE(a[1][2] == a12);
2140 1 : CATCH_REQUIRE(a[1][3] == a13);
2141 1 : CATCH_REQUIRE(a[2][0] == a20);
2142 1 : CATCH_REQUIRE(a[2][1] == a21);
2143 1 : CATCH_REQUIRE(a[2][2] == a22);
2144 1 : CATCH_REQUIRE(a[2][3] == a23);
2145 1 : CATCH_REQUIRE(a[3][0] == a30);
2146 1 : CATCH_REQUIRE(a[3][1] == a31);
2147 1 : CATCH_REQUIRE(a[3][2] == a32);
2148 1 : CATCH_REQUIRE(a[3][3] == a33);
2149 :
2150 2 : snapdev::matrix<double> t = a.transpose();
2151 :
2152 1 : CATCH_REQUIRE(t.rows() == 4);
2153 1 : CATCH_REQUIRE(t.columns() == 4);
2154 :
2155 1 : CATCH_REQUIRE(t[0][0] == a00);
2156 1 : CATCH_REQUIRE(t[0][1] == a10);
2157 1 : CATCH_REQUIRE(t[0][2] == a20);
2158 1 : CATCH_REQUIRE(t[0][3] == a30);
2159 1 : CATCH_REQUIRE(t[1][0] == a01);
2160 1 : CATCH_REQUIRE(t[1][1] == a11);
2161 1 : CATCH_REQUIRE(t[1][2] == a21);
2162 1 : CATCH_REQUIRE(t[1][3] == a31);
2163 1 : CATCH_REQUIRE(t[2][0] == a02);
2164 1 : CATCH_REQUIRE(t[2][1] == a12);
2165 1 : CATCH_REQUIRE(t[2][2] == a22);
2166 1 : CATCH_REQUIRE(t[2][3] == a32);
2167 1 : CATCH_REQUIRE(t[3][0] == a03);
2168 1 : CATCH_REQUIRE(t[3][1] == a13);
2169 1 : CATCH_REQUIRE(t[3][2] == a23);
2170 1 : CATCH_REQUIRE(t[3][3] == a33);
2171 : }
2172 : CATCH_END_SECTION()
2173 :
2174 10 : CATCH_START_SECTION("matrix: transpose 6x2")
2175 : {
2176 : // setup A
2177 : //
2178 2 : snapdev::matrix<double> a(6, 2);
2179 :
2180 1 : CATCH_REQUIRE(a.rows() == 6);
2181 1 : CATCH_REQUIRE(a.columns() == 2);
2182 :
2183 1 : double const a00(frand());
2184 1 : double const a01(frand());
2185 1 : double const a10(frand());
2186 1 : double const a11(frand());
2187 1 : double const a20(frand());
2188 1 : double const a21(frand());
2189 1 : double const a30(frand());
2190 1 : double const a31(frand());
2191 1 : double const a40(frand());
2192 1 : double const a41(frand());
2193 1 : double const a50(frand());
2194 1 : double const a51(frand());
2195 :
2196 1 : a[0][0] = a00;
2197 1 : a[0][1] = a01;
2198 1 : a[1][0] = a10;
2199 1 : a[1][1] = a11;
2200 1 : a[2][0] = a20;
2201 1 : a[2][1] = a21;
2202 1 : a[3][0] = a30;
2203 1 : a[3][1] = a31;
2204 1 : a[4][0] = a40;
2205 1 : a[4][1] = a41;
2206 1 : a[5][0] = a50;
2207 1 : a[5][1] = a51;
2208 :
2209 1 : CATCH_REQUIRE(a[0][0] == a00);
2210 1 : CATCH_REQUIRE(a[0][1] == a01);
2211 1 : CATCH_REQUIRE(a[1][0] == a10);
2212 1 : CATCH_REQUIRE(a[1][1] == a11);
2213 1 : CATCH_REQUIRE(a[2][0] == a20);
2214 1 : CATCH_REQUIRE(a[2][1] == a21);
2215 1 : CATCH_REQUIRE(a[3][0] == a30);
2216 1 : CATCH_REQUIRE(a[3][1] == a31);
2217 1 : CATCH_REQUIRE(a[4][0] == a40);
2218 1 : CATCH_REQUIRE(a[4][1] == a41);
2219 1 : CATCH_REQUIRE(a[5][0] == a50);
2220 1 : CATCH_REQUIRE(a[5][1] == a51);
2221 :
2222 2 : snapdev::matrix<double> t = a.transpose();
2223 :
2224 1 : CATCH_REQUIRE(t.rows() == 2);
2225 1 : CATCH_REQUIRE(t.columns() == 6);
2226 :
2227 : // original
2228 : // a00 a01
2229 : // a10 a11
2230 : // a20 a21
2231 : // a30 a31
2232 : // a40 a41
2233 : // a50 a51
2234 : //
2235 : // transposed
2236 : // a00 a10 a20 a30 a40 a50
2237 : // a01 a11 a21 a31 a41 a51
2238 : //
2239 :
2240 1 : CATCH_REQUIRE(t[0][0] == a00);
2241 1 : CATCH_REQUIRE(t[0][1] == a10);
2242 1 : CATCH_REQUIRE(t[0][2] == a20);
2243 1 : CATCH_REQUIRE(t[0][3] == a30);
2244 1 : CATCH_REQUIRE(t[0][4] == a40);
2245 1 : CATCH_REQUIRE(t[0][5] == a50);
2246 1 : CATCH_REQUIRE(t[1][0] == a01);
2247 1 : CATCH_REQUIRE(t[1][1] == a11);
2248 1 : CATCH_REQUIRE(t[1][2] == a21);
2249 1 : CATCH_REQUIRE(t[1][3] == a31);
2250 1 : CATCH_REQUIRE(t[1][4] == a41);
2251 1 : CATCH_REQUIRE(t[1][5] == a51);
2252 : }
2253 : CATCH_END_SECTION()
2254 :
2255 10 : CATCH_START_SECTION("matrix: adjugate 2x2")
2256 : {
2257 : // setup A
2258 : //
2259 2 : snapdev::matrix<double> a(2, 2);
2260 :
2261 1 : CATCH_REQUIRE(a.rows() == 2);
2262 1 : CATCH_REQUIRE(a.columns() == 2);
2263 :
2264 1 : double const a00(frand());
2265 1 : double const a01(frand());
2266 1 : double const a10(frand());
2267 1 : double const a11(frand());
2268 :
2269 1 : a[0][0] = a00;
2270 1 : a[0][1] = a01;
2271 1 : a[1][0] = a10;
2272 1 : a[1][1] = a11;
2273 :
2274 1 : CATCH_REQUIRE(a[0][0] == a00);
2275 1 : CATCH_REQUIRE(a[0][1] == a01);
2276 1 : CATCH_REQUIRE(a[1][0] == a10);
2277 1 : CATCH_REQUIRE(a[1][1] == a11);
2278 :
2279 2 : snapdev::matrix<double> m = a.adjugate();
2280 :
2281 1 : CATCH_REQUIRE(m.rows() == 2);
2282 1 : CATCH_REQUIRE(m.columns() == 2);
2283 :
2284 1 : CATCH_REQUIRE(m[0][0] == a11);
2285 1 : CATCH_REQUIRE(m[0][1] == -a01);
2286 1 : CATCH_REQUIRE(m[1][0] == -a10);
2287 1 : CATCH_REQUIRE(m[1][1] == a00);
2288 : }
2289 : CATCH_END_SECTION()
2290 :
2291 10 : CATCH_START_SECTION("matrix: adjugate 3x3")
2292 : {
2293 : // setup A
2294 : //
2295 2 : snapdev::matrix<double> a(3, 3);
2296 :
2297 1 : CATCH_REQUIRE(a.rows() == 3);
2298 1 : CATCH_REQUIRE(a.columns() == 3);
2299 :
2300 1 : double const a00(frand());
2301 1 : double const a01(frand());
2302 1 : double const a02(frand());
2303 1 : double const a10(frand());
2304 1 : double const a11(frand());
2305 1 : double const a12(frand());
2306 1 : double const a20(frand());
2307 1 : double const a21(frand());
2308 1 : double const a22(frand());
2309 :
2310 1 : a[0][0] = a00;
2311 1 : a[0][1] = a01;
2312 1 : a[0][2] = a02;
2313 1 : a[1][0] = a10;
2314 1 : a[1][1] = a11;
2315 1 : a[1][2] = a12;
2316 1 : a[2][0] = a20;
2317 1 : a[2][1] = a21;
2318 1 : a[2][2] = a22;
2319 :
2320 1 : CATCH_REQUIRE(a[0][0] == a00);
2321 1 : CATCH_REQUIRE(a[0][1] == a01);
2322 1 : CATCH_REQUIRE(a[0][2] == a02);
2323 1 : CATCH_REQUIRE(a[1][0] == a10);
2324 1 : CATCH_REQUIRE(a[1][1] == a11);
2325 1 : CATCH_REQUIRE(a[1][2] == a12);
2326 1 : CATCH_REQUIRE(a[2][0] == a20);
2327 1 : CATCH_REQUIRE(a[2][1] == a21);
2328 1 : CATCH_REQUIRE(a[2][2] == a22);
2329 :
2330 2 : snapdev::matrix<double> m = a.adjugate();
2331 :
2332 1 : CATCH_REQUIRE(m.rows() == 3);
2333 1 : CATCH_REQUIRE(m.columns() == 3);
2334 :
2335 : {
2336 2 : snapdev::matrix<double> p(a.minor_matrix(0, 0));
2337 1 : double const e(m[0][0] - p.determinant());
2338 1 : CATCH_REQUIRE(fabs(e) < 0.0001);
2339 : }
2340 :
2341 : {
2342 2 : snapdev::matrix<double> p(a.minor_matrix(0, 1));
2343 1 : double const e(m[1][0] + p.determinant());
2344 1 : CATCH_REQUIRE(fabs(e) < 0.0001);
2345 : }
2346 :
2347 : {
2348 2 : snapdev::matrix<double> p(a.minor_matrix(0, 2));
2349 1 : double const e(m[2][0] - p.determinant());
2350 1 : CATCH_REQUIRE(fabs(e) < 0.0001);
2351 : }
2352 :
2353 : {
2354 2 : snapdev::matrix<double> p(a.minor_matrix(1, 0));
2355 1 : double const e(m[0][1] + p.determinant());
2356 1 : CATCH_REQUIRE(fabs(e) < 0.0001);
2357 : }
2358 :
2359 : {
2360 2 : snapdev::matrix<double> p(a.minor_matrix(1, 1));
2361 1 : double const e(m[1][1] - p.determinant());
2362 1 : CATCH_REQUIRE(fabs(e) < 0.0001);
2363 : }
2364 :
2365 : {
2366 2 : snapdev::matrix<double> p(a.minor_matrix(1, 2));
2367 1 : double const e(m[2][1] + p.determinant());
2368 1 : CATCH_REQUIRE(fabs(e) < 0.0001);
2369 : }
2370 :
2371 : {
2372 2 : snapdev::matrix<double> p(a.minor_matrix(2, 0));
2373 1 : double const e(m[0][2] - p.determinant());
2374 1 : CATCH_REQUIRE(fabs(e) < 0.0001);
2375 : }
2376 :
2377 : {
2378 2 : snapdev::matrix<double> p(a.minor_matrix(2, 1));
2379 1 : double const e(m[1][2] + p.determinant());
2380 1 : CATCH_REQUIRE(fabs(e) < 0.0001);
2381 : }
2382 :
2383 : {
2384 2 : snapdev::matrix<double> p(a.minor_matrix(2, 2));
2385 1 : double const e(m[2][2] - p.determinant());
2386 1 : CATCH_REQUIRE(fabs(e) < 0.0001);
2387 : }
2388 : }
2389 : CATCH_END_SECTION()
2390 : }
2391 :
2392 : // calculations such as the determinant
2393 : //
2394 16 : CATCH_GIVEN("calculations")
2395 : {
2396 6 : CATCH_START_SECTION("matrix: 2x2 determinant")
2397 : {
2398 : // setup A
2399 : //
2400 2 : snapdev::matrix<double> a(2, 2);
2401 :
2402 1 : CATCH_REQUIRE(a.rows() == 2);
2403 1 : CATCH_REQUIRE(a.columns() == 2);
2404 :
2405 1 : double const a00(frand());
2406 1 : double const a01(frand());
2407 1 : double const a10(frand());
2408 1 : double const a11(frand());
2409 :
2410 1 : a[0][0] = a00;
2411 1 : a[0][1] = a01;
2412 1 : a[1][0] = a10;
2413 1 : a[1][1] = a11;
2414 :
2415 1 : CATCH_REQUIRE(a[0][0] == a00);
2416 1 : CATCH_REQUIRE(a[0][1] == a01);
2417 1 : CATCH_REQUIRE(a[1][0] == a10);
2418 1 : CATCH_REQUIRE(a[1][1] == a11);
2419 :
2420 1 : double d(a.determinant());
2421 :
2422 : //return f_vector[0 + 0 * 2] * f_vector[1 + 1 * 2]
2423 : // - f_vector[1 + 0 * 2] * f_vector[0 + 1 * 2];
2424 1 : double e = a00 * a11 - a10 * a01;
2425 :
2426 1 : CATCH_REQUIRE(fabs(d - e) < 0.0001);
2427 : }
2428 : CATCH_END_SECTION()
2429 :
2430 6 : CATCH_START_SECTION("matrix: 3x3 determinant with specific data")
2431 : {
2432 : // setup A
2433 : //
2434 2 : snapdev::matrix<double> a(3, 3);
2435 :
2436 1 : CATCH_REQUIRE(a.rows() == 3);
2437 1 : CATCH_REQUIRE(a.columns() == 3);
2438 :
2439 1 : a[0][0] = 5.0;
2440 1 : a[0][1] = -2.0;
2441 1 : a[0][2] = 1.0;
2442 1 : a[1][0] = 0.0;
2443 1 : a[1][1] = 3.0;
2444 1 : a[1][2] = -1.0;
2445 1 : a[2][0] = 2.0;
2446 1 : a[2][1] = 0.0;
2447 1 : a[2][2] = 7.0;
2448 :
2449 1 : double d(a.determinant());
2450 :
2451 : // we know that the exact answer for this one is 103
2452 : //
2453 1 : CATCH_REQUIRE(fabs(d - 103.0) < 0.0001);
2454 : }
2455 : CATCH_END_SECTION()
2456 :
2457 6 : CATCH_START_SECTION("matrix: 3x3 determinant with random data")
2458 : {
2459 : // verify 10 times (with different values)
2460 : //
2461 11 : for(int repeat(0); repeat < 10; ++repeat)
2462 : {
2463 : // setup A
2464 : //
2465 20 : snapdev::matrix<double> a(3, 3);
2466 :
2467 10 : CATCH_REQUIRE(a.rows() == 3);
2468 10 : CATCH_REQUIRE(a.columns() == 3);
2469 :
2470 10 : double const a00(frand());
2471 10 : double const a01(frand());
2472 10 : double const a02(frand());
2473 10 : double const a10(frand());
2474 10 : double const a11(frand());
2475 10 : double const a12(frand());
2476 10 : double const a20(frand());
2477 10 : double const a21(frand());
2478 10 : double const a22(frand());
2479 :
2480 10 : a[0][0] = a00;
2481 10 : a[0][1] = a01;
2482 10 : a[0][2] = a02;
2483 10 : a[1][0] = a10;
2484 10 : a[1][1] = a11;
2485 10 : a[1][2] = a12;
2486 10 : a[2][0] = a20;
2487 10 : a[2][1] = a21;
2488 10 : a[2][2] = a22;
2489 :
2490 10 : CATCH_REQUIRE(a[0][0] == a00);
2491 10 : CATCH_REQUIRE(a[0][1] == a01);
2492 10 : CATCH_REQUIRE(a[0][2] == a02);
2493 10 : CATCH_REQUIRE(a[1][0] == a10);
2494 10 : CATCH_REQUIRE(a[1][1] == a11);
2495 10 : CATCH_REQUIRE(a[1][2] == a12);
2496 10 : CATCH_REQUIRE(a[2][0] == a20);
2497 10 : CATCH_REQUIRE(a[2][1] == a21);
2498 10 : CATCH_REQUIRE(a[2][2] == a22);
2499 :
2500 10 : double d(a.determinant());
2501 :
2502 : //return f_vector[0 + 0 * 2] * f_vector[1 + 1 * 2]
2503 : // - f_vector[1 + 0 * 2] * f_vector[0 + 1 * 2];
2504 20 : double e = a00 * a11 * a22 + a01 * a12 * a20 + a02 * a10 * a21
2505 10 : - a00 * a12 * a21 - a01 * a10 * a22 - a02 * a11 * a20;
2506 :
2507 10 : CATCH_REQUIRE(fabs(d - e) < 0.0001);
2508 : }
2509 : }
2510 : CATCH_END_SECTION()
2511 : }
2512 8 : }
2513 :
2514 :
2515 10 : CATCH_TEST_CASE("matrix_multiplicative", "[matrix]")
2516 : {
2517 : // create two random 4x4 matrices and make sure the add works
2518 : //
2519 16 : CATCH_GIVEN("mul")
2520 : {
2521 8 : CATCH_START_SECTION("matrix: b=a*<scalar>")
2522 : {
2523 : // setup A
2524 : //
2525 2 : snapdev::matrix<double> a(4, 4);
2526 :
2527 1 : CATCH_REQUIRE(a.rows() == 4);
2528 1 : CATCH_REQUIRE(a.columns() == 4);
2529 :
2530 1 : double const a00(frand());
2531 1 : double const a01(frand());
2532 1 : double const a02(frand());
2533 1 : double const a03(frand());
2534 1 : double const a10(frand());
2535 1 : double const a11(frand());
2536 1 : double const a12(frand());
2537 1 : double const a13(frand());
2538 1 : double const a20(frand());
2539 1 : double const a21(frand());
2540 1 : double const a22(frand());
2541 1 : double const a23(frand());
2542 1 : double const a30(frand());
2543 1 : double const a31(frand());
2544 1 : double const a32(frand());
2545 1 : double const a33(frand());
2546 :
2547 1 : a[0][0] = a00;
2548 1 : a[0][1] = a01;
2549 1 : a[0][2] = a02;
2550 1 : a[0][3] = a03;
2551 1 : a[1][0] = a10;
2552 1 : a[1][1] = a11;
2553 1 : a[1][2] = a12;
2554 1 : a[1][3] = a13;
2555 1 : a[2][0] = a20;
2556 1 : a[2][1] = a21;
2557 1 : a[2][2] = a22;
2558 1 : a[2][3] = a23;
2559 1 : a[3][0] = a30;
2560 1 : a[3][1] = a31;
2561 1 : a[3][2] = a32;
2562 1 : a[3][3] = a33;
2563 :
2564 1 : CATCH_REQUIRE(a[0][0] == a00);
2565 1 : CATCH_REQUIRE(a[0][1] == a01);
2566 1 : CATCH_REQUIRE(a[0][2] == a02);
2567 1 : CATCH_REQUIRE(a[0][3] == a03);
2568 1 : CATCH_REQUIRE(a[1][0] == a10);
2569 1 : CATCH_REQUIRE(a[1][1] == a11);
2570 1 : CATCH_REQUIRE(a[1][2] == a12);
2571 1 : CATCH_REQUIRE(a[1][3] == a13);
2572 1 : CATCH_REQUIRE(a[2][0] == a20);
2573 1 : CATCH_REQUIRE(a[2][1] == a21);
2574 1 : CATCH_REQUIRE(a[2][2] == a22);
2575 1 : CATCH_REQUIRE(a[2][3] == a23);
2576 1 : CATCH_REQUIRE(a[3][0] == a30);
2577 1 : CATCH_REQUIRE(a[3][1] == a31);
2578 1 : CATCH_REQUIRE(a[3][2] == a32);
2579 1 : CATCH_REQUIRE(a[3][3] == a33);
2580 :
2581 : // setup B
2582 : //
2583 2 : snapdev::matrix<double> b(4, 4);
2584 :
2585 1 : CATCH_REQUIRE(b.rows() == 4);
2586 1 : CATCH_REQUIRE(b.columns() == 4);
2587 :
2588 1 : double const b00(frand());
2589 1 : double const b01(frand());
2590 1 : double const b02(frand());
2591 1 : double const b03(frand());
2592 1 : double const b10(frand());
2593 1 : double const b11(frand());
2594 1 : double const b12(frand());
2595 1 : double const b13(frand());
2596 1 : double const b20(frand());
2597 1 : double const b21(frand());
2598 1 : double const b22(frand());
2599 1 : double const b23(frand());
2600 1 : double const b30(frand());
2601 1 : double const b31(frand());
2602 1 : double const b32(frand());
2603 1 : double const b33(frand());
2604 :
2605 1 : b[0][0] = b00;
2606 1 : b[0][1] = b01;
2607 1 : b[0][2] = b02;
2608 1 : b[0][3] = b03;
2609 1 : b[1][0] = b10;
2610 1 : b[1][1] = b11;
2611 1 : b[1][2] = b12;
2612 1 : b[1][3] = b13;
2613 1 : b[2][0] = b20;
2614 1 : b[2][1] = b21;
2615 1 : b[2][2] = b22;
2616 1 : b[2][3] = b23;
2617 1 : b[3][0] = b30;
2618 1 : b[3][1] = b31;
2619 1 : b[3][2] = b32;
2620 1 : b[3][3] = b33;
2621 :
2622 1 : CATCH_REQUIRE(b[0][0] == b00);
2623 1 : CATCH_REQUIRE(b[0][1] == b01);
2624 1 : CATCH_REQUIRE(b[0][2] == b02);
2625 1 : CATCH_REQUIRE(b[0][3] == b03);
2626 1 : CATCH_REQUIRE(b[1][0] == b10);
2627 1 : CATCH_REQUIRE(b[1][1] == b11);
2628 1 : CATCH_REQUIRE(b[1][2] == b12);
2629 1 : CATCH_REQUIRE(b[1][3] == b13);
2630 1 : CATCH_REQUIRE(b[2][0] == b20);
2631 1 : CATCH_REQUIRE(b[2][1] == b21);
2632 1 : CATCH_REQUIRE(b[2][2] == b22);
2633 1 : CATCH_REQUIRE(b[2][3] == b23);
2634 1 : CATCH_REQUIRE(b[3][0] == b30);
2635 1 : CATCH_REQUIRE(b[3][1] == b31);
2636 1 : CATCH_REQUIRE(b[3][2] == b32);
2637 1 : CATCH_REQUIRE(b[3][3] == b33);
2638 :
2639 : // create a scalar for our test
2640 : //
2641 1 : double scalar(frand());
2642 1 : while(scalar == 0.0 || scalar == 1.0)
2643 : {
2644 0 : scalar = frand();
2645 : }
2646 :
2647 : // run operation B = A + <scalar>
2648 : //
2649 1 : b = a * scalar;
2650 :
2651 1 : CATCH_REQUIRE(a[0][0] == a00);
2652 1 : CATCH_REQUIRE(a[0][1] == a01);
2653 1 : CATCH_REQUIRE(a[0][2] == a02);
2654 1 : CATCH_REQUIRE(a[0][3] == a03);
2655 1 : CATCH_REQUIRE(a[1][0] == a10);
2656 1 : CATCH_REQUIRE(a[1][1] == a11);
2657 1 : CATCH_REQUIRE(a[1][2] == a12);
2658 1 : CATCH_REQUIRE(a[1][3] == a13);
2659 1 : CATCH_REQUIRE(a[2][0] == a20);
2660 1 : CATCH_REQUIRE(a[2][1] == a21);
2661 1 : CATCH_REQUIRE(a[2][2] == a22);
2662 1 : CATCH_REQUIRE(a[2][3] == a23);
2663 1 : CATCH_REQUIRE(a[3][0] == a30);
2664 1 : CATCH_REQUIRE(a[3][1] == a31);
2665 1 : CATCH_REQUIRE(a[3][2] == a32);
2666 1 : CATCH_REQUIRE(a[3][3] == a33);
2667 :
2668 : // this can't fail because we ensure scalar != 0 or 1
2669 : //
2670 1 : CATCH_REQUIRE_FALSE(b[0][0] == b00);
2671 1 : CATCH_REQUIRE_FALSE(b[0][1] == b01);
2672 1 : CATCH_REQUIRE_FALSE(b[0][2] == b02);
2673 1 : CATCH_REQUIRE_FALSE(b[0][3] == b03);
2674 1 : CATCH_REQUIRE_FALSE(b[1][0] == b10);
2675 1 : CATCH_REQUIRE_FALSE(b[1][1] == b11);
2676 1 : CATCH_REQUIRE_FALSE(b[1][2] == b12);
2677 1 : CATCH_REQUIRE_FALSE(b[1][3] == b13);
2678 1 : CATCH_REQUIRE_FALSE(b[2][0] == b20);
2679 1 : CATCH_REQUIRE_FALSE(b[2][1] == b21);
2680 1 : CATCH_REQUIRE_FALSE(b[2][2] == b22);
2681 1 : CATCH_REQUIRE_FALSE(b[2][3] == b23);
2682 1 : CATCH_REQUIRE_FALSE(b[3][0] == b30);
2683 1 : CATCH_REQUIRE_FALSE(b[3][1] == b31);
2684 1 : CATCH_REQUIRE_FALSE(b[3][2] == b32);
2685 1 : CATCH_REQUIRE_FALSE(b[3][3] == b33);
2686 :
2687 1 : CATCH_REQUIRE(b[0][0] == a00 * scalar);
2688 1 : CATCH_REQUIRE(b[0][1] == a01 * scalar);
2689 1 : CATCH_REQUIRE(b[0][2] == a02 * scalar);
2690 1 : CATCH_REQUIRE(b[0][3] == a03 * scalar);
2691 1 : CATCH_REQUIRE(b[1][0] == a10 * scalar);
2692 1 : CATCH_REQUIRE(b[1][1] == a11 * scalar);
2693 1 : CATCH_REQUIRE(b[1][2] == a12 * scalar);
2694 1 : CATCH_REQUIRE(b[1][3] == a13 * scalar);
2695 1 : CATCH_REQUIRE(b[2][0] == a20 * scalar);
2696 1 : CATCH_REQUIRE(b[2][1] == a21 * scalar);
2697 1 : CATCH_REQUIRE(b[2][2] == a22 * scalar);
2698 1 : CATCH_REQUIRE(b[2][3] == a23 * scalar);
2699 1 : CATCH_REQUIRE(b[3][0] == a30 * scalar);
2700 1 : CATCH_REQUIRE(b[3][1] == a31 * scalar);
2701 1 : CATCH_REQUIRE(b[3][2] == a32 * scalar);
2702 1 : CATCH_REQUIRE(b[3][3] == a33 * scalar);
2703 : }
2704 : CATCH_END_SECTION()
2705 :
2706 8 : CATCH_START_SECTION("matrix: a*=<scalar>")
2707 : {
2708 : // setup A
2709 : //
2710 2 : snapdev::matrix<double> a(4, 4);
2711 :
2712 1 : CATCH_REQUIRE(a.rows() == 4);
2713 1 : CATCH_REQUIRE(a.columns() == 4);
2714 :
2715 1 : double const a00(frand());
2716 1 : double const a01(frand());
2717 1 : double const a02(frand());
2718 1 : double const a03(frand());
2719 1 : double const a10(frand());
2720 1 : double const a11(frand());
2721 1 : double const a12(frand());
2722 1 : double const a13(frand());
2723 1 : double const a20(frand());
2724 1 : double const a21(frand());
2725 1 : double const a22(frand());
2726 1 : double const a23(frand());
2727 1 : double const a30(frand());
2728 1 : double const a31(frand());
2729 1 : double const a32(frand());
2730 1 : double const a33(frand());
2731 :
2732 1 : a[0][0] = a00;
2733 1 : a[0][1] = a01;
2734 1 : a[0][2] = a02;
2735 1 : a[0][3] = a03;
2736 1 : a[1][0] = a10;
2737 1 : a[1][1] = a11;
2738 1 : a[1][2] = a12;
2739 1 : a[1][3] = a13;
2740 1 : a[2][0] = a20;
2741 1 : a[2][1] = a21;
2742 1 : a[2][2] = a22;
2743 1 : a[2][3] = a23;
2744 1 : a[3][0] = a30;
2745 1 : a[3][1] = a31;
2746 1 : a[3][2] = a32;
2747 1 : a[3][3] = a33;
2748 :
2749 1 : CATCH_REQUIRE(a[0][0] == a00);
2750 1 : CATCH_REQUIRE(a[0][1] == a01);
2751 1 : CATCH_REQUIRE(a[0][2] == a02);
2752 1 : CATCH_REQUIRE(a[0][3] == a03);
2753 1 : CATCH_REQUIRE(a[1][0] == a10);
2754 1 : CATCH_REQUIRE(a[1][1] == a11);
2755 1 : CATCH_REQUIRE(a[1][2] == a12);
2756 1 : CATCH_REQUIRE(a[1][3] == a13);
2757 1 : CATCH_REQUIRE(a[2][0] == a20);
2758 1 : CATCH_REQUIRE(a[2][1] == a21);
2759 1 : CATCH_REQUIRE(a[2][2] == a22);
2760 1 : CATCH_REQUIRE(a[2][3] == a23);
2761 1 : CATCH_REQUIRE(a[3][0] == a30);
2762 1 : CATCH_REQUIRE(a[3][1] == a31);
2763 1 : CATCH_REQUIRE(a[3][2] == a32);
2764 1 : CATCH_REQUIRE(a[3][3] == a33);
2765 :
2766 : // create a scalar for our test
2767 : //
2768 1 : double scalar(frand());
2769 1 : while(scalar == 0.0 || scalar == 1.0)
2770 : {
2771 0 : scalar = frand();
2772 : }
2773 :
2774 : // run operation A *= <scalar>
2775 : //
2776 1 : a *= scalar;
2777 :
2778 : // this can't fail because we ensure scalar != 0 or 1
2779 : //
2780 1 : CATCH_REQUIRE_FALSE(a[0][0] == a00);
2781 1 : CATCH_REQUIRE_FALSE(a[0][1] == a01);
2782 1 : CATCH_REQUIRE_FALSE(a[0][2] == a02);
2783 1 : CATCH_REQUIRE_FALSE(a[0][3] == a03);
2784 1 : CATCH_REQUIRE_FALSE(a[1][0] == a10);
2785 1 : CATCH_REQUIRE_FALSE(a[1][1] == a11);
2786 1 : CATCH_REQUIRE_FALSE(a[1][2] == a12);
2787 1 : CATCH_REQUIRE_FALSE(a[1][3] == a13);
2788 1 : CATCH_REQUIRE_FALSE(a[2][0] == a20);
2789 1 : CATCH_REQUIRE_FALSE(a[2][1] == a21);
2790 1 : CATCH_REQUIRE_FALSE(a[2][2] == a22);
2791 1 : CATCH_REQUIRE_FALSE(a[2][3] == a23);
2792 1 : CATCH_REQUIRE_FALSE(a[3][0] == a30);
2793 1 : CATCH_REQUIRE_FALSE(a[3][1] == a31);
2794 1 : CATCH_REQUIRE_FALSE(a[3][2] == a32);
2795 1 : CATCH_REQUIRE_FALSE(a[3][3] == a33);
2796 :
2797 1 : CATCH_REQUIRE(a[0][0] == a00 * scalar);
2798 1 : CATCH_REQUIRE(a[0][1] == a01 * scalar);
2799 1 : CATCH_REQUIRE(a[0][2] == a02 * scalar);
2800 1 : CATCH_REQUIRE(a[0][3] == a03 * scalar);
2801 1 : CATCH_REQUIRE(a[1][0] == a10 * scalar);
2802 1 : CATCH_REQUIRE(a[1][1] == a11 * scalar);
2803 1 : CATCH_REQUIRE(a[1][2] == a12 * scalar);
2804 1 : CATCH_REQUIRE(a[1][3] == a13 * scalar);
2805 1 : CATCH_REQUIRE(a[2][0] == a20 * scalar);
2806 1 : CATCH_REQUIRE(a[2][1] == a21 * scalar);
2807 1 : CATCH_REQUIRE(a[2][2] == a22 * scalar);
2808 1 : CATCH_REQUIRE(a[2][3] == a23 * scalar);
2809 1 : CATCH_REQUIRE(a[3][0] == a30 * scalar);
2810 1 : CATCH_REQUIRE(a[3][1] == a31 * scalar);
2811 1 : CATCH_REQUIRE(a[3][2] == a32 * scalar);
2812 1 : CATCH_REQUIRE(a[3][3] == a33 * scalar);
2813 : }
2814 : CATCH_END_SECTION()
2815 :
2816 8 : CATCH_START_SECTION("matrix: c=a*b")
2817 : {
2818 : // setup A
2819 : //
2820 2 : snapdev::matrix<double> a(4, 4);
2821 :
2822 1 : CATCH_REQUIRE(a.rows() == 4);
2823 1 : CATCH_REQUIRE(a.columns() == 4);
2824 :
2825 1 : double const a00(frand());
2826 1 : double const a01(frand());
2827 1 : double const a02(frand());
2828 1 : double const a03(frand());
2829 1 : double const a10(frand());
2830 1 : double const a11(frand());
2831 1 : double const a12(frand());
2832 1 : double const a13(frand());
2833 1 : double const a20(frand());
2834 1 : double const a21(frand());
2835 1 : double const a22(frand());
2836 1 : double const a23(frand());
2837 1 : double const a30(frand());
2838 1 : double const a31(frand());
2839 1 : double const a32(frand());
2840 1 : double const a33(frand());
2841 :
2842 1 : a[0][0] = a00;
2843 1 : a[0][1] = a01;
2844 1 : a[0][2] = a02;
2845 1 : a[0][3] = a03;
2846 1 : a[1][0] = a10;
2847 1 : a[1][1] = a11;
2848 1 : a[1][2] = a12;
2849 1 : a[1][3] = a13;
2850 1 : a[2][0] = a20;
2851 1 : a[2][1] = a21;
2852 1 : a[2][2] = a22;
2853 1 : a[2][3] = a23;
2854 1 : a[3][0] = a30;
2855 1 : a[3][1] = a31;
2856 1 : a[3][2] = a32;
2857 1 : a[3][3] = a33;
2858 :
2859 1 : CATCH_REQUIRE(a[0][0] == a00);
2860 1 : CATCH_REQUIRE(a[0][1] == a01);
2861 1 : CATCH_REQUIRE(a[0][2] == a02);
2862 1 : CATCH_REQUIRE(a[0][3] == a03);
2863 1 : CATCH_REQUIRE(a[1][0] == a10);
2864 1 : CATCH_REQUIRE(a[1][1] == a11);
2865 1 : CATCH_REQUIRE(a[1][2] == a12);
2866 1 : CATCH_REQUIRE(a[1][3] == a13);
2867 1 : CATCH_REQUIRE(a[2][0] == a20);
2868 1 : CATCH_REQUIRE(a[2][1] == a21);
2869 1 : CATCH_REQUIRE(a[2][2] == a22);
2870 1 : CATCH_REQUIRE(a[2][3] == a23);
2871 1 : CATCH_REQUIRE(a[3][0] == a30);
2872 1 : CATCH_REQUIRE(a[3][1] == a31);
2873 1 : CATCH_REQUIRE(a[3][2] == a32);
2874 1 : CATCH_REQUIRE(a[3][3] == a33);
2875 :
2876 : // setup B
2877 : //
2878 2 : snapdev::matrix<double> b(4, 4);
2879 :
2880 1 : CATCH_REQUIRE(b.rows() == 4);
2881 1 : CATCH_REQUIRE(b.columns() == 4);
2882 :
2883 1 : double const b00(frand());
2884 1 : double const b01(frand());
2885 1 : double const b02(frand());
2886 1 : double const b03(frand());
2887 1 : double const b10(frand());
2888 1 : double const b11(frand());
2889 1 : double const b12(frand());
2890 1 : double const b13(frand());
2891 1 : double const b20(frand());
2892 1 : double const b21(frand());
2893 1 : double const b22(frand());
2894 1 : double const b23(frand());
2895 1 : double const b30(frand());
2896 1 : double const b31(frand());
2897 1 : double const b32(frand());
2898 1 : double const b33(frand());
2899 :
2900 1 : b[0][0] = b00;
2901 1 : b[0][1] = b01;
2902 1 : b[0][2] = b02;
2903 1 : b[0][3] = b03;
2904 1 : b[1][0] = b10;
2905 1 : b[1][1] = b11;
2906 1 : b[1][2] = b12;
2907 1 : b[1][3] = b13;
2908 1 : b[2][0] = b20;
2909 1 : b[2][1] = b21;
2910 1 : b[2][2] = b22;
2911 1 : b[2][3] = b23;
2912 1 : b[3][0] = b30;
2913 1 : b[3][1] = b31;
2914 1 : b[3][2] = b32;
2915 1 : b[3][3] = b33;
2916 :
2917 1 : CATCH_REQUIRE(b[0][0] == b00);
2918 1 : CATCH_REQUIRE(b[0][1] == b01);
2919 1 : CATCH_REQUIRE(b[0][2] == b02);
2920 1 : CATCH_REQUIRE(b[0][3] == b03);
2921 1 : CATCH_REQUIRE(b[1][0] == b10);
2922 1 : CATCH_REQUIRE(b[1][1] == b11);
2923 1 : CATCH_REQUIRE(b[1][2] == b12);
2924 1 : CATCH_REQUIRE(b[1][3] == b13);
2925 1 : CATCH_REQUIRE(b[2][0] == b20);
2926 1 : CATCH_REQUIRE(b[2][1] == b21);
2927 1 : CATCH_REQUIRE(b[2][2] == b22);
2928 1 : CATCH_REQUIRE(b[2][3] == b23);
2929 1 : CATCH_REQUIRE(b[3][0] == b30);
2930 1 : CATCH_REQUIRE(b[3][1] == b31);
2931 1 : CATCH_REQUIRE(b[3][2] == b32);
2932 1 : CATCH_REQUIRE(b[3][3] == b33);
2933 :
2934 : // setup C
2935 : //
2936 2 : snapdev::matrix<double> c(4, 4);
2937 :
2938 1 : CATCH_REQUIRE(c.rows() == 4);
2939 1 : CATCH_REQUIRE(c.columns() == 4);
2940 :
2941 1 : double const c00(frand());
2942 1 : double const c01(frand());
2943 1 : double const c02(frand());
2944 1 : double const c03(frand());
2945 1 : double const c10(frand());
2946 1 : double const c11(frand());
2947 1 : double const c12(frand());
2948 1 : double const c13(frand());
2949 1 : double const c20(frand());
2950 1 : double const c21(frand());
2951 1 : double const c22(frand());
2952 1 : double const c23(frand());
2953 1 : double const c30(frand());
2954 1 : double const c31(frand());
2955 1 : double const c32(frand());
2956 1 : double const c33(frand());
2957 :
2958 1 : c[0][0] = c00;
2959 1 : c[0][1] = c01;
2960 1 : c[0][2] = c02;
2961 1 : c[0][3] = c03;
2962 1 : c[1][0] = c10;
2963 1 : c[1][1] = c11;
2964 1 : c[1][2] = c12;
2965 1 : c[1][3] = c13;
2966 1 : c[2][0] = c20;
2967 1 : c[2][1] = c21;
2968 1 : c[2][2] = c22;
2969 1 : c[2][3] = c23;
2970 1 : c[3][0] = c30;
2971 1 : c[3][1] = c31;
2972 1 : c[3][2] = c32;
2973 1 : c[3][3] = c33;
2974 :
2975 1 : CATCH_REQUIRE(c[0][0] == c00);
2976 1 : CATCH_REQUIRE(c[0][1] == c01);
2977 1 : CATCH_REQUIRE(c[0][2] == c02);
2978 1 : CATCH_REQUIRE(c[0][3] == c03);
2979 1 : CATCH_REQUIRE(c[1][0] == c10);
2980 1 : CATCH_REQUIRE(c[1][1] == c11);
2981 1 : CATCH_REQUIRE(c[1][2] == c12);
2982 1 : CATCH_REQUIRE(c[1][3] == c13);
2983 1 : CATCH_REQUIRE(c[2][0] == c20);
2984 1 : CATCH_REQUIRE(c[2][1] == c21);
2985 1 : CATCH_REQUIRE(c[2][2] == c22);
2986 1 : CATCH_REQUIRE(c[2][3] == c23);
2987 1 : CATCH_REQUIRE(c[3][0] == c30);
2988 1 : CATCH_REQUIRE(c[3][1] == c31);
2989 1 : CATCH_REQUIRE(c[3][2] == c32);
2990 1 : CATCH_REQUIRE(c[3][3] == c33);
2991 :
2992 : // run operation C = A * B
2993 : //
2994 1 : c = a * b;
2995 :
2996 1 : CATCH_REQUIRE(a[0][0] == a00);
2997 1 : CATCH_REQUIRE(a[0][1] == a01);
2998 1 : CATCH_REQUIRE(a[0][2] == a02);
2999 1 : CATCH_REQUIRE(a[0][3] == a03);
3000 1 : CATCH_REQUIRE(a[1][0] == a10);
3001 1 : CATCH_REQUIRE(a[1][1] == a11);
3002 1 : CATCH_REQUIRE(a[1][2] == a12);
3003 1 : CATCH_REQUIRE(a[1][3] == a13);
3004 1 : CATCH_REQUIRE(a[2][0] == a20);
3005 1 : CATCH_REQUIRE(a[2][1] == a21);
3006 1 : CATCH_REQUIRE(a[2][2] == a22);
3007 1 : CATCH_REQUIRE(a[2][3] == a23);
3008 1 : CATCH_REQUIRE(a[3][0] == a30);
3009 1 : CATCH_REQUIRE(a[3][1] == a31);
3010 1 : CATCH_REQUIRE(a[3][2] == a32);
3011 1 : CATCH_REQUIRE(a[3][3] == a33);
3012 :
3013 1 : CATCH_REQUIRE(b[0][0] == b00);
3014 1 : CATCH_REQUIRE(b[0][1] == b01);
3015 1 : CATCH_REQUIRE(b[0][2] == b02);
3016 1 : CATCH_REQUIRE(b[0][3] == b03);
3017 1 : CATCH_REQUIRE(b[1][0] == b10);
3018 1 : CATCH_REQUIRE(b[1][1] == b11);
3019 1 : CATCH_REQUIRE(b[1][2] == b12);
3020 1 : CATCH_REQUIRE(b[1][3] == b13);
3021 1 : CATCH_REQUIRE(b[2][0] == b20);
3022 1 : CATCH_REQUIRE(b[2][1] == b21);
3023 1 : CATCH_REQUIRE(b[2][2] == b22);
3024 1 : CATCH_REQUIRE(b[2][3] == b23);
3025 1 : CATCH_REQUIRE(b[3][0] == b30);
3026 1 : CATCH_REQUIRE(b[3][1] == b31);
3027 1 : CATCH_REQUIRE(b[3][2] == b32);
3028 1 : CATCH_REQUIRE(b[3][3] == b33);
3029 :
3030 : // this could fail (albeit rather unlikely)
3031 : //CATCH_REQUIRE_FALSE(c[0][0] == c00);
3032 : //CATCH_REQUIRE_FALSE(c[0][1] == c01);
3033 : //CATCH_REQUIRE_FALSE(c[0][2] == c02);
3034 : //CATCH_REQUIRE_FALSE(c[0][3] == c03);
3035 : //CATCH_REQUIRE_FALSE(c[1][0] == c10);
3036 : //CATCH_REQUIRE_FALSE(c[1][1] == c11);
3037 : //CATCH_REQUIRE_FALSE(c[1][2] == c12);
3038 : //CATCH_REQUIRE_FALSE(c[1][3] == c13);
3039 : //CATCH_REQUIRE_FALSE(c[2][0] == c20);
3040 : //CATCH_REQUIRE_FALSE(c[2][1] == c21);
3041 : //CATCH_REQUIRE_FALSE(c[2][2] == c22);
3042 : //CATCH_REQUIRE_FALSE(c[2][3] == c23);
3043 : //CATCH_REQUIRE_FALSE(c[3][0] == c30);
3044 : //CATCH_REQUIRE_FALSE(c[3][1] == c31);
3045 : //CATCH_REQUIRE_FALSE(c[3][2] == c32);
3046 : //CATCH_REQUIRE_FALSE(c[3][3] == c33);
3047 :
3048 : //std::cout << "A = " << a << std::endl;
3049 : //std::cout << "B = " << b << std::endl;
3050 : //std::cout << "C = " << c << std::endl;
3051 :
3052 1 : CATCH_REQUIRE(c[0][0] == a00 * b00 + a01 * b10 + a02 * b20 + a03 * b30);
3053 1 : CATCH_REQUIRE(c[0][1] == a00 * b01 + a01 * b11 + a02 * b21 + a03 * b31);
3054 1 : CATCH_REQUIRE(c[0][2] == a00 * b02 + a01 * b12 + a02 * b22 + a03 * b32);
3055 1 : CATCH_REQUIRE(c[0][3] == a00 * b03 + a01 * b13 + a02 * b23 + a03 * b33);
3056 :
3057 1 : CATCH_REQUIRE(c[1][0] == a10 * b00 + a11 * b10 + a12 * b20 + a13 * b30);
3058 1 : CATCH_REQUIRE(c[1][1] == a10 * b01 + a11 * b11 + a12 * b21 + a13 * b31);
3059 1 : CATCH_REQUIRE(c[1][2] == a10 * b02 + a11 * b12 + a12 * b22 + a13 * b32);
3060 1 : CATCH_REQUIRE(c[1][3] == a10 * b03 + a11 * b13 + a12 * b23 + a13 * b33);
3061 :
3062 1 : CATCH_REQUIRE(c[2][0] == a20 * b00 + a21 * b10 + a22 * b20 + a23 * b30);
3063 1 : CATCH_REQUIRE(c[2][1] == a20 * b01 + a21 * b11 + a22 * b21 + a23 * b31);
3064 1 : CATCH_REQUIRE(c[2][2] == a20 * b02 + a21 * b12 + a22 * b22 + a23 * b32);
3065 1 : CATCH_REQUIRE(c[2][3] == a20 * b03 + a21 * b13 + a22 * b23 + a23 * b33);
3066 :
3067 1 : CATCH_REQUIRE(c[3][0] == a30 * b00 + a31 * b10 + a32 * b20 + a33 * b30);
3068 1 : CATCH_REQUIRE(c[3][1] == a30 * b01 + a31 * b11 + a32 * b21 + a33 * b31);
3069 1 : CATCH_REQUIRE(c[3][2] == a30 * b02 + a31 * b12 + a32 * b22 + a33 * b32);
3070 1 : CATCH_REQUIRE(c[3][3] == a30 * b03 + a31 * b13 + a32 * b23 + a33 * b33);
3071 :
3072 : // if swapped it would equal this instead
3073 : //CATCH_REQUIRE(c[0][0] == a00 * b00 + a10 * b01 + a20 * b02 + a30 * b03);
3074 : //CATCH_REQUIRE(c[0][1] == a01 * b00 + a11 * b01 + a21 * b02 + a31 * b03);
3075 : //CATCH_REQUIRE(c[0][2] == a02 * b00 + a12 * b01 + a22 * b02 + a32 * b03);
3076 : //CATCH_REQUIRE(c[0][3] == a03 * b00 + a13 * b01 + a23 * b02 + a33 * b03);
3077 :
3078 : //CATCH_REQUIRE(c[1][0] == a00 * b10 + a10 * b11 + a20 * b12 + a30 * b13);
3079 : //CATCH_REQUIRE(c[1][1] == a01 * b10 + a11 * b11 + a21 * b12 + a31 * b13);
3080 : //CATCH_REQUIRE(c[1][2] == a02 * b10 + a12 * b11 + a22 * b12 + a32 * b13);
3081 : //CATCH_REQUIRE(c[1][3] == a03 * b10 + a13 * b11 + a23 * b12 + a33 * b13);
3082 :
3083 : //CATCH_REQUIRE(c[2][0] == a00 * b20 + a10 * b21 + a20 * b22 + a30 * b23);
3084 : //CATCH_REQUIRE(c[2][1] == a01 * b20 + a11 * b21 + a21 * b22 + a31 * b23);
3085 : //CATCH_REQUIRE(c[2][2] == a02 * b20 + a12 * b21 + a22 * b22 + a32 * b23);
3086 : //CATCH_REQUIRE(c[2][3] == a03 * b20 + a13 * b21 + a23 * b22 + a33 * b23);
3087 :
3088 : //CATCH_REQUIRE(c[3][0] == a00 * b30 + a10 * b31 + a20 * b32 + a30 * b33);
3089 : //CATCH_REQUIRE(c[3][1] == a01 * b30 + a11 * b31 + a21 * b32 + a31 * b33);
3090 : //CATCH_REQUIRE(c[3][2] == a02 * b30 + a12 * b31 + a22 * b32 + a32 * b33);
3091 : //CATCH_REQUIRE(c[3][3] == a03 * b30 + a13 * b31 + a23 * b32 + a33 * b33);
3092 : }
3093 : CATCH_END_SECTION()
3094 :
3095 8 : CATCH_START_SECTION("matrix: a*=b")
3096 : {
3097 : // setup A
3098 : //
3099 2 : snapdev::matrix<double> a(4, 4);
3100 :
3101 1 : CATCH_REQUIRE(a.rows() == 4);
3102 1 : CATCH_REQUIRE(a.columns() == 4);
3103 :
3104 1 : double const a00(frand());
3105 1 : double const a01(frand());
3106 1 : double const a02(frand());
3107 1 : double const a03(frand());
3108 1 : double const a10(frand());
3109 1 : double const a11(frand());
3110 1 : double const a12(frand());
3111 1 : double const a13(frand());
3112 1 : double const a20(frand());
3113 1 : double const a21(frand());
3114 1 : double const a22(frand());
3115 1 : double const a23(frand());
3116 1 : double const a30(frand());
3117 1 : double const a31(frand());
3118 1 : double const a32(frand());
3119 1 : double const a33(frand());
3120 :
3121 1 : a[0][0] = a00;
3122 1 : a[0][1] = a01;
3123 1 : a[0][2] = a02;
3124 1 : a[0][3] = a03;
3125 1 : a[1][0] = a10;
3126 1 : a[1][1] = a11;
3127 1 : a[1][2] = a12;
3128 1 : a[1][3] = a13;
3129 1 : a[2][0] = a20;
3130 1 : a[2][1] = a21;
3131 1 : a[2][2] = a22;
3132 1 : a[2][3] = a23;
3133 1 : a[3][0] = a30;
3134 1 : a[3][1] = a31;
3135 1 : a[3][2] = a32;
3136 1 : a[3][3] = a33;
3137 :
3138 1 : CATCH_REQUIRE(a[0][0] == a00);
3139 1 : CATCH_REQUIRE(a[0][1] == a01);
3140 1 : CATCH_REQUIRE(a[0][2] == a02);
3141 1 : CATCH_REQUIRE(a[0][3] == a03);
3142 1 : CATCH_REQUIRE(a[1][0] == a10);
3143 1 : CATCH_REQUIRE(a[1][1] == a11);
3144 1 : CATCH_REQUIRE(a[1][2] == a12);
3145 1 : CATCH_REQUIRE(a[1][3] == a13);
3146 1 : CATCH_REQUIRE(a[2][0] == a20);
3147 1 : CATCH_REQUIRE(a[2][1] == a21);
3148 1 : CATCH_REQUIRE(a[2][2] == a22);
3149 1 : CATCH_REQUIRE(a[2][3] == a23);
3150 1 : CATCH_REQUIRE(a[3][0] == a30);
3151 1 : CATCH_REQUIRE(a[3][1] == a31);
3152 1 : CATCH_REQUIRE(a[3][2] == a32);
3153 1 : CATCH_REQUIRE(a[3][3] == a33);
3154 :
3155 : // setup B
3156 : //
3157 2 : snapdev::matrix<double> b(4, 4);
3158 :
3159 1 : CATCH_REQUIRE(b.rows() == 4);
3160 1 : CATCH_REQUIRE(b.columns() == 4);
3161 :
3162 1 : double const b00(frand());
3163 1 : double const b01(frand());
3164 1 : double const b02(frand());
3165 1 : double const b03(frand());
3166 1 : double const b10(frand());
3167 1 : double const b11(frand());
3168 1 : double const b12(frand());
3169 1 : double const b13(frand());
3170 1 : double const b20(frand());
3171 1 : double const b21(frand());
3172 1 : double const b22(frand());
3173 1 : double const b23(frand());
3174 1 : double const b30(frand());
3175 1 : double const b31(frand());
3176 1 : double const b32(frand());
3177 1 : double const b33(frand());
3178 :
3179 1 : b[0][0] = b00;
3180 1 : b[0][1] = b01;
3181 1 : b[0][2] = b02;
3182 1 : b[0][3] = b03;
3183 1 : b[1][0] = b10;
3184 1 : b[1][1] = b11;
3185 1 : b[1][2] = b12;
3186 1 : b[1][3] = b13;
3187 1 : b[2][0] = b20;
3188 1 : b[2][1] = b21;
3189 1 : b[2][2] = b22;
3190 1 : b[2][3] = b23;
3191 1 : b[3][0] = b30;
3192 1 : b[3][1] = b31;
3193 1 : b[3][2] = b32;
3194 1 : b[3][3] = b33;
3195 :
3196 1 : CATCH_REQUIRE(b[0][0] == b00);
3197 1 : CATCH_REQUIRE(b[0][1] == b01);
3198 1 : CATCH_REQUIRE(b[0][2] == b02);
3199 1 : CATCH_REQUIRE(b[0][3] == b03);
3200 1 : CATCH_REQUIRE(b[1][0] == b10);
3201 1 : CATCH_REQUIRE(b[1][1] == b11);
3202 1 : CATCH_REQUIRE(b[1][2] == b12);
3203 1 : CATCH_REQUIRE(b[1][3] == b13);
3204 1 : CATCH_REQUIRE(b[2][0] == b20);
3205 1 : CATCH_REQUIRE(b[2][1] == b21);
3206 1 : CATCH_REQUIRE(b[2][2] == b22);
3207 1 : CATCH_REQUIRE(b[2][3] == b23);
3208 1 : CATCH_REQUIRE(b[3][0] == b30);
3209 1 : CATCH_REQUIRE(b[3][1] == b31);
3210 1 : CATCH_REQUIRE(b[3][2] == b32);
3211 1 : CATCH_REQUIRE(b[3][3] == b33);
3212 :
3213 : // run operation A *= B
3214 : //
3215 1 : a *= b;
3216 :
3217 : // this could fail
3218 : //
3219 : //CATCH_REQUIRE_FALSE(a[0][0] == a00);
3220 : //CATCH_REQUIRE_FALSE(a[0][1] == a01);
3221 : //CATCH_REQUIRE_FALSE(a[0][2] == a02);
3222 : //CATCH_REQUIRE_FALSE(a[0][3] == a03);
3223 : //CATCH_REQUIRE_FALSE(a[1][0] == a10);
3224 : //CATCH_REQUIRE_FALSE(a[1][1] == a11);
3225 : //CATCH_REQUIRE_FALSE(a[1][2] == a12);
3226 : //CATCH_REQUIRE_FALSE(a[1][3] == a13);
3227 : //CATCH_REQUIRE_FALSE(a[2][0] == a20);
3228 : //CATCH_REQUIRE_FALSE(a[2][1] == a21);
3229 : //CATCH_REQUIRE_FALSE(a[2][2] == a22);
3230 : //CATCH_REQUIRE_FALSE(a[2][3] == a23);
3231 : //CATCH_REQUIRE_FALSE(a[3][0] == a30);
3232 : //CATCH_REQUIRE_FALSE(a[3][1] == a31);
3233 : //CATCH_REQUIRE_FALSE(a[3][2] == a32);
3234 : //CATCH_REQUIRE_FALSE(a[3][3] == a33);
3235 :
3236 1 : CATCH_REQUIRE(b[0][0] == b00);
3237 1 : CATCH_REQUIRE(b[0][1] == b01);
3238 1 : CATCH_REQUIRE(b[0][2] == b02);
3239 1 : CATCH_REQUIRE(b[0][3] == b03);
3240 1 : CATCH_REQUIRE(b[1][0] == b10);
3241 1 : CATCH_REQUIRE(b[1][1] == b11);
3242 1 : CATCH_REQUIRE(b[1][2] == b12);
3243 1 : CATCH_REQUIRE(b[1][3] == b13);
3244 1 : CATCH_REQUIRE(b[2][0] == b20);
3245 1 : CATCH_REQUIRE(b[2][1] == b21);
3246 1 : CATCH_REQUIRE(b[2][2] == b22);
3247 1 : CATCH_REQUIRE(b[2][3] == b23);
3248 1 : CATCH_REQUIRE(b[3][0] == b30);
3249 1 : CATCH_REQUIRE(b[3][1] == b31);
3250 1 : CATCH_REQUIRE(b[3][2] == b32);
3251 1 : CATCH_REQUIRE(b[3][3] == b33);
3252 :
3253 1 : CATCH_REQUIRE(a[0][0] == a00 * b00 + a01 * b10 + a02 * b20 + a03 * b30);
3254 1 : CATCH_REQUIRE(a[0][1] == a00 * b01 + a01 * b11 + a02 * b21 + a03 * b31);
3255 1 : CATCH_REQUIRE(a[0][2] == a00 * b02 + a01 * b12 + a02 * b22 + a03 * b32);
3256 1 : CATCH_REQUIRE(a[0][3] == a00 * b03 + a01 * b13 + a02 * b23 + a03 * b33);
3257 :
3258 1 : CATCH_REQUIRE(a[1][0] == a10 * b00 + a11 * b10 + a12 * b20 + a13 * b30);
3259 1 : CATCH_REQUIRE(a[1][1] == a10 * b01 + a11 * b11 + a12 * b21 + a13 * b31);
3260 1 : CATCH_REQUIRE(a[1][2] == a10 * b02 + a11 * b12 + a12 * b22 + a13 * b32);
3261 1 : CATCH_REQUIRE(a[1][3] == a10 * b03 + a11 * b13 + a12 * b23 + a13 * b33);
3262 :
3263 1 : CATCH_REQUIRE(a[2][0] == a20 * b00 + a21 * b10 + a22 * b20 + a23 * b30);
3264 1 : CATCH_REQUIRE(a[2][1] == a20 * b01 + a21 * b11 + a22 * b21 + a23 * b31);
3265 1 : CATCH_REQUIRE(a[2][2] == a20 * b02 + a21 * b12 + a22 * b22 + a23 * b32);
3266 1 : CATCH_REQUIRE(a[2][3] == a20 * b03 + a21 * b13 + a22 * b23 + a23 * b33);
3267 :
3268 1 : CATCH_REQUIRE(a[3][0] == a30 * b00 + a31 * b10 + a32 * b20 + a33 * b30);
3269 1 : CATCH_REQUIRE(a[3][1] == a30 * b01 + a31 * b11 + a32 * b21 + a33 * b31);
3270 1 : CATCH_REQUIRE(a[3][2] == a30 * b02 + a31 * b12 + a32 * b22 + a33 * b32);
3271 1 : CATCH_REQUIRE(a[3][3] == a30 * b03 + a31 * b13 + a32 * b23 + a33 * b33);
3272 :
3273 : //CATCH_REQUIRE(a[0][0] == a00 * b00 + a10 * b01 + a20 * b02 + a30 * b03);
3274 : //CATCH_REQUIRE(a[0][1] == a01 * b00 + a11 * b01 + a21 * b02 + a31 * b03);
3275 : //CATCH_REQUIRE(a[0][2] == a02 * b00 + a12 * b01 + a22 * b02 + a32 * b03);
3276 : //CATCH_REQUIRE(a[0][3] == a03 * b00 + a13 * b01 + a23 * b02 + a33 * b03);
3277 :
3278 : //CATCH_REQUIRE(a[1][0] == a00 * b10 + a10 * b11 + a20 * b12 + a30 * b13);
3279 : //CATCH_REQUIRE(a[1][1] == a01 * b10 + a11 * b11 + a21 * b12 + a31 * b13);
3280 : //CATCH_REQUIRE(a[1][2] == a02 * b10 + a12 * b11 + a22 * b12 + a32 * b13);
3281 : //CATCH_REQUIRE(a[1][3] == a03 * b10 + a13 * b11 + a23 * b12 + a33 * b13);
3282 :
3283 : //CATCH_REQUIRE(a[2][0] == a00 * b20 + a10 * b21 + a20 * b22 + a30 * b23);
3284 : //CATCH_REQUIRE(a[2][1] == a01 * b20 + a11 * b21 + a21 * b22 + a31 * b23);
3285 : //CATCH_REQUIRE(a[2][2] == a02 * b20 + a12 * b21 + a22 * b22 + a32 * b23);
3286 : //CATCH_REQUIRE(a[2][3] == a03 * b20 + a13 * b21 + a23 * b22 + a33 * b23);
3287 :
3288 : //CATCH_REQUIRE(a[3][0] == a00 * b30 + a10 * b31 + a20 * b32 + a30 * b33);
3289 : //CATCH_REQUIRE(a[3][1] == a01 * b30 + a11 * b31 + a21 * b32 + a31 * b33);
3290 : //CATCH_REQUIRE(a[3][2] == a02 * b30 + a12 * b31 + a22 * b32 + a32 * b33);
3291 : //CATCH_REQUIRE(a[3][3] == a03 * b30 + a13 * b31 + a23 * b32 + a33 * b33);
3292 : }
3293 : CATCH_END_SECTION()
3294 : }
3295 :
3296 : // create two random 4x4 matrices and make sure the add works
3297 : //
3298 16 : CATCH_GIVEN("div")
3299 : {
3300 8 : CATCH_START_SECTION("matrix: b=a/<scalar>")
3301 : {
3302 : // setup A
3303 : //
3304 2 : snapdev::matrix<double> a(4, 4);
3305 :
3306 1 : CATCH_REQUIRE(a.rows() == 4);
3307 1 : CATCH_REQUIRE(a.columns() == 4);
3308 :
3309 1 : double const a00(frand());
3310 1 : double const a01(frand());
3311 1 : double const a02(frand());
3312 1 : double const a03(frand());
3313 1 : double const a10(frand());
3314 1 : double const a11(frand());
3315 1 : double const a12(frand());
3316 1 : double const a13(frand());
3317 1 : double const a20(frand());
3318 1 : double const a21(frand());
3319 1 : double const a22(frand());
3320 1 : double const a23(frand());
3321 1 : double const a30(frand());
3322 1 : double const a31(frand());
3323 1 : double const a32(frand());
3324 1 : double const a33(frand());
3325 :
3326 1 : a[0][0] = a00;
3327 1 : a[0][1] = a01;
3328 1 : a[0][2] = a02;
3329 1 : a[0][3] = a03;
3330 1 : a[1][0] = a10;
3331 1 : a[1][1] = a11;
3332 1 : a[1][2] = a12;
3333 1 : a[1][3] = a13;
3334 1 : a[2][0] = a20;
3335 1 : a[2][1] = a21;
3336 1 : a[2][2] = a22;
3337 1 : a[2][3] = a23;
3338 1 : a[3][0] = a30;
3339 1 : a[3][1] = a31;
3340 1 : a[3][2] = a32;
3341 1 : a[3][3] = a33;
3342 :
3343 1 : CATCH_REQUIRE(a[0][0] == a00);
3344 1 : CATCH_REQUIRE(a[0][1] == a01);
3345 1 : CATCH_REQUIRE(a[0][2] == a02);
3346 1 : CATCH_REQUIRE(a[0][3] == a03);
3347 1 : CATCH_REQUIRE(a[1][0] == a10);
3348 1 : CATCH_REQUIRE(a[1][1] == a11);
3349 1 : CATCH_REQUIRE(a[1][2] == a12);
3350 1 : CATCH_REQUIRE(a[1][3] == a13);
3351 1 : CATCH_REQUIRE(a[2][0] == a20);
3352 1 : CATCH_REQUIRE(a[2][1] == a21);
3353 1 : CATCH_REQUIRE(a[2][2] == a22);
3354 1 : CATCH_REQUIRE(a[2][3] == a23);
3355 1 : CATCH_REQUIRE(a[3][0] == a30);
3356 1 : CATCH_REQUIRE(a[3][1] == a31);
3357 1 : CATCH_REQUIRE(a[3][2] == a32);
3358 1 : CATCH_REQUIRE(a[3][3] == a33);
3359 :
3360 : // setup B
3361 : //
3362 2 : snapdev::matrix<double> b(4, 4);
3363 :
3364 1 : CATCH_REQUIRE(b.rows() == 4);
3365 1 : CATCH_REQUIRE(b.columns() == 4);
3366 :
3367 1 : double const b00(frand());
3368 1 : double const b01(frand());
3369 1 : double const b02(frand());
3370 1 : double const b03(frand());
3371 1 : double const b10(frand());
3372 1 : double const b11(frand());
3373 1 : double const b12(frand());
3374 1 : double const b13(frand());
3375 1 : double const b20(frand());
3376 1 : double const b21(frand());
3377 1 : double const b22(frand());
3378 1 : double const b23(frand());
3379 1 : double const b30(frand());
3380 1 : double const b31(frand());
3381 1 : double const b32(frand());
3382 1 : double const b33(frand());
3383 :
3384 1 : b[0][0] = b00;
3385 1 : b[0][1] = b01;
3386 1 : b[0][2] = b02;
3387 1 : b[0][3] = b03;
3388 1 : b[1][0] = b10;
3389 1 : b[1][1] = b11;
3390 1 : b[1][2] = b12;
3391 1 : b[1][3] = b13;
3392 1 : b[2][0] = b20;
3393 1 : b[2][1] = b21;
3394 1 : b[2][2] = b22;
3395 1 : b[2][3] = b23;
3396 1 : b[3][0] = b30;
3397 1 : b[3][1] = b31;
3398 1 : b[3][2] = b32;
3399 1 : b[3][3] = b33;
3400 :
3401 1 : CATCH_REQUIRE(b[0][0] == b00);
3402 1 : CATCH_REQUIRE(b[0][1] == b01);
3403 1 : CATCH_REQUIRE(b[0][2] == b02);
3404 1 : CATCH_REQUIRE(b[0][3] == b03);
3405 1 : CATCH_REQUIRE(b[1][0] == b10);
3406 1 : CATCH_REQUIRE(b[1][1] == b11);
3407 1 : CATCH_REQUIRE(b[1][2] == b12);
3408 1 : CATCH_REQUIRE(b[1][3] == b13);
3409 1 : CATCH_REQUIRE(b[2][0] == b20);
3410 1 : CATCH_REQUIRE(b[2][1] == b21);
3411 1 : CATCH_REQUIRE(b[2][2] == b22);
3412 1 : CATCH_REQUIRE(b[2][3] == b23);
3413 1 : CATCH_REQUIRE(b[3][0] == b30);
3414 1 : CATCH_REQUIRE(b[3][1] == b31);
3415 1 : CATCH_REQUIRE(b[3][2] == b32);
3416 1 : CATCH_REQUIRE(b[3][3] == b33);
3417 :
3418 : // create a scalar for our test
3419 : //
3420 1 : double scalar(frand());
3421 1 : while(scalar == 0.0 || scalar == 1.0)
3422 : {
3423 0 : scalar = frand();
3424 : }
3425 :
3426 : // run operation B = A / <scalar>
3427 : //
3428 1 : b = a / scalar;
3429 :
3430 1 : CATCH_REQUIRE(a[0][0] == a00);
3431 1 : CATCH_REQUIRE(a[0][1] == a01);
3432 1 : CATCH_REQUIRE(a[0][2] == a02);
3433 1 : CATCH_REQUIRE(a[0][3] == a03);
3434 1 : CATCH_REQUIRE(a[1][0] == a10);
3435 1 : CATCH_REQUIRE(a[1][1] == a11);
3436 1 : CATCH_REQUIRE(a[1][2] == a12);
3437 1 : CATCH_REQUIRE(a[1][3] == a13);
3438 1 : CATCH_REQUIRE(a[2][0] == a20);
3439 1 : CATCH_REQUIRE(a[2][1] == a21);
3440 1 : CATCH_REQUIRE(a[2][2] == a22);
3441 1 : CATCH_REQUIRE(a[2][3] == a23);
3442 1 : CATCH_REQUIRE(a[3][0] == a30);
3443 1 : CATCH_REQUIRE(a[3][1] == a31);
3444 1 : CATCH_REQUIRE(a[3][2] == a32);
3445 1 : CATCH_REQUIRE(a[3][3] == a33);
3446 :
3447 : // this can't fail because we ensure scalar != 0 or 1
3448 : //
3449 1 : CATCH_REQUIRE_FALSE(b[0][0] == b00);
3450 1 : CATCH_REQUIRE_FALSE(b[0][1] == b01);
3451 1 : CATCH_REQUIRE_FALSE(b[0][2] == b02);
3452 1 : CATCH_REQUIRE_FALSE(b[0][3] == b03);
3453 1 : CATCH_REQUIRE_FALSE(b[1][0] == b10);
3454 1 : CATCH_REQUIRE_FALSE(b[1][1] == b11);
3455 1 : CATCH_REQUIRE_FALSE(b[1][2] == b12);
3456 1 : CATCH_REQUIRE_FALSE(b[1][3] == b13);
3457 1 : CATCH_REQUIRE_FALSE(b[2][0] == b20);
3458 1 : CATCH_REQUIRE_FALSE(b[2][1] == b21);
3459 1 : CATCH_REQUIRE_FALSE(b[2][2] == b22);
3460 1 : CATCH_REQUIRE_FALSE(b[2][3] == b23);
3461 1 : CATCH_REQUIRE_FALSE(b[3][0] == b30);
3462 1 : CATCH_REQUIRE_FALSE(b[3][1] == b31);
3463 1 : CATCH_REQUIRE_FALSE(b[3][2] == b32);
3464 1 : CATCH_REQUIRE_FALSE(b[3][3] == b33);
3465 :
3466 1 : CATCH_REQUIRE(b[0][0] == a00 / scalar);
3467 1 : CATCH_REQUIRE(b[0][1] == a01 / scalar);
3468 1 : CATCH_REQUIRE(b[0][2] == a02 / scalar);
3469 1 : CATCH_REQUIRE(b[0][3] == a03 / scalar);
3470 1 : CATCH_REQUIRE(b[1][0] == a10 / scalar);
3471 1 : CATCH_REQUIRE(b[1][1] == a11 / scalar);
3472 1 : CATCH_REQUIRE(b[1][2] == a12 / scalar);
3473 1 : CATCH_REQUIRE(b[1][3] == a13 / scalar);
3474 1 : CATCH_REQUIRE(b[2][0] == a20 / scalar);
3475 1 : CATCH_REQUIRE(b[2][1] == a21 / scalar);
3476 1 : CATCH_REQUIRE(b[2][2] == a22 / scalar);
3477 1 : CATCH_REQUIRE(b[2][3] == a23 / scalar);
3478 1 : CATCH_REQUIRE(b[3][0] == a30 / scalar);
3479 1 : CATCH_REQUIRE(b[3][1] == a31 / scalar);
3480 1 : CATCH_REQUIRE(b[3][2] == a32 / scalar);
3481 1 : CATCH_REQUIRE(b[3][3] == a33 / scalar);
3482 : }
3483 : CATCH_END_SECTION()
3484 :
3485 8 : CATCH_START_SECTION("matrix: a/=<scalar>")
3486 : {
3487 : // setup A
3488 : //
3489 2 : snapdev::matrix<double> a(4, 4);
3490 :
3491 1 : CATCH_REQUIRE(a.rows() == 4);
3492 1 : CATCH_REQUIRE(a.columns() == 4);
3493 :
3494 1 : double const a00(frand());
3495 1 : double const a01(frand());
3496 1 : double const a02(frand());
3497 1 : double const a03(frand());
3498 1 : double const a10(frand());
3499 1 : double const a11(frand());
3500 1 : double const a12(frand());
3501 1 : double const a13(frand());
3502 1 : double const a20(frand());
3503 1 : double const a21(frand());
3504 1 : double const a22(frand());
3505 1 : double const a23(frand());
3506 1 : double const a30(frand());
3507 1 : double const a31(frand());
3508 1 : double const a32(frand());
3509 1 : double const a33(frand());
3510 :
3511 1 : a[0][0] = a00;
3512 1 : a[0][1] = a01;
3513 1 : a[0][2] = a02;
3514 1 : a[0][3] = a03;
3515 1 : a[1][0] = a10;
3516 1 : a[1][1] = a11;
3517 1 : a[1][2] = a12;
3518 1 : a[1][3] = a13;
3519 1 : a[2][0] = a20;
3520 1 : a[2][1] = a21;
3521 1 : a[2][2] = a22;
3522 1 : a[2][3] = a23;
3523 1 : a[3][0] = a30;
3524 1 : a[3][1] = a31;
3525 1 : a[3][2] = a32;
3526 1 : a[3][3] = a33;
3527 :
3528 1 : CATCH_REQUIRE(a[0][0] == a00);
3529 1 : CATCH_REQUIRE(a[0][1] == a01);
3530 1 : CATCH_REQUIRE(a[0][2] == a02);
3531 1 : CATCH_REQUIRE(a[0][3] == a03);
3532 1 : CATCH_REQUIRE(a[1][0] == a10);
3533 1 : CATCH_REQUIRE(a[1][1] == a11);
3534 1 : CATCH_REQUIRE(a[1][2] == a12);
3535 1 : CATCH_REQUIRE(a[1][3] == a13);
3536 1 : CATCH_REQUIRE(a[2][0] == a20);
3537 1 : CATCH_REQUIRE(a[2][1] == a21);
3538 1 : CATCH_REQUIRE(a[2][2] == a22);
3539 1 : CATCH_REQUIRE(a[2][3] == a23);
3540 1 : CATCH_REQUIRE(a[3][0] == a30);
3541 1 : CATCH_REQUIRE(a[3][1] == a31);
3542 1 : CATCH_REQUIRE(a[3][2] == a32);
3543 1 : CATCH_REQUIRE(a[3][3] == a33);
3544 :
3545 : // create a scalar for our test
3546 : //
3547 1 : double scalar(frand());
3548 1 : while(scalar == 0.0)
3549 : {
3550 0 : scalar = frand();
3551 : }
3552 :
3553 : // run operation A /= <scalar>
3554 : //
3555 1 : a /= scalar;
3556 :
3557 : // this can't fail because we ensure scalar != 0 or 1
3558 : //
3559 1 : CATCH_REQUIRE_FALSE(a[0][0] == a00);
3560 1 : CATCH_REQUIRE_FALSE(a[0][1] == a01);
3561 1 : CATCH_REQUIRE_FALSE(a[0][2] == a02);
3562 1 : CATCH_REQUIRE_FALSE(a[0][3] == a03);
3563 1 : CATCH_REQUIRE_FALSE(a[1][0] == a10);
3564 1 : CATCH_REQUIRE_FALSE(a[1][1] == a11);
3565 1 : CATCH_REQUIRE_FALSE(a[1][2] == a12);
3566 1 : CATCH_REQUIRE_FALSE(a[1][3] == a13);
3567 1 : CATCH_REQUIRE_FALSE(a[2][0] == a20);
3568 1 : CATCH_REQUIRE_FALSE(a[2][1] == a21);
3569 1 : CATCH_REQUIRE_FALSE(a[2][2] == a22);
3570 1 : CATCH_REQUIRE_FALSE(a[2][3] == a23);
3571 1 : CATCH_REQUIRE_FALSE(a[3][0] == a30);
3572 1 : CATCH_REQUIRE_FALSE(a[3][1] == a31);
3573 1 : CATCH_REQUIRE_FALSE(a[3][2] == a32);
3574 1 : CATCH_REQUIRE_FALSE(a[3][3] == a33);
3575 :
3576 1 : CATCH_REQUIRE(a[0][0] == a00 / scalar);
3577 1 : CATCH_REQUIRE(a[0][1] == a01 / scalar);
3578 1 : CATCH_REQUIRE(a[0][2] == a02 / scalar);
3579 1 : CATCH_REQUIRE(a[0][3] == a03 / scalar);
3580 1 : CATCH_REQUIRE(a[1][0] == a10 / scalar);
3581 1 : CATCH_REQUIRE(a[1][1] == a11 / scalar);
3582 1 : CATCH_REQUIRE(a[1][2] == a12 / scalar);
3583 1 : CATCH_REQUIRE(a[1][3] == a13 / scalar);
3584 1 : CATCH_REQUIRE(a[2][0] == a20 / scalar);
3585 1 : CATCH_REQUIRE(a[2][1] == a21 / scalar);
3586 1 : CATCH_REQUIRE(a[2][2] == a22 / scalar);
3587 1 : CATCH_REQUIRE(a[2][3] == a23 / scalar);
3588 1 : CATCH_REQUIRE(a[3][0] == a30 / scalar);
3589 1 : CATCH_REQUIRE(a[3][1] == a31 / scalar);
3590 1 : CATCH_REQUIRE(a[3][2] == a32 / scalar);
3591 1 : CATCH_REQUIRE(a[3][3] == a33 / scalar);
3592 : }
3593 : CATCH_END_SECTION()
3594 :
3595 8 : CATCH_START_SECTION("matrix: c=a/b")
3596 : {
3597 : // setup A
3598 : //
3599 2 : snapdev::matrix<double> a(4, 4);
3600 :
3601 1 : CATCH_REQUIRE(a.rows() == 4);
3602 1 : CATCH_REQUIRE(a.columns() == 4);
3603 :
3604 1 : double const a00(frand());
3605 1 : double const a01(frand());
3606 1 : double const a02(frand());
3607 1 : double const a03(frand());
3608 1 : double const a10(frand());
3609 1 : double const a11(frand());
3610 1 : double const a12(frand());
3611 1 : double const a13(frand());
3612 1 : double const a20(frand());
3613 1 : double const a21(frand());
3614 1 : double const a22(frand());
3615 1 : double const a23(frand());
3616 1 : double const a30(frand());
3617 1 : double const a31(frand());
3618 1 : double const a32(frand());
3619 1 : double const a33(frand());
3620 :
3621 1 : a[0][0] = a00;
3622 1 : a[0][1] = a01;
3623 1 : a[0][2] = a02;
3624 1 : a[0][3] = a03;
3625 1 : a[1][0] = a10;
3626 1 : a[1][1] = a11;
3627 1 : a[1][2] = a12;
3628 1 : a[1][3] = a13;
3629 1 : a[2][0] = a20;
3630 1 : a[2][1] = a21;
3631 1 : a[2][2] = a22;
3632 1 : a[2][3] = a23;
3633 1 : a[3][0] = a30;
3634 1 : a[3][1] = a31;
3635 1 : a[3][2] = a32;
3636 1 : a[3][3] = a33;
3637 :
3638 1 : CATCH_REQUIRE(a[0][0] == a00);
3639 1 : CATCH_REQUIRE(a[0][1] == a01);
3640 1 : CATCH_REQUIRE(a[0][2] == a02);
3641 1 : CATCH_REQUIRE(a[0][3] == a03);
3642 1 : CATCH_REQUIRE(a[1][0] == a10);
3643 1 : CATCH_REQUIRE(a[1][1] == a11);
3644 1 : CATCH_REQUIRE(a[1][2] == a12);
3645 1 : CATCH_REQUIRE(a[1][3] == a13);
3646 1 : CATCH_REQUIRE(a[2][0] == a20);
3647 1 : CATCH_REQUIRE(a[2][1] == a21);
3648 1 : CATCH_REQUIRE(a[2][2] == a22);
3649 1 : CATCH_REQUIRE(a[2][3] == a23);
3650 1 : CATCH_REQUIRE(a[3][0] == a30);
3651 1 : CATCH_REQUIRE(a[3][1] == a31);
3652 1 : CATCH_REQUIRE(a[3][2] == a32);
3653 1 : CATCH_REQUIRE(a[3][3] == a33);
3654 :
3655 : // setup B
3656 : //
3657 2 : snapdev::matrix<double> b(4, 4);
3658 :
3659 1 : CATCH_REQUIRE(b.rows() == 4);
3660 1 : CATCH_REQUIRE(b.columns() == 4);
3661 :
3662 : double b00;
3663 : double b01;
3664 : double b02;
3665 : double b03;
3666 : double b10;
3667 : double b11;
3668 : double b12;
3669 : double b13;
3670 : double b20;
3671 : double b21;
3672 : double b22;
3673 : double b23;
3674 : double b30;
3675 : double b31;
3676 : double b32;
3677 : double b33;
3678 :
3679 : // create an inversible matrix (most are with random numbers)
3680 2 : snapdev::matrix<double> t(4, 4);
3681 0 : do
3682 : {
3683 1 : b00 = frand();
3684 1 : b01 = frand();
3685 1 : b02 = frand();
3686 1 : b03 = frand();
3687 1 : b10 = frand();
3688 1 : b11 = frand();
3689 1 : b12 = frand();
3690 1 : b13 = frand();
3691 1 : b20 = frand();
3692 1 : b21 = frand();
3693 1 : b22 = frand();
3694 1 : b23 = frand();
3695 1 : b30 = frand();
3696 1 : b31 = frand();
3697 1 : b32 = frand();
3698 1 : b33 = frand();
3699 :
3700 1 : b[0][0] = b00;
3701 1 : b[0][1] = b01;
3702 1 : b[0][2] = b02;
3703 1 : b[0][3] = b03;
3704 1 : b[1][0] = b10;
3705 1 : b[1][1] = b11;
3706 1 : b[1][2] = b12;
3707 1 : b[1][3] = b13;
3708 1 : b[2][0] = b20;
3709 1 : b[2][1] = b21;
3710 1 : b[2][2] = b22;
3711 1 : b[2][3] = b23;
3712 1 : b[3][0] = b30;
3713 1 : b[3][1] = b31;
3714 1 : b[3][2] = b32;
3715 1 : b[3][3] = b33;
3716 :
3717 1 : CATCH_REQUIRE(b[0][0] == b00);
3718 1 : CATCH_REQUIRE(b[0][1] == b01);
3719 1 : CATCH_REQUIRE(b[0][2] == b02);
3720 1 : CATCH_REQUIRE(b[0][3] == b03);
3721 1 : CATCH_REQUIRE(b[1][0] == b10);
3722 1 : CATCH_REQUIRE(b[1][1] == b11);
3723 1 : CATCH_REQUIRE(b[1][2] == b12);
3724 1 : CATCH_REQUIRE(b[1][3] == b13);
3725 1 : CATCH_REQUIRE(b[2][0] == b20);
3726 1 : CATCH_REQUIRE(b[2][1] == b21);
3727 1 : CATCH_REQUIRE(b[2][2] == b22);
3728 1 : CATCH_REQUIRE(b[2][3] == b23);
3729 1 : CATCH_REQUIRE(b[3][0] == b30);
3730 1 : CATCH_REQUIRE(b[3][1] == b31);
3731 1 : CATCH_REQUIRE(b[3][2] == b32);
3732 1 : CATCH_REQUIRE(b[3][3] == b33);
3733 :
3734 1 : t = b;
3735 : }
3736 1 : while(!t.inverse());
3737 :
3738 : // setup C
3739 : //
3740 2 : snapdev::matrix<double> c(4, 4);
3741 :
3742 1 : CATCH_REQUIRE(c.rows() == 4);
3743 1 : CATCH_REQUIRE(c.columns() == 4);
3744 :
3745 1 : double const c00(frand());
3746 1 : double const c01(frand());
3747 1 : double const c02(frand());
3748 1 : double const c03(frand());
3749 1 : double const c10(frand());
3750 1 : double const c11(frand());
3751 1 : double const c12(frand());
3752 1 : double const c13(frand());
3753 1 : double const c20(frand());
3754 1 : double const c21(frand());
3755 1 : double const c22(frand());
3756 1 : double const c23(frand());
3757 1 : double const c30(frand());
3758 1 : double const c31(frand());
3759 1 : double const c32(frand());
3760 1 : double const c33(frand());
3761 :
3762 1 : c[0][0] = c00;
3763 1 : c[0][1] = c01;
3764 1 : c[0][2] = c02;
3765 1 : c[0][3] = c03;
3766 1 : c[1][0] = c10;
3767 1 : c[1][1] = c11;
3768 1 : c[1][2] = c12;
3769 1 : c[1][3] = c13;
3770 1 : c[2][0] = c20;
3771 1 : c[2][1] = c21;
3772 1 : c[2][2] = c22;
3773 1 : c[2][3] = c23;
3774 1 : c[3][0] = c30;
3775 1 : c[3][1] = c31;
3776 1 : c[3][2] = c32;
3777 1 : c[3][3] = c33;
3778 :
3779 1 : CATCH_REQUIRE(c[0][0] == c00);
3780 1 : CATCH_REQUIRE(c[0][1] == c01);
3781 1 : CATCH_REQUIRE(c[0][2] == c02);
3782 1 : CATCH_REQUIRE(c[0][3] == c03);
3783 1 : CATCH_REQUIRE(c[1][0] == c10);
3784 1 : CATCH_REQUIRE(c[1][1] == c11);
3785 1 : CATCH_REQUIRE(c[1][2] == c12);
3786 1 : CATCH_REQUIRE(c[1][3] == c13);
3787 1 : CATCH_REQUIRE(c[2][0] == c20);
3788 1 : CATCH_REQUIRE(c[2][1] == c21);
3789 1 : CATCH_REQUIRE(c[2][2] == c22);
3790 1 : CATCH_REQUIRE(c[2][3] == c23);
3791 1 : CATCH_REQUIRE(c[3][0] == c30);
3792 1 : CATCH_REQUIRE(c[3][1] == c31);
3793 1 : CATCH_REQUIRE(c[3][2] == c32);
3794 1 : CATCH_REQUIRE(c[3][3] == c33);
3795 :
3796 : // run operation C = A / B
3797 : //
3798 1 : c = a / b;
3799 :
3800 1 : CATCH_REQUIRE(a[0][0] == a00);
3801 1 : CATCH_REQUIRE(a[0][1] == a01);
3802 1 : CATCH_REQUIRE(a[0][2] == a02);
3803 1 : CATCH_REQUIRE(a[0][3] == a03);
3804 1 : CATCH_REQUIRE(a[1][0] == a10);
3805 1 : CATCH_REQUIRE(a[1][1] == a11);
3806 1 : CATCH_REQUIRE(a[1][2] == a12);
3807 1 : CATCH_REQUIRE(a[1][3] == a13);
3808 1 : CATCH_REQUIRE(a[2][0] == a20);
3809 1 : CATCH_REQUIRE(a[2][1] == a21);
3810 1 : CATCH_REQUIRE(a[2][2] == a22);
3811 1 : CATCH_REQUIRE(a[2][3] == a23);
3812 1 : CATCH_REQUIRE(a[3][0] == a30);
3813 1 : CATCH_REQUIRE(a[3][1] == a31);
3814 1 : CATCH_REQUIRE(a[3][2] == a32);
3815 1 : CATCH_REQUIRE(a[3][3] == a33);
3816 :
3817 1 : CATCH_REQUIRE(b[0][0] == b00);
3818 1 : CATCH_REQUIRE(b[0][1] == b01);
3819 1 : CATCH_REQUIRE(b[0][2] == b02);
3820 1 : CATCH_REQUIRE(b[0][3] == b03);
3821 1 : CATCH_REQUIRE(b[1][0] == b10);
3822 1 : CATCH_REQUIRE(b[1][1] == b11);
3823 1 : CATCH_REQUIRE(b[1][2] == b12);
3824 1 : CATCH_REQUIRE(b[1][3] == b13);
3825 1 : CATCH_REQUIRE(b[2][0] == b20);
3826 1 : CATCH_REQUIRE(b[2][1] == b21);
3827 1 : CATCH_REQUIRE(b[2][2] == b22);
3828 1 : CATCH_REQUIRE(b[2][3] == b23);
3829 1 : CATCH_REQUIRE(b[3][0] == b30);
3830 1 : CATCH_REQUIRE(b[3][1] == b31);
3831 1 : CATCH_REQUIRE(b[3][2] == b32);
3832 1 : CATCH_REQUIRE(b[3][3] == b33);
3833 :
3834 : // this could fail
3835 : //CATCH_REQUIRE_FALSE(c[0][0] == c00);
3836 : //CATCH_REQUIRE_FALSE(c[0][1] == c01);
3837 : //CATCH_REQUIRE_FALSE(c[0][2] == c02);
3838 : //CATCH_REQUIRE_FALSE(c[0][3] == c03);
3839 : //CATCH_REQUIRE_FALSE(c[1][0] == c10);
3840 : //CATCH_REQUIRE_FALSE(c[1][1] == c11);
3841 : //CATCH_REQUIRE_FALSE(c[1][2] == c12);
3842 : //CATCH_REQUIRE_FALSE(c[1][3] == c13);
3843 : //CATCH_REQUIRE_FALSE(c[2][0] == c20);
3844 : //CATCH_REQUIRE_FALSE(c[2][1] == c21);
3845 : //CATCH_REQUIRE_FALSE(c[2][2] == c22);
3846 : //CATCH_REQUIRE_FALSE(c[2][3] == c23);
3847 : //CATCH_REQUIRE_FALSE(c[3][0] == c30);
3848 : //CATCH_REQUIRE_FALSE(c[3][1] == c31);
3849 : //CATCH_REQUIRE_FALSE(c[3][2] == c32);
3850 : //CATCH_REQUIRE_FALSE(c[3][3] == c33);
3851 :
3852 : // first level verification, which is exactly what the
3853 : // division function does, it does not mean that the
3854 : // division works properly, though
3855 : //
3856 2 : snapdev::matrix<double> r(4, 4);
3857 1 : r = a * t;
3858 :
3859 1 : CATCH_REQUIRE(c[0][0] == r[0][0]);
3860 1 : CATCH_REQUIRE(c[0][1] == r[0][1]);
3861 1 : CATCH_REQUIRE(c[0][2] == r[0][2]);
3862 1 : CATCH_REQUIRE(c[0][3] == r[0][3]);
3863 1 : CATCH_REQUIRE(c[1][0] == r[1][0]);
3864 1 : CATCH_REQUIRE(c[1][1] == r[1][1]);
3865 1 : CATCH_REQUIRE(c[1][2] == r[1][2]);
3866 1 : CATCH_REQUIRE(c[1][3] == r[1][3]);
3867 1 : CATCH_REQUIRE(c[2][0] == r[2][0]);
3868 1 : CATCH_REQUIRE(c[2][1] == r[2][1]);
3869 1 : CATCH_REQUIRE(c[2][2] == r[2][2]);
3870 1 : CATCH_REQUIRE(c[2][3] == r[2][3]);
3871 1 : CATCH_REQUIRE(c[3][0] == r[3][0]);
3872 1 : CATCH_REQUIRE(c[3][1] == r[3][1]);
3873 1 : CATCH_REQUIRE(c[3][2] == r[3][2]);
3874 1 : CATCH_REQUIRE(c[3][3] == r[3][3]);
3875 :
3876 1 : double const determinant(b.determinant());
3877 2 : snapdev::matrix<double> adjugate(b.adjugate());
3878 :
3879 2 : snapdev::matrix<double> inv(adjugate * (1.0 / determinant));
3880 2 : snapdev::matrix<double> div(a * inv);
3881 :
3882 1 : CATCH_REQUIRE(fabs(c[0][0] - div[0][0]) < 0.0001);
3883 1 : CATCH_REQUIRE(fabs(c[0][1] - div[0][1]) < 0.0001);
3884 1 : CATCH_REQUIRE(fabs(c[0][2] - div[0][2]) < 0.0001);
3885 1 : CATCH_REQUIRE(fabs(c[0][3] - div[0][3]) < 0.0001);
3886 1 : CATCH_REQUIRE(fabs(c[1][0] - div[1][0]) < 0.0001);
3887 1 : CATCH_REQUIRE(fabs(c[1][1] - div[1][1]) < 0.0001);
3888 1 : CATCH_REQUIRE(fabs(c[1][2] - div[1][2]) < 0.0001);
3889 1 : CATCH_REQUIRE(fabs(c[1][3] - div[1][3]) < 0.0001);
3890 1 : CATCH_REQUIRE(fabs(c[2][0] - div[2][0]) < 0.0001);
3891 1 : CATCH_REQUIRE(fabs(c[2][1] - div[2][1]) < 0.0001);
3892 1 : CATCH_REQUIRE(fabs(c[2][2] - div[2][2]) < 0.0001);
3893 1 : CATCH_REQUIRE(fabs(c[2][3] - div[2][3]) < 0.0001);
3894 1 : CATCH_REQUIRE(fabs(c[3][0] - div[3][0]) < 0.0001);
3895 1 : CATCH_REQUIRE(fabs(c[3][1] - div[3][1]) < 0.0001);
3896 1 : CATCH_REQUIRE(fabs(c[3][2] - div[3][2]) < 0.0001);
3897 1 : CATCH_REQUIRE(fabs(c[3][3] - div[3][3]) < 0.0001);
3898 : }
3899 : CATCH_END_SECTION()
3900 :
3901 8 : CATCH_START_SECTION("matrix: a/=b")
3902 : {
3903 : // setup A
3904 : //
3905 2 : snapdev::matrix<double> a(4, 4);
3906 :
3907 1 : CATCH_REQUIRE(a.rows() == 4);
3908 1 : CATCH_REQUIRE(a.columns() == 4);
3909 :
3910 1 : double const a00(frand());
3911 1 : double const a01(frand());
3912 1 : double const a02(frand());
3913 1 : double const a03(frand());
3914 1 : double const a10(frand());
3915 1 : double const a11(frand());
3916 1 : double const a12(frand());
3917 1 : double const a13(frand());
3918 1 : double const a20(frand());
3919 1 : double const a21(frand());
3920 1 : double const a22(frand());
3921 1 : double const a23(frand());
3922 1 : double const a30(frand());
3923 1 : double const a31(frand());
3924 1 : double const a32(frand());
3925 1 : double const a33(frand());
3926 :
3927 1 : a[0][0] = a00;
3928 1 : a[0][1] = a01;
3929 1 : a[0][2] = a02;
3930 1 : a[0][3] = a03;
3931 1 : a[1][0] = a10;
3932 1 : a[1][1] = a11;
3933 1 : a[1][2] = a12;
3934 1 : a[1][3] = a13;
3935 1 : a[2][0] = a20;
3936 1 : a[2][1] = a21;
3937 1 : a[2][2] = a22;
3938 1 : a[2][3] = a23;
3939 1 : a[3][0] = a30;
3940 1 : a[3][1] = a31;
3941 1 : a[3][2] = a32;
3942 1 : a[3][3] = a33;
3943 :
3944 1 : CATCH_REQUIRE(a[0][0] == a00);
3945 1 : CATCH_REQUIRE(a[0][1] == a01);
3946 1 : CATCH_REQUIRE(a[0][2] == a02);
3947 1 : CATCH_REQUIRE(a[0][3] == a03);
3948 1 : CATCH_REQUIRE(a[1][0] == a10);
3949 1 : CATCH_REQUIRE(a[1][1] == a11);
3950 1 : CATCH_REQUIRE(a[1][2] == a12);
3951 1 : CATCH_REQUIRE(a[1][3] == a13);
3952 1 : CATCH_REQUIRE(a[2][0] == a20);
3953 1 : CATCH_REQUIRE(a[2][1] == a21);
3954 1 : CATCH_REQUIRE(a[2][2] == a22);
3955 1 : CATCH_REQUIRE(a[2][3] == a23);
3956 1 : CATCH_REQUIRE(a[3][0] == a30);
3957 1 : CATCH_REQUIRE(a[3][1] == a31);
3958 1 : CATCH_REQUIRE(a[3][2] == a32);
3959 1 : CATCH_REQUIRE(a[3][3] == a33);
3960 :
3961 : // setup B
3962 : //
3963 2 : snapdev::matrix<double> b(4, 4);
3964 :
3965 1 : CATCH_REQUIRE(b.rows() == 4);
3966 1 : CATCH_REQUIRE(b.columns() == 4);
3967 :
3968 1 : double const b00(frand());
3969 1 : double const b01(frand());
3970 1 : double const b02(frand());
3971 1 : double const b03(frand());
3972 1 : double const b10(frand());
3973 1 : double const b11(frand());
3974 1 : double const b12(frand());
3975 1 : double const b13(frand());
3976 1 : double const b20(frand());
3977 1 : double const b21(frand());
3978 1 : double const b22(frand());
3979 1 : double const b23(frand());
3980 1 : double const b30(frand());
3981 1 : double const b31(frand());
3982 1 : double const b32(frand());
3983 1 : double const b33(frand());
3984 :
3985 1 : b[0][0] = b00;
3986 1 : b[0][1] = b01;
3987 1 : b[0][2] = b02;
3988 1 : b[0][3] = b03;
3989 1 : b[1][0] = b10;
3990 1 : b[1][1] = b11;
3991 1 : b[1][2] = b12;
3992 1 : b[1][3] = b13;
3993 1 : b[2][0] = b20;
3994 1 : b[2][1] = b21;
3995 1 : b[2][2] = b22;
3996 1 : b[2][3] = b23;
3997 1 : b[3][0] = b30;
3998 1 : b[3][1] = b31;
3999 1 : b[3][2] = b32;
4000 1 : b[3][3] = b33;
4001 :
4002 1 : CATCH_REQUIRE(b[0][0] == b00);
4003 1 : CATCH_REQUIRE(b[0][1] == b01);
4004 1 : CATCH_REQUIRE(b[0][2] == b02);
4005 1 : CATCH_REQUIRE(b[0][3] == b03);
4006 1 : CATCH_REQUIRE(b[1][0] == b10);
4007 1 : CATCH_REQUIRE(b[1][1] == b11);
4008 1 : CATCH_REQUIRE(b[1][2] == b12);
4009 1 : CATCH_REQUIRE(b[1][3] == b13);
4010 1 : CATCH_REQUIRE(b[2][0] == b20);
4011 1 : CATCH_REQUIRE(b[2][1] == b21);
4012 1 : CATCH_REQUIRE(b[2][2] == b22);
4013 1 : CATCH_REQUIRE(b[2][3] == b23);
4014 1 : CATCH_REQUIRE(b[3][0] == b30);
4015 1 : CATCH_REQUIRE(b[3][1] == b31);
4016 1 : CATCH_REQUIRE(b[3][2] == b32);
4017 1 : CATCH_REQUIRE(b[3][3] == b33);
4018 :
4019 1 : double const determinant(b.determinant());
4020 2 : snapdev::matrix<double> adjugate(b.adjugate());
4021 :
4022 2 : snapdev::matrix<double> inv(adjugate * (1.0 / determinant));
4023 2 : snapdev::matrix<double> div(a * inv);
4024 :
4025 : // run operation A /= B
4026 : //
4027 1 : a /= b;
4028 :
4029 : // this could fail if one of bXX is 0.0
4030 : //
4031 : //CATCH_REQUIRE(a[0][0] == a00);
4032 : //CATCH_REQUIRE(a[0][1] == a01);
4033 : //CATCH_REQUIRE(a[0][2] == a02);
4034 : //CATCH_REQUIRE(a[0][3] == a03);
4035 : //CATCH_REQUIRE(a[1][0] == a10);
4036 : //CATCH_REQUIRE(a[1][1] == a11);
4037 : //CATCH_REQUIRE(a[1][2] == a12);
4038 : //CATCH_REQUIRE(a[1][3] == a13);
4039 : //CATCH_REQUIRE(a[2][0] == a20);
4040 : //CATCH_REQUIRE(a[2][1] == a21);
4041 : //CATCH_REQUIRE(a[2][2] == a22);
4042 : //CATCH_REQUIRE(a[2][3] == a23);
4043 : //CATCH_REQUIRE(a[3][0] == a30);
4044 : //CATCH_REQUIRE(a[3][1] == a31);
4045 : //CATCH_REQUIRE(a[3][2] == a32);
4046 : //CATCH_REQUIRE(a[3][3] == a33);
4047 :
4048 1 : CATCH_REQUIRE(b[0][0] == b00);
4049 1 : CATCH_REQUIRE(b[0][1] == b01);
4050 1 : CATCH_REQUIRE(b[0][2] == b02);
4051 1 : CATCH_REQUIRE(b[0][3] == b03);
4052 1 : CATCH_REQUIRE(b[1][0] == b10);
4053 1 : CATCH_REQUIRE(b[1][1] == b11);
4054 1 : CATCH_REQUIRE(b[1][2] == b12);
4055 1 : CATCH_REQUIRE(b[1][3] == b13);
4056 1 : CATCH_REQUIRE(b[2][0] == b20);
4057 1 : CATCH_REQUIRE(b[2][1] == b21);
4058 1 : CATCH_REQUIRE(b[2][2] == b22);
4059 1 : CATCH_REQUIRE(b[2][3] == b23);
4060 1 : CATCH_REQUIRE(b[3][0] == b30);
4061 1 : CATCH_REQUIRE(b[3][1] == b31);
4062 1 : CATCH_REQUIRE(b[3][2] == b32);
4063 1 : CATCH_REQUIRE(b[3][3] == b33);
4064 :
4065 1 : CATCH_REQUIRE(fabs(a[0][0] - div[0][0]) < 0.0001);
4066 1 : CATCH_REQUIRE(fabs(a[0][1] - div[0][1]) < 0.0001);
4067 1 : CATCH_REQUIRE(fabs(a[0][2] - div[0][2]) < 0.0001);
4068 1 : CATCH_REQUIRE(fabs(a[0][3] - div[0][3]) < 0.0001);
4069 1 : CATCH_REQUIRE(fabs(a[1][0] - div[1][0]) < 0.0001);
4070 1 : CATCH_REQUIRE(fabs(a[1][1] - div[1][1]) < 0.0001);
4071 1 : CATCH_REQUIRE(fabs(a[1][2] - div[1][2]) < 0.0001);
4072 1 : CATCH_REQUIRE(fabs(a[1][3] - div[1][3]) < 0.0001);
4073 1 : CATCH_REQUIRE(fabs(a[2][0] - div[2][0]) < 0.0001);
4074 1 : CATCH_REQUIRE(fabs(a[2][1] - div[2][1]) < 0.0001);
4075 1 : CATCH_REQUIRE(fabs(a[2][2] - div[2][2]) < 0.0001);
4076 1 : CATCH_REQUIRE(fabs(a[2][3] - div[2][3]) < 0.0001);
4077 1 : CATCH_REQUIRE(fabs(a[3][0] - div[3][0]) < 0.0001);
4078 1 : CATCH_REQUIRE(fabs(a[3][1] - div[3][1]) < 0.0001);
4079 1 : CATCH_REQUIRE(fabs(a[3][2] - div[3][2]) < 0.0001);
4080 1 : CATCH_REQUIRE(fabs(a[3][3] - div[3][3]) < 0.0001);
4081 : }
4082 : CATCH_END_SECTION()
4083 : }
4084 8 : }
4085 :
4086 :
4087 5 : CATCH_TEST_CASE("matrix_color", "[matrix]")
4088 : {
4089 : // generate a few brightness matrices
4090 : //
4091 6 : CATCH_GIVEN("brightness")
4092 : {
4093 2 : CATCH_START_SECTION("matrix: b=a*<brightness> (a is identity)")
4094 : {
4095 : // setup A
4096 : //
4097 2 : snapdev::matrix<double> a(4, 4);
4098 :
4099 1 : CATCH_REQUIRE(a.rows() == 4);
4100 1 : CATCH_REQUIRE(a.columns() == 4);
4101 :
4102 : // setup B
4103 : //
4104 2 : snapdev::matrix<double> b(4, 4);
4105 :
4106 1 : CATCH_REQUIRE(b.rows() == 4);
4107 1 : CATCH_REQUIRE(b.columns() == 4);
4108 :
4109 1001 : for(int idx(0); idx < 1000; ++idx)
4110 : {
4111 1000 : double const brightness(static_cast<double>(idx) / 1000.0);
4112 1000 : b = a.brightness(brightness);
4113 :
4114 1000 : CATCH_REQUIRE(fabs(b[0][0] - brightness) < 0.0001);
4115 1000 : CATCH_REQUIRE(fabs(b[0][1] - 0.0) < 0.0001);
4116 1000 : CATCH_REQUIRE(fabs(b[0][2] - 0.0) < 0.0001);
4117 1000 : CATCH_REQUIRE(fabs(b[0][3] - 0.0) < 0.0001);
4118 1000 : CATCH_REQUIRE(fabs(b[1][0] - 0.0) < 0.0001);
4119 1000 : CATCH_REQUIRE(fabs(b[1][1] - brightness) < 0.0001);
4120 1000 : CATCH_REQUIRE(fabs(b[1][2] - 0.0) < 0.0001);
4121 1000 : CATCH_REQUIRE(fabs(b[1][3] - 0.0) < 0.0001);
4122 1000 : CATCH_REQUIRE(fabs(b[2][0] - 0.0) < 0.0001);
4123 1000 : CATCH_REQUIRE(fabs(b[2][1] - 0.0) < 0.0001);
4124 1000 : CATCH_REQUIRE(fabs(b[2][2] - brightness) < 0.0001);
4125 1000 : CATCH_REQUIRE(fabs(b[2][3] - 0.0) < 0.0001);
4126 1000 : CATCH_REQUIRE(fabs(b[3][0] - 0.0) < 0.0001);
4127 1000 : CATCH_REQUIRE(fabs(b[3][1] - 0.0) < 0.0001);
4128 1000 : CATCH_REQUIRE(fabs(b[3][2] - 0.0) < 0.0001);
4129 1000 : CATCH_REQUIRE(fabs(b[3][3] - 1.0) < 0.0001);
4130 :
4131 : // * brightness_r & 0 & 0 & 0
4132 : // * \\ 0 & brightness_g & 0 & 0
4133 : // * \\ 0 & 0 & brightness_b & 0
4134 : // * \\ 0 & 0 & 0 & 1
4135 :
4136 : }
4137 : }
4138 : CATCH_END_SECTION()
4139 : }
4140 :
4141 : // generate a few saturation matrices
4142 : //
4143 6 : CATCH_GIVEN("saturation")
4144 : {
4145 2 : CATCH_START_SECTION("matrix: b=a*<saturation> (a is identity)")
4146 : {
4147 : // setup A
4148 : //
4149 2 : snapdev::matrix<double> a(4, 4);
4150 :
4151 1 : CATCH_REQUIRE(a.rows() == 4);
4152 1 : CATCH_REQUIRE(a.columns() == 4);
4153 :
4154 : // setup B
4155 : //
4156 2 : snapdev::matrix<double> b(4, 4);
4157 :
4158 1 : CATCH_REQUIRE(b.rows() == 4);
4159 1 : CATCH_REQUIRE(b.columns() == 4);
4160 :
4161 1 : double luma[6][3];
4162 :
4163 1 : luma[0][0] = a.HDTV_LUMA_RED ;
4164 1 : luma[0][1] = a.HDTV_LUMA_GREEN;
4165 1 : luma[0][2] = a.HDTV_LUMA_BLUE ;
4166 :
4167 1 : luma[1][0] = a.LED_LUMA_RED ;
4168 1 : luma[1][1] = a.LED_LUMA_GREEN;
4169 1 : luma[1][2] = a.LED_LUMA_BLUE ;
4170 :
4171 1 : luma[2][0] = a.CRT_LUMA_RED ;
4172 1 : luma[2][1] = a.CRT_LUMA_GREEN;
4173 1 : luma[2][2] = a.CRT_LUMA_BLUE ;
4174 :
4175 1 : luma[3][0] = a.NTSC_LUMA_RED ;
4176 1 : luma[3][1] = a.NTSC_LUMA_GREEN;
4177 1 : luma[3][2] = a.NTSC_LUMA_BLUE ;
4178 :
4179 1 : luma[4][0] = a.AVERAGE_LUMA_RED ;
4180 1 : luma[4][1] = a.AVERAGE_LUMA_GREEN;
4181 1 : luma[4][2] = a.AVERAGE_LUMA_BLUE ;
4182 :
4183 1 : luma[5][0] = 0.2; // to test any dynamic version if such exists for this algorithm
4184 1 : luma[5][1] = 0.7;
4185 1 : luma[5][2] = 0.1;
4186 :
4187 7 : for(int l(0); l < 6; ++l)
4188 : {
4189 : // set the luma weights
4190 : //
4191 6 : a.set_luma_vector(luma[l][0], luma[l][1], luma[l][2]);
4192 :
4193 6006 : for(int idx(0); idx < 1000; ++idx)
4194 : {
4195 6000 : double const s(static_cast<double>(idx) / 1000.0);
4196 6000 : b = a.saturation(s);
4197 :
4198 6000 : CATCH_REQUIRE(fabs(b[0][0] - (luma[l][0] * (1.0 - s) + s)) < 0.0001);
4199 6000 : CATCH_REQUIRE(fabs(b[0][1] - (luma[l][0] * (1.0 - s))) < 0.0001);
4200 6000 : CATCH_REQUIRE(fabs(b[0][2] - (luma[l][0] * (1.0 - s))) < 0.0001);
4201 6000 : CATCH_REQUIRE(fabs(b[0][3] - (0.0)) < 0.0001);
4202 6000 : CATCH_REQUIRE(fabs(b[1][0] - (luma[l][1] * (1.0 - s))) < 0.0001);
4203 6000 : CATCH_REQUIRE(fabs(b[1][1] - (luma[l][1] * (1.0 - s) + s)) < 0.0001);
4204 6000 : CATCH_REQUIRE(fabs(b[1][2] - (luma[l][1] * (1.0 - s))) < 0.0001);
4205 6000 : CATCH_REQUIRE(fabs(b[1][3] - (0.0)) < 0.0001);
4206 6000 : CATCH_REQUIRE(fabs(b[2][0] - (luma[l][2] * (1.0 - s))) < 0.0001);
4207 6000 : CATCH_REQUIRE(fabs(b[2][1] - (luma[l][2] * (1.0 - s))) < 0.0001);
4208 6000 : CATCH_REQUIRE(fabs(b[2][2] - (luma[l][2] * (1.0 - s) + s)) < 0.0001);
4209 6000 : CATCH_REQUIRE(fabs(b[2][3] - (0.0)) < 0.0001);
4210 6000 : CATCH_REQUIRE(fabs(b[3][0] - (0.0)) < 0.0001);
4211 6000 : CATCH_REQUIRE(fabs(b[3][1] - (0.0)) < 0.0001);
4212 6000 : CATCH_REQUIRE(fabs(b[3][2] - (0.0)) < 0.0001);
4213 6000 : CATCH_REQUIRE(fabs(b[3][3] - (1.0)) < 0.0001);
4214 : }
4215 : }
4216 : }
4217 : CATCH_END_SECTION()
4218 : }
4219 :
4220 : // generate a few brightness matrices
4221 : //
4222 6 : CATCH_GIVEN("hue")
4223 : {
4224 2 : CATCH_START_SECTION("matrix: b=a*<hue> (a is identity)")
4225 : {
4226 : // setup A
4227 : //
4228 2 : snapdev::matrix<double> a(4, 4);
4229 :
4230 1 : CATCH_REQUIRE(a.rows() == 4);
4231 1 : CATCH_REQUIRE(a.columns() == 4);
4232 :
4233 : // setup B
4234 : //
4235 2 : snapdev::matrix<double> b(4, 4);
4236 :
4237 1 : CATCH_REQUIRE(b.rows() == 4);
4238 1 : CATCH_REQUIRE(b.columns() == 4);
4239 :
4240 1 : double luma[6][3];
4241 :
4242 1 : luma[0][0] = a.HDTV_LUMA_RED ;
4243 1 : luma[0][1] = a.HDTV_LUMA_GREEN;
4244 1 : luma[0][2] = a.HDTV_LUMA_BLUE ;
4245 :
4246 1 : luma[1][0] = a.LED_LUMA_RED ;
4247 1 : luma[1][1] = a.LED_LUMA_GREEN;
4248 1 : luma[1][2] = a.LED_LUMA_BLUE ;
4249 :
4250 1 : luma[2][0] = a.CRT_LUMA_RED ;
4251 1 : luma[2][1] = a.CRT_LUMA_GREEN;
4252 1 : luma[2][2] = a.CRT_LUMA_BLUE ;
4253 :
4254 1 : luma[3][0] = a.NTSC_LUMA_RED ;
4255 1 : luma[3][1] = a.NTSC_LUMA_GREEN;
4256 1 : luma[3][2] = a.NTSC_LUMA_BLUE ;
4257 :
4258 1 : luma[4][0] = a.AVERAGE_LUMA_RED ;
4259 1 : luma[4][1] = a.AVERAGE_LUMA_GREEN;
4260 1 : luma[4][2] = a.AVERAGE_LUMA_BLUE ;
4261 :
4262 1 : luma[5][0] = 0.2; // to test any dynamic version if such exists for this algorithm
4263 1 : luma[5][1] = 0.7;
4264 1 : luma[5][2] = 0.1;
4265 :
4266 1 : std::vector<std::string> last_hue_matrix{
4267 : "hdtv",
4268 : "led",
4269 : "crt",
4270 : "ntsc",
4271 : "average",
4272 : "dynamic"
4273 2 : };
4274 :
4275 7 : for(int l(0); l < 6; ++l)
4276 : {
4277 : // set the luma weights
4278 : //
4279 6 : a.set_luma_vector(luma[l][0], luma[l][1], luma[l][2]);
4280 :
4281 : //std::cerr << "*** working on luma " << l << "\n";
4282 6006 : for(int idx(0); idx < 1000; ++idx)
4283 : {
4284 : // hue is an angle from 0 to 360 (albeit in radian)
4285 : //
4286 6000 : double const hue(static_cast<double>(idx) / 1000.0 * M_PI * 2.0);
4287 6000 : b = a.hue(hue);
4288 :
4289 : #ifdef _DEBUG
4290 6000 : CATCH_REQUIRE(a.get_last_hue_matrix() == last_hue_matrix[l]);
4291 : #endif
4292 :
4293 : // this one requires us to work!
4294 : // here we calculate the hue matrix by hand instead of using
4295 : // pre-calculated numbers; that way we know we have it right
4296 : //
4297 :
4298 : // $R_r$ (red rotation)
4299 : //
4300 12000 : snapdev::matrix<double> r_r(4, 4);
4301 6000 : r_r[1][1] = 1.0 / sqrt(2.0);
4302 6000 : r_r[1][2] = 1.0 / sqrt(2.0);
4303 6000 : r_r[2][1] = -1.0 / sqrt(2.0);
4304 6000 : r_r[2][2] = 1.0 / sqrt(2.0);
4305 :
4306 : // $R_g$ (green rotation)
4307 : //
4308 12000 : snapdev::matrix<double> r_g(4, 4);
4309 6000 : r_g[0][0] = sqrt(2.0) / sqrt(3.0);
4310 6000 : r_g[0][2] = 1.0 / sqrt(3.0);
4311 6000 : r_g[2][0] = -1.0 / sqrt(3.0);
4312 6000 : r_g[2][2] = sqrt(2.0) / sqrt(3.0);
4313 :
4314 : // $R_rg = R_r R_g$ (red and green rotations combined)
4315 : //
4316 12000 : snapdev::matrix<double> r_rg(r_r * r_g);
4317 :
4318 : // $w$ (luma vector, a.k.a. color weights to calcalute the perfect grayscale for your monitor)
4319 : //
4320 12000 : snapdev::matrix<double> w(a.get_luma_vector());
4321 : //w[0][0] = a.RED_WEIGHT;
4322 : //w[1][0] = a.GREEN_WEIGHT;
4323 : //w[2][0] = a.BLUE_WEIGHT;
4324 : //w[3][0] = 0;
4325 :
4326 : // sm (skew matrix)
4327 : //
4328 12000 : snapdev::matrix<double> sm(r_rg * w);
4329 :
4330 : // s (scaling with skew matrix)
4331 : //
4332 12000 : snapdev::matrix<double> s(4, 4);
4333 6000 : s[0][2] = sm[0][0] / sm[2][0];
4334 6000 : s[1][2] = sm[1][0] / sm[2][0];
4335 :
4336 : // p (skewed red & green rotation)
4337 : //
4338 12000 : snapdev::matrix<double> p(r_rg);
4339 6000 : p *= s;
4340 :
4341 : // r_b (blue rotation)
4342 : //
4343 12000 : snapdev::matrix<double> r_b(4, 4);
4344 6000 : double const rot_cos(cos(hue));
4345 6000 : double const rot_sin(sin(hue));
4346 6000 : r_b[0][0] = rot_cos;
4347 6000 : r_b[0][1] = rot_sin;
4348 6000 : r_b[1][0] = -rot_sin;
4349 6000 : r_b[1][1] = rot_cos;
4350 :
4351 : // hue_matrix (the resulting matrix depending on angle)
4352 : //
4353 12000 : snapdev::matrix<double> hue_matrix(p * r_b / p);
4354 :
4355 : //std::cout << "Got b = " << b << "\n";
4356 : //std::cout << "Got our hue matrix = " << hue_matrix << "\n";
4357 :
4358 6000 : CATCH_REQUIRE(fabs(b[0][0] - hue_matrix[0][0]) < 0.0001);
4359 6000 : CATCH_REQUIRE(fabs(b[0][1] - hue_matrix[0][1]) < 0.0001);
4360 6000 : CATCH_REQUIRE(fabs(b[0][2] - hue_matrix[0][2]) < 0.0001);
4361 6000 : CATCH_REQUIRE(fabs(b[0][3] - hue_matrix[0][3]) < 0.0001);
4362 6000 : CATCH_REQUIRE(fabs(b[1][0] - hue_matrix[1][0]) < 0.0001);
4363 6000 : CATCH_REQUIRE(fabs(b[1][1] - hue_matrix[1][1]) < 0.0001);
4364 6000 : CATCH_REQUIRE(fabs(b[1][2] - hue_matrix[1][2]) < 0.0001);
4365 6000 : CATCH_REQUIRE(fabs(b[1][3] - hue_matrix[1][3]) < 0.0001);
4366 6000 : CATCH_REQUIRE(fabs(b[2][0] - hue_matrix[2][0]) < 0.0001);
4367 6000 : CATCH_REQUIRE(fabs(b[2][1] - hue_matrix[2][1]) < 0.0001);
4368 6000 : CATCH_REQUIRE(fabs(b[2][2] - hue_matrix[2][2]) < 0.0001);
4369 6000 : CATCH_REQUIRE(fabs(b[2][3] - hue_matrix[2][3]) < 0.0001);
4370 6000 : CATCH_REQUIRE(fabs(b[3][0] - hue_matrix[3][0]) < 0.0001);
4371 6000 : CATCH_REQUIRE(fabs(b[3][1] - hue_matrix[3][1]) < 0.0001);
4372 6000 : CATCH_REQUIRE(fabs(b[3][2] - hue_matrix[3][2]) < 0.0001);
4373 6000 : CATCH_REQUIRE(fabs(b[3][3] - hue_matrix[3][3]) < 0.0001);
4374 :
4375 : // * $$
4376 : // * R_r =
4377 : // * \begin{bmatrix}
4378 : // * 1 & 0 & 0 & 0
4379 : // * \\ 0 & {1 \over \sqrt 2 } & {1 \over \sqrt 2 } & 0
4380 : // * \\ 0 & -{1 \over \sqrt 2 } & {1 \over \sqrt 2 } & 0
4381 : // * \\ 0 & 0 & 0 & 1
4382 : // * \end{bmatrix}
4383 : // * $$
4384 : // *
4385 : // * $$
4386 : // * R_g =
4387 : // * \begin{bmatrix}
4388 : // * 1 & 0 & 0 & 0
4389 : // * \\ 0 & {\sqrt 2 \over \sqrt 3} & -{1 \over \sqrt 3} & 0
4390 : // * \\ 0 & {1 \over \sqrt 3} & {\sqrt 2 \over \sqrt 3} & 0
4391 : // * \\ 0 & 0 & 0 & 1
4392 : // * \end{bmatrix}
4393 : // * $$
4394 : // *
4395 : // * $$
4396 : // * R_{rg} = R_r R_g
4397 : // * $$
4398 : // *
4399 : // * $$
4400 : // * \begin{bmatrix}
4401 : // * L_r
4402 : // * \\ L_g
4403 : // * \\ L_b
4404 : // * \\ 0
4405 : // * \end{bmatrix}
4406 : // * =
4407 : // * R_{rg}
4408 : // * \begin{bmatrix}
4409 : // * W_r
4410 : // * \\ W_g
4411 : // * \\ W_b
4412 : // * \\ 0
4413 : // * \end{bmatrix}
4414 : // * $$
4415 : // *
4416 : // * Now we define a rotation matrix for the blue axis. This one uses
4417 : // * a variable as the angle. This is the actual rotation angle which
4418 : // * can go from 0 to $2 \pi$:
4419 : // *
4420 : // * $$
4421 : // * R_b =
4422 : // * \begin{bmatrix}
4423 : // * cos( \theta ) & sin( \theta ) & 0 & 0
4424 : // * \\ -sin( \theta ) & cos( \theta ) & 0 & 0
4425 : // * \\ 0 & 0 & 1 & 0
4426 : // * \\ 0 & 0 & 0 & 1
4427 : // * \end{bmatrix}
4428 : // * $$
4429 : // *
4430 : // * Now we have all the matrices to caculate the hue rotation
4431 : // * of all the components of an image:
4432 : // *
4433 : // * $$
4434 : // * H = R_{rg} S R_b S^{-1} R_{rg}^{-1}
4435 : // * $$
4436 : // *
4437 : // * $H$ can then be used as in:
4438 : // *
4439 : // * $$
4440 : // * \begin{bmatrix}
4441 : // * R'
4442 : // * \\ G'
4443 : // * \\ B'
4444 : // * \end{bmatrix}
4445 : // * =
4446 : // * H
4447 : // * \begin{bmatrix}
4448 : // * R
4449 : // * \\ G
4450 : // * \\ B
4451 : // * \end{bmatrix}
4452 : // * $$
4453 : // *
4454 : }
4455 : }
4456 : }
4457 : CATCH_END_SECTION()
4458 : }
4459 9 : }
4460 :
4461 :
4462 : // vim: ts=4 sw=4 et
|