Line data Source code
1 : /*
2 : * File: TestMatrixExtras.cpp
3 : * Author: Pantelis Sopasakis
4 : *
5 : * Created on Nov 8, 2015, 4:33:46 PM
6 : *
7 : * ForBES is free software: you can redistribute it and/or modify
8 : * it under the terms of the GNU Lesser General Public License as published by
9 : * the Free Software Foundation, either version 3 of the License, or
10 : * (at your option) any later version.
11 : *
12 : * ForBES is distributed in the hope that it will be useful,
13 : * but WITHOUT ANY WARRANTY; without even the implied warranty of
14 : * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 : * GNU Lesser General Public License for more details.
16 : *
17 : * You should have received a copy of the GNU Lesser General Public License
18 : * along with ForBES. If not, see <http://www.gnu.org/licenses/>.
19 : */
20 :
21 : #include "TestMatrixExtras.h"
22 : #include "MatrixFactory.h"
23 : #include "ForBES.h"
24 : #include <cmath>
25 :
26 1 : CPPUNIT_TEST_SUITE_REGISTRATION(TestMatrixExtras);
27 :
28 23 : TestMatrixExtras::TestMatrixExtras() {
29 23 : }
30 :
31 46 : TestMatrixExtras::~TestMatrixExtras() {
32 46 : }
33 :
34 23 : void TestMatrixExtras::setUp() {
35 23 : }
36 :
37 23 : void TestMatrixExtras::tearDown() {
38 23 : Matrix::destroy_handle();
39 23 : }
40 :
41 1 : void TestMatrixExtras::test_add_wrong_args() {
42 1 : Matrix A(5, 6);
43 2 : Matrix B(5, 7);
44 2 : _ASSERT_EXCEPTION(Matrix::add(A, 1.3, B, 0.8), std::invalid_argument);
45 1 : }
46 :
47 1 : void TestMatrixExtras::test_mult_wrong_args() {
48 : {
49 1 : Matrix A(5, 6);
50 2 : Matrix B(6, 7);
51 2 : Matrix C(5, 4);
52 2 : _ASSERT_EXCEPTION(Matrix::mult(C, 1.3, A, B, 0.8), std::invalid_argument);
53 : }
54 : {
55 1 : Matrix A(5, 6);
56 2 : Matrix B(3, 7);
57 2 : Matrix C(5, 7);
58 2 : _ASSERT_EXCEPTION(Matrix::mult(C, 1.3, A, B, 0.8), std::invalid_argument);
59 : }
60 1 : }
61 :
62 : /***** ADDITION ******/
63 :
64 1 : void TestMatrixExtras::test_add_DD() {
65 1 : size_t n = 10;
66 1 : size_t m = 15;
67 1 : size_t repetitions = 20;
68 :
69 21 : for (size_t r = 0; r < repetitions; r++) {
70 20 : Matrix A = MatrixFactory::MakeRandomMatrix(n, m, 2.0, 1.0);
71 40 : Matrix B = MatrixFactory::MakeRandomMatrix(n, m, 2.0, 1.0);
72 :
73 20 : double alpha = 2.0 * static_cast<double> (std::rand()) / static_cast<double> (RAND_MAX);
74 20 : double gamma = -0.5 + static_cast<double> (std::rand()) / static_cast<double> (RAND_MAX);
75 :
76 40 : Matrix A_copy(A);
77 40 : Matrix B_copy(B);
78 :
79 : // A = gamma*A + alpha*B
80 20 : int status = Matrix::add(A, alpha, B, gamma);
81 20 : _ASSERT_EQ(ForBESUtils::STATUS_OK, status);
82 20 : _ASSERT_EQ(B_copy, B);
83 :
84 20 : A_copy *= gamma;
85 20 : B_copy *= alpha;
86 40 : Matrix C = A_copy + B_copy;
87 :
88 20 : _ASSERT_EQ(C, A);
89 20 : }
90 :
91 1 : }
92 :
93 1 : void TestMatrixExtras::test_add_DS() {
94 1 : size_t n = 10;
95 1 : size_t m = 15;
96 1 : size_t nnz = 100;
97 1 : Matrix A = MatrixFactory::MakeRandomMatrix(n, m, 2.0, 1.0);
98 2 : Matrix B = MatrixFactory::MakeRandomSparse(n, m, nnz, 2.0, 1.0);
99 :
100 1 : double alpha = -1.0 + 2.0 * static_cast<double> (std::rand()) / static_cast<double> (RAND_MAX);
101 1 : double gamma = -0.5 + static_cast<double> (std::rand()) / static_cast<double> (RAND_MAX);
102 :
103 :
104 2 : Matrix R(A);
105 1 : R *= gamma;
106 2 : Matrix T(B);
107 1 : T *= alpha;
108 2 : Matrix S = R + T;
109 :
110 1 : int status = Matrix::add(A, alpha, B, gamma); // A = gamma * A + alpha * B
111 1 : _ASSERT_EQ(ForBESUtils::STATUS_OK, status);
112 2 : _ASSERT_EQ(S, A);
113 :
114 1 : }
115 :
116 1 : void TestMatrixExtras::test_add_DX() {
117 1 : size_t n = 10;
118 :
119 1 : Matrix D = MatrixFactory::MakeRandomMatrix(n, n, 2.0, 6.0);
120 2 : Matrix Dorig(D);
121 :
122 2 : Matrix X = MatrixFactory::MakeRandomMatrix(n, 1, 4.0, 10.0);
123 1 : X.toggle_diagonal();
124 :
125 1 : double alpha = -4.23512323;
126 1 : double gamma = 0.24325345;
127 1 : Matrix::add(D, alpha, X, gamma);
128 :
129 1 : const double tol = 1e-7;
130 :
131 : // check diagonal entries
132 11 : for (size_t i = 0; i < n; i++) {
133 10 : _ASSERT_NUM_EQ(gamma * Dorig.get(i, i) + alpha * X[i], D.get(i, i), tol);
134 : }
135 :
136 : // check all non-diagonal entries
137 11 : for (size_t i = 0; i < n; i++) {
138 55 : for (size_t j = 0; j < i; j++) {
139 45 : _ASSERT_NUM_EQ(gamma * Dorig.get(i, j), D.get(i, j), tol);
140 : }
141 55 : for (size_t j = i + 1; j < n; j++) {
142 45 : _ASSERT_NUM_EQ(gamma * Dorig.get(i, j), D.get(i, j), tol);
143 : }
144 1 : }
145 :
146 1 : }
147 :
148 1 : void TestMatrixExtras::test_add_DST() {
149 7 : for (size_t n = 11; n < 40; n += 5) {
150 42 : for (size_t m = 13; m < 30; m += 3) {
151 36 : size_t nnz = 100;
152 36 : Matrix A = MatrixFactory::MakeRandomMatrix(n, m, 2.0, 1.0);
153 72 : Matrix B = MatrixFactory::MakeRandomSparse(m, n, nnz, 2.0, 1.0);
154 :
155 72 : Matrix A_copy(A);
156 72 : Matrix B_copy(B);
157 :
158 36 : B.transpose();
159 :
160 36 : double alpha = -1.0 + 2.0 * static_cast<double> (std::rand()) / static_cast<double> (RAND_MAX);
161 36 : double gamma = -0.5 + static_cast<double> (std::rand()) / static_cast<double> (RAND_MAX);
162 36 : const double tol = 1e-8;
163 :
164 36 : int status = Matrix::add(A, alpha, B, gamma); // A = gamma * A + alpha * B
165 36 : _ASSERT_EQ(ForBESUtils::STATUS_OK, status);
166 :
167 36 : _ASSERT_EQ(n, A.getNrows());
168 882 : for (size_t i = 0; i < n; i++) {
169 18189 : for (size_t j = 0; j < m; j++) {
170 17343 : _ASSERT_NUM_EQ(gamma * A_copy.get(i, j) + alpha * B_copy.get(j, i), A.get(i, j), tol);
171 : }
172 : }
173 36 : }
174 : }
175 1 : }
176 :
177 1 : void TestMatrixExtras::test_add_SS() {
178 1 : size_t n = 10;
179 1 : size_t m = 15;
180 1 : size_t nnz = 100;
181 1 : size_t repetitions = 300;
182 :
183 1 : Matrix A = MatrixFactory::MakeRandomSparse(n, m, nnz, 2.0, 1.0);
184 301 : for (size_t r = 0; r < repetitions; r++) {
185 300 : Matrix B = MatrixFactory::MakeRandomSparse(n, m, nnz, 2.0, 1.0);
186 :
187 300 : double alpha = 2.0 * static_cast<double> (std::rand()) / RAND_MAX;
188 300 : double gamma = -0.5 + static_cast<double> (std::rand()) / RAND_MAX;
189 :
190 600 : Matrix A_copy(A);
191 600 : Matrix B_copy(B);
192 :
193 : // A = gamma*A + alpha*B
194 300 : int status = Matrix::add(A, alpha, B, gamma);
195 300 : _ASSERT_EQ(ForBESUtils::STATUS_OK, status);
196 300 : _ASSERT_EQ(B_copy, B);
197 :
198 300 : A_copy *= gamma;
199 300 : B_copy *= alpha;
200 600 : Matrix C = A_copy + B_copy;
201 :
202 300 : _ASSERT_EQ(C, A);
203 301 : }
204 1 : }
205 :
206 1 : void TestMatrixExtras::test_add_DDT() {
207 1 : size_t n = 10;
208 1 : size_t m = 15;
209 1 : size_t repetitions = 20;
210 :
211 21 : for (size_t r = 0; r < repetitions; r++) {
212 20 : Matrix A = MatrixFactory::MakeRandomMatrix(m, n, 2.0, 1.0);
213 40 : Matrix B = MatrixFactory::MakeRandomMatrix(n, m, 2.0, 1.0);
214 :
215 20 : double alpha = 2.0 * static_cast<float> (std::rand()) / static_cast<float> (RAND_MAX);
216 20 : double gamma = -0.5 + static_cast<float> (std::rand()) / static_cast<float> (RAND_MAX);
217 :
218 20 : A.transpose();
219 :
220 40 : Matrix A_copy(A);
221 40 : Matrix B_copy(B);
222 :
223 : // A = gamma*A + alpha*B
224 20 : int status = Matrix::add(A, alpha, B, gamma);
225 20 : _ASSERT_EQ(ForBESUtils::STATUS_OK, status);
226 20 : _ASSERT_EQ(B_copy, B);
227 :
228 20 : _ASSERT_EQ(n, A.getNrows());
229 20 : _ASSERT_EQ(m, A.getNcols());
230 :
231 20 : A_copy *= gamma;
232 20 : B_copy *= alpha;
233 40 : Matrix C = A_copy + B_copy;
234 :
235 20 : _ASSERT_EQ(C, A);
236 20 : }
237 1 : }
238 :
239 1 : void TestMatrixExtras::test_add_DTD() {
240 1 : size_t n = 10;
241 1 : size_t m = 15;
242 1 : size_t repetitions = 20;
243 :
244 21 : for (size_t r = 0; r < repetitions; r++) {
245 20 : Matrix A = MatrixFactory::MakeRandomMatrix(n, m, 2.0, 1.0);
246 40 : Matrix B = MatrixFactory::MakeRandomMatrix(m, n, 2.0, 1.0);
247 :
248 20 : double alpha = 2.0 * static_cast<float> (std::rand()) / static_cast<float> (RAND_MAX);
249 20 : double gamma = -0.5 + static_cast<float> (std::rand()) / static_cast<float> (RAND_MAX);
250 :
251 20 : B.transpose();
252 :
253 40 : Matrix A_copy(A);
254 40 : Matrix B_copy(B);
255 :
256 20 : int status = Matrix::add(A, alpha, B, gamma); // A = gamma*A + alpha*B'
257 20 : _ASSERT_EQ(ForBESUtils::STATUS_OK, status);
258 20 : _ASSERT_EQ(B_copy, B);
259 :
260 20 : _ASSERT_EQ(n, A.getNrows());
261 20 : _ASSERT_EQ(m, A.getNcols());
262 :
263 20 : A_copy *= gamma;
264 20 : B_copy *= alpha;
265 40 : Matrix C = A_copy + B_copy;
266 :
267 20 : _ASSERT_EQ(C, A);
268 20 : }
269 1 : }
270 :
271 1 : void TestMatrixExtras::test_add_DTDT() {
272 1 : size_t n = 10;
273 1 : size_t m = 15;
274 1 : size_t repetitions = 20;
275 :
276 21 : for (size_t r = 0; r < repetitions; r++) {
277 20 : Matrix A = MatrixFactory::MakeRandomMatrix(n, m, 2.0, 1.0);
278 40 : Matrix B = MatrixFactory::MakeRandomMatrix(n, m, 2.0, 1.0);
279 :
280 20 : double alpha = 2.0 * static_cast<float> (std::rand()) / static_cast<float> (RAND_MAX);
281 20 : double gamma = -0.5 + static_cast<float> (std::rand()) / static_cast<float> (RAND_MAX);
282 :
283 40 : Matrix A_copy(A);
284 40 : Matrix B_copy(B);
285 :
286 20 : A.transpose();
287 20 : B.transpose();
288 :
289 : // A := gamma*A + alpha*B
290 20 : int status = Matrix::add(A, alpha, B, gamma);
291 20 : _ASSERT_EQ(ForBESUtils::STATUS_OK, status);
292 20 : _ASSERT_EQ(m, A.getNrows());
293 20 : _ASSERT_EQ(n, A.getNcols());
294 :
295 320 : for (size_t i = 0; i < m; i++) {
296 3300 : for (size_t j = 0; j < n; j++) {
297 3000 : _ASSERT_EQ(gamma * A_copy.get(j, i) + alpha * B_copy.get(j, i), A.get(i, j));
298 : }
299 : }
300 :
301 :
302 20 : }
303 1 : }
304 :
305 1 : void TestMatrixExtras::test_add_SST() {
306 1 : size_t n = 10;
307 1 : size_t m = 15;
308 1 : size_t nnz = std::ceil(n * m / 2);
309 1 : size_t repetitions = 300;
310 :
311 301 : for (size_t r = 0; r < repetitions; r++) {
312 300 : Matrix A = MatrixFactory::MakeRandomSparse(n, m, nnz, 2.0, 1.0);
313 600 : Matrix B = MatrixFactory::MakeRandomSparse(m, n, nnz, 2.0, 1.0);
314 :
315 300 : double alpha = 2.0 * static_cast<float> (std::rand()) / static_cast<float> (RAND_MAX);
316 300 : double gamma = -0.5 + static_cast<float> (std::rand()) / static_cast<float> (RAND_MAX);
317 :
318 600 : Matrix B_copy(B);
319 300 : B.transpose();
320 600 : Matrix A_copy(A);
321 :
322 300 : int status = Matrix::add(A, alpha, B, gamma); // A = gamma*A + alpha*B'
323 300 : _ASSERT_EQ(ForBESUtils::STATUS_OK, status);
324 :
325 300 : _ASSERT_EQ(n, A.getNrows());
326 300 : _ASSERT_EQ(m, A.getNcols());
327 :
328 3300 : for (size_t i = 0; i < n; i++) {
329 48000 : for (size_t j = 0; j < m; j++) {
330 45000 : _ASSERT_NUM_EQ(gamma * A_copy.get(i, j) + alpha * B_copy.get(j, i), A.get(i, j), 1e-8);
331 : }
332 : }
333 :
334 600 : Matrix S(A_copy);
335 300 : S *= gamma;
336 600 : Matrix T(B);
337 300 : T *= alpha;
338 :
339 300 : S += T;
340 :
341 300 : _ASSERT_EQ(S, A);
342 :
343 300 : }
344 1 : }
345 :
346 1 : void TestMatrixExtras::test_add_STS() {
347 1 : size_t n = 10;
348 1 : size_t m = 15;
349 1 : size_t nnz = 100;
350 1 : size_t repetitions = 300;
351 :
352 301 : for (size_t r = 0; r < repetitions; r++) {
353 300 : Matrix B = MatrixFactory::MakeRandomSparse(n, m, nnz, 2.0, 1.0);
354 600 : Matrix A = MatrixFactory::MakeRandomSparse(m, n, nnz, 2.0, 1.0);
355 :
356 300 : double alpha = 2.0 * static_cast<float> (std::rand()) / static_cast<float> (RAND_MAX);
357 300 : double gamma = -0.5 + static_cast<float> (std::rand()) / static_cast<float> (RAND_MAX);
358 :
359 :
360 600 : Matrix A_copy(A);
361 300 : A.transpose();
362 :
363 300 : int status = Matrix::add(A, alpha, B, gamma); // A' = gamma*A' + alpha*B
364 300 : _ASSERT_EQ(ForBESUtils::STATUS_OK, status);
365 :
366 300 : _ASSERT_EQ(n, A.getNrows());
367 300 : _ASSERT_EQ(m, A.getNcols());
368 :
369 3300 : for (size_t i = 0; i < n; i++) {
370 48000 : for (size_t j = 0; j < m; j++) {
371 45000 : _ASSERT_EQ(gamma * A_copy.get(j, i) + alpha * B.get(i, j), A.get(i, j));
372 : }
373 : }
374 :
375 300 : }
376 1 : }
377 :
378 1 : void TestMatrixExtras::test_add_STST() {
379 1 : size_t repetitions = 300;
380 :
381 301 : for (size_t r = 0; r < repetitions; r++) {
382 300 : size_t n = 10;
383 300 : size_t m = 15;
384 300 : size_t nnz = 70;
385 300 : Matrix A = MatrixFactory::MakeRandomSparse(n, m, nnz, 2.0, 1.0);
386 600 : Matrix B = MatrixFactory::MakeRandomSparse(n, m, nnz, 2.0, 1.0);
387 :
388 300 : double alpha = 2.0 * static_cast<float> (std::rand()) / static_cast<float> (RAND_MAX);
389 300 : double gamma = -0.5 + static_cast<float> (std::rand()) / static_cast<float> (RAND_MAX);
390 :
391 :
392 600 : Matrix A_copy(A);
393 600 : Matrix B_copy(B);
394 300 : A.transpose();
395 300 : B.transpose();
396 :
397 300 : int status = Matrix::add(A, alpha, B, gamma); // A' = gamma*A' + alpha*B
398 300 : _ASSERT_EQ(ForBESUtils::STATUS_OK, status);
399 :
400 300 : _ASSERT_EQ(m, A.getNrows());
401 300 : _ASSERT_EQ(n, A.getNcols());
402 4800 : for (size_t i = 0; i < m; i++) {
403 49500 : for (size_t j = 0; j < n; j++) {
404 45000 : _ASSERT_EQ(gamma * A_copy.get(j, i) + alpha * B_copy.get(j, i), A.get(i, j));
405 : }
406 : }
407 :
408 300 : }
409 1 : }
410 :
411 : /***** MULTIPLICATION ******/
412 :
413 1 : void TestMatrixExtras::test_mult_DD() {
414 1 : size_t n = 8;
415 1 : size_t k = 6;
416 1 : size_t m = 5;
417 1 : size_t repetitions = 300;
418 :
419 301 : for (size_t r = 0; r < repetitions; r++) {
420 300 : Matrix A = MatrixFactory::MakeRandomMatrix(n, k, 0.0, 1.0);
421 600 : Matrix B = MatrixFactory::MakeRandomMatrix(k, m, 0.0, 1.0);
422 600 : Matrix C = MatrixFactory::MakeRandomMatrix(n, m, 0.0, 1.0);
423 600 : Matrix C_copy(C);
424 :
425 300 : double alpha = 2.0 * static_cast<float> (std::rand()) / static_cast<float> (RAND_MAX);
426 300 : double gamma = -0.5 + static_cast<float> (std::rand()) / static_cast<float> (RAND_MAX);
427 :
428 300 : int status = Matrix::mult(C, alpha, A, B, gamma);
429 300 : _ASSERT_EQ(ForBESUtils::STATUS_OK, status);
430 :
431 600 : Matrix AB = A*B;
432 300 : AB *= alpha;
433 300 : C_copy *= gamma;
434 300 : C_copy += AB;
435 :
436 300 : _ASSERT_EQ(C_copy, C);
437 300 : }
438 1 : }
439 :
440 1 : void TestMatrixExtras::test_mult_DH() {
441 1 : size_t n = 8;
442 1 : size_t k = 6;
443 1 : size_t repetitions = 100;
444 :
445 101 : for (size_t r = 0; r < repetitions; r++) {
446 100 : Matrix A = MatrixFactory::MakeRandomMatrix(n, k, 0.0, 1.0);
447 200 : Matrix B = MatrixFactory::MakeRandomMatrix(k, k, 0.0, 1.0, Matrix::MATRIX_SYMMETRIC);
448 200 : Matrix C = MatrixFactory::MakeRandomMatrix(n, k, 0.0, 1.0);
449 200 : Matrix C_copy(C);
450 :
451 100 : double alpha = 2.0 * static_cast<float> (std::rand()) / static_cast<float> (RAND_MAX);
452 100 : double gamma = -0.5 + static_cast<float> (std::rand()) / static_cast<float> (RAND_MAX);
453 :
454 100 : int status = Matrix::mult(C, alpha, A, B, gamma);
455 100 : _ASSERT_EQ(ForBESUtils::STATUS_OK, status);
456 :
457 200 : Matrix AB = A*B;
458 100 : AB *= alpha;
459 100 : C_copy *= gamma;
460 100 : _ASSERT_NOT(C_copy == C);
461 100 : C_copy += AB;
462 :
463 100 : _ASSERT_EQ(C_copy, C);
464 100 : }
465 1 : }
466 :
467 1 : void TestMatrixExtras::test_mult_DDT() {
468 1 : size_t n = 8;
469 1 : size_t k = 6;
470 1 : size_t m = 5;
471 1 : size_t repetitions = 300;
472 :
473 301 : for (size_t r = 0; r < repetitions; r++) {
474 300 : Matrix A = MatrixFactory::MakeRandomMatrix(n, k, 0.0, 1.0);
475 600 : Matrix B = MatrixFactory::MakeRandomMatrix(m, k, 0.0, 1.0);
476 600 : Matrix C = MatrixFactory::MakeRandomMatrix(n, m, 0.0, 1.0);
477 600 : Matrix C_copy(C);
478 :
479 300 : B.transpose();
480 :
481 300 : double alpha = 2.0 * static_cast<float> (std::rand()) / static_cast<float> (RAND_MAX);
482 300 : double gamma = -0.5 + static_cast<float> (std::rand()) / static_cast<float> (RAND_MAX);
483 :
484 300 : int status = Matrix::mult(C, alpha, A, B, gamma);
485 300 : _ASSERT_EQ(ForBESUtils::STATUS_OK, status);
486 :
487 600 : Matrix AB = A*B;
488 300 : AB *= alpha;
489 300 : C_copy *= gamma;
490 300 : _ASSERT_NOT(C_copy == C);
491 300 : C_copy += AB;
492 :
493 300 : _ASSERT_EQ(C_copy, C);
494 300 : }
495 1 : }
496 :
497 1 : void TestMatrixExtras::test_mult_DX() {
498 1 : size_t n = 8;
499 1 : size_t k = 6;
500 1 : size_t repetitions = 300;
501 :
502 301 : for (size_t r = 0; r < repetitions; r++) {
503 300 : Matrix A = MatrixFactory::MakeRandomMatrix(n, k, 0.0, 1.0);
504 600 : Matrix B = MatrixFactory::MakeRandomMatrix(k, k, 0.0, 1.0, Matrix::MATRIX_DIAGONAL);
505 600 : Matrix C = MatrixFactory::MakeRandomMatrix(n, k, 0.0, 1.0);
506 600 : Matrix C_copy(C);
507 :
508 300 : double alpha = 2.0 * static_cast<float> (std::rand()) / static_cast<float> (RAND_MAX);
509 300 : double gamma = -0.5 + static_cast<float> (std::rand()) / static_cast<float> (RAND_MAX);
510 :
511 300 : int status = Matrix::mult(C, alpha, A, B, gamma);
512 300 : _ASSERT_EQ(ForBESUtils::STATUS_OK, status);
513 :
514 600 : Matrix AB = A*B;
515 300 : AB *= alpha;
516 300 : C_copy *= gamma;
517 300 : C_copy += AB;
518 :
519 300 : _ASSERT_EQ(C_copy, C);
520 300 : }
521 1 : }
522 :
523 1 : void TestMatrixExtras::test_mult_SS() {
524 1 : size_t n = 10;
525 1 : size_t k = 8;
526 1 : size_t m = 9;
527 :
528 1 : size_t repetitions = 300;
529 :
530 301 : for (size_t r = 0; r < repetitions; r++) {
531 :
532 300 : Matrix A = MatrixFactory::MakeRandomSparse(n, k, std::ceil(n * k / 2), 2.0, 1.0);
533 600 : Matrix B = MatrixFactory::MakeRandomSparse(k, m, std::ceil(k * m / 2), 2.0, 1.0);
534 600 : Matrix C = MatrixFactory::MakeRandomSparse(n, m, std::ceil(n * m / 2), 2.0, 1.0);
535 :
536 300 : double alpha = -2.0 + 2.0 * static_cast<double> (std::rand()) / static_cast<double> (RAND_MAX);
537 300 : double gamma = -2.0 + 3.0 * static_cast<double> (std::rand()) / static_cast<double> (RAND_MAX);
538 :
539 :
540 600 : Matrix aAB = A*B;
541 300 : aAB *= alpha; /* aAB = a * A * B */
542 :
543 600 : Matrix gC(C);
544 300 : gC *= gamma; /* gC = g * C */
545 :
546 600 : Matrix R = aAB + gC;
547 :
548 600 : Matrix C_copy(C);
549 300 : int status = Matrix::mult(C_copy, alpha, A, B, gamma);
550 300 : _ASSERT_EQ(ForBESUtils::STATUS_OK, status);
551 :
552 300 : _ASSERT_EQ(R, C_copy);
553 300 : }
554 :
555 1 : }
556 :
557 1 : void TestMatrixExtras::test_mult_SS2() {
558 : // testing with gamma = 0.0;
559 1 : size_t repetitions = 300;
560 301 : for (size_t r = 0; r < repetitions; r++) {
561 300 : size_t n = 10;
562 300 : size_t k = 8;
563 300 : size_t m = 9;
564 300 : Matrix A = MatrixFactory::MakeRandomSparse(n, k, std::ceil(n * k / 2), 2.0, 1.0);
565 600 : Matrix B = MatrixFactory::MakeRandomSparse(k, m, std::ceil(k * m / 2), 2.0, 1.0);
566 600 : Matrix C = MatrixFactory::MakeRandomSparse(n, m, std::ceil(n * m / 2), 2.0, 1.0);
567 :
568 300 : double alpha = -2.0 + 2.0 * static_cast<double> (std::rand()) / static_cast<double> (RAND_MAX);
569 300 : double gamma = 0.0;
570 :
571 :
572 600 : Matrix aAB = A*B;
573 300 : aAB *= alpha; /* aAB = a * A * B */
574 :
575 :
576 600 : Matrix C_copy(C);
577 300 : int status = Matrix::mult(C_copy, alpha, A, B, gamma); /* C_copy = a * A * B */
578 300 : _ASSERT_EQ(ForBESUtils::STATUS_OK, status);
579 300 : _ASSERT_EQ(aAB, C_copy);
580 300 : }
581 1 : }
582 :
583 1 : void TestMatrixExtras::test_mult_SS3() {
584 : // testing with C : DENSE
585 :
586 1 : size_t repetitions = 300;
587 :
588 301 : for (size_t r = 0; r < repetitions; r++) {
589 300 : size_t n = 10;
590 300 : size_t k = 8;
591 300 : size_t m = 9;
592 300 : Matrix A = MatrixFactory::MakeRandomSparse(n, k, std::ceil(n * k / 2), 2.0, 1.0);
593 600 : Matrix B = MatrixFactory::MakeRandomSparse(k, m, std::ceil(m * k / 2), 2.0, 1.0);
594 600 : Matrix C = MatrixFactory::MakeRandomMatrix(n, m, 2.0, 1.0);
595 :
596 300 : double alpha = -2.0 + 2.0 * static_cast<double> (std::rand()) / static_cast<double> (RAND_MAX);
597 300 : double gamma = -2.0 + 3.0 * static_cast<double> (std::rand()) / static_cast<double> (RAND_MAX);
598 :
599 :
600 600 : Matrix aAB = A*B;
601 300 : aAB *= alpha; /* aAB = a * A * B */
602 :
603 600 : Matrix gC(C);
604 300 : gC *= gamma; /* gC = g * C */
605 :
606 600 : Matrix R = aAB + gC;
607 :
608 600 : Matrix C_copy(C);
609 300 : int status = Matrix::mult(C_copy, alpha, A, B, gamma);
610 300 : _ASSERT_EQ(ForBESUtils::STATUS_OK, status);
611 :
612 300 : _ASSERT_EQ(R, C_copy);
613 300 : }
614 1 : }
615 :
616 1 : void TestMatrixExtras::test_mult_SST() {
617 1 : size_t n = 10;
618 1 : size_t k = 8;
619 1 : size_t m = 9;
620 :
621 1 : size_t repetitions = 300;
622 :
623 301 : for (size_t r = 0; r < repetitions; r++) {
624 300 : Matrix A = MatrixFactory::MakeRandomSparse(n, k, std::ceil(n * k / 2), 2.0, 1.0);
625 600 : Matrix B = MatrixFactory::MakeRandomSparse(m, k, std::ceil(m * k / 2), 2.0, 1.0);
626 600 : Matrix C = MatrixFactory::MakeRandomSparse(n, m, std::ceil(n * m / 2), 2.0, 1.0);
627 :
628 300 : B.transpose();
629 :
630 300 : double alpha = -2.0 + 2.0 * static_cast<double> (std::rand()) / static_cast<double> (RAND_MAX);
631 300 : double gamma = -2.0 + 3.0 * static_cast<double> (std::rand()) / static_cast<double> (RAND_MAX);
632 :
633 600 : Matrix C_orig(C);
634 :
635 300 : int status = Matrix::mult(C, alpha, A, B, gamma);
636 300 : _ASSERT_EQ(ForBESUtils::STATUS_OK, status);
637 :
638 300 : C_orig *= gamma;
639 600 : Matrix AB = A * B;
640 300 : AB *= alpha;
641 300 : C_orig += AB;
642 :
643 300 : _ASSERT_EQ(C_orig, C);
644 300 : }
645 1 : }
646 :
647 1 : void TestMatrixExtras::test_mult_SX() {
648 :
649 1 : size_t repetitions = 300;
650 301 : for (size_t r = 0; r < repetitions; r++) {
651 300 : size_t n = 10;
652 300 : size_t k = 8;
653 300 : size_t nnz = std::ceil(n * k / 2);
654 300 : Matrix A = MatrixFactory::MakeRandomSparse(n, k, nnz, 2.0, 1.0);
655 600 : Matrix B = MatrixFactory::MakeRandomMatrix(k, k, 0.0, 1.0, Matrix::MATRIX_DIAGONAL);
656 600 : Matrix C = MatrixFactory::MakeRandomSparse(n, k, nnz, 2.0, 1.0);
657 600 : Matrix D = MatrixFactory::MakeRandomMatrix(n, k, 0.0, 1.0);
658 :
659 300 : double alpha = -2.0 + 2.0 * static_cast<double> (std::rand()) / static_cast<double> (RAND_MAX);
660 300 : double gamma = -2.0 + 3.0 * static_cast<double> (std::rand()) / static_cast<double> (RAND_MAX);
661 :
662 600 : Matrix D_orig(D);
663 600 : Matrix C_orig(C);
664 :
665 300 : int status = Matrix::mult(C, alpha, A, B, gamma);
666 300 : _ASSERT_EQ(ForBESUtils::STATUS_OK, status);
667 :
668 300 : status = Matrix::mult(D, alpha, A, B, gamma);
669 300 : _ASSERT_EQ(ForBESUtils::STATUS_OK, status);
670 :
671 300 : C_orig *= gamma;
672 300 : D_orig *= gamma;
673 600 : Matrix AB = A * B;
674 300 : AB *= alpha;
675 300 : C_orig += AB;
676 300 : D_orig += AB;
677 :
678 300 : _ASSERT_EQ(C_orig, C);
679 300 : _ASSERT_EQ(D_orig, D);
680 :
681 300 : }
682 1 : }
683 :
684 1 : void TestMatrixExtras::test_mult_Hv() {
685 1 : size_t n = 10;
686 1 : Matrix x = MatrixFactory::MakeRandomMatrix(n, 1, 0.0, 1.0);
687 2 : Matrix A = MatrixFactory::MakeRandomMatrix(n, n, 0.0, 1.0, Matrix::MATRIX_SYMMETRIC);
688 2 : Matrix y = A*x;
689 :
690 2 : Matrix y2(n, 1);
691 1 : Matrix::mult(y2, 1.0, A, x, 0.0);
692 :
693 2 : _ASSERT_EQ(y2, y);
694 4 : }
|