ARTS  2.4.0(git:4fb77825)
test_matpack.cc
Go to the documentation of this file.
1 /* Copyright (C) 2001-2012 Stefan Buehler <sbuehler@ltu.se>
2 
3  This program is free software; you can redistribute it and/or modify it
4  under the terms of the GNU General Public License as published by the
5  Free Software Foundation; either version 2, or (at your option) any
6  later version.
7 
8  This program is distributed in the hope that it will be useful,
9  but WITHOUT ANY WARRANTY; without even the implied warranty of
10  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
11  GNU General Public License for more details.
12 
13  You should have received a copy of the GNU General Public License
14  along with this program; if not, write to the Free Software
15  Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307,
16  USA. */
17 
18 #include <algorithm>
19 #include <cmath>
20 #include <cstdlib>
21 #include <iostream>
22 #include "array.h"
23 #include "describe.h"
24 #include "exceptions.h"
25 #include "logic.h"
26 #include "math_funcs.h"
27 #include "matpackII.h"
28 #include "matpackVII.h"
29 #include "mystring.h"
30 #include "rational.h"
31 #include "test_utils.h"
32 #include "time.h"
33 #include "wigner_functions.h"
34 
35 using std::cout;
36 using std::endl;
37 using std::runtime_error;
38 
39 Numeric by_reference(const Numeric& x) { return x + 1; }
40 
41 Numeric by_value(Numeric x) { return x + 1; }
42 
43 void fill_with_junk(VectorView x) { x = 999; }
44 
45 void fill_with_junk(MatrixView x) { x = 888; }
46 
47 int test1() {
48  Vector v(20);
49 
50  cout << "v.nelem() = " << v.nelem() << "\n";
51 
52  for (Index i = 0; i < v.nelem(); ++i) v[i] = (Numeric)i;
53 
54  cout << "v.begin() = " << *v.begin() << "\n";
55 
56  cout << "v = \n" << v << "\n";
57 
58  fill_with_junk(v[Range(1, 8, 2)][Range(2, joker)]);
59  // fill_with_junk(v);
60 
61  Vector v2 = v[Range(2, 4)];
62 
63  cout << "v2 = \n" << v2 << "\n";
64 
65  for (Index i = 0; i < 1000; ++i) {
66  Vector v3(1000);
67  v3 = (Numeric)i;
68  }
69 
70  v2[Range(joker)] = 88;
71 
72  v2[Range(0, 2)] = 77;
73 
74  cout << "v = \n" << v << "\n";
75  cout << "v2 = \n" << v2 << "\n";
76  cout << "v2.nelem() = \n" << v2.nelem() << "\n";
77 
78  Vector v3;
79  v3.resize(v2.nelem());
80  v3 = v2;
81 
82  cout << "\nv3 = \n" << v3 << "\n";
84  cout << "\nv3 after junking v2 = \n" << v3 << "\n";
85  v3 *= 2;
86  cout << "\nv3 after *2 = \n" << v3 << "\n";
87 
88  Matrix M(10, 15);
89  {
90  Numeric n = 0;
91  for (Index i = 0; i < M.nrows(); ++i)
92  for (Index j = 0; j < M.ncols(); ++j) M(i, j) = ++n;
93  }
94 
95  cout << "\nM =\n" << M << "\n";
96 
97  cout << "\nM(Range(2,4),Range(2,4)) =\n"
98  << M(Range(2, 4), Range(2, 4)) << "\n";
99 
100  cout << "\nM(Range(2,4),Range(2,4))(Range(1,2),Range(1,2)) =\n"
101  << M(Range(2, 4), Range(2, 4))(Range(1, 2), Range(1, 2)) << "\n";
102 
103  cout << "\nM(1,Range(joker)) =\n" << M(1, Range(joker)) << "\n";
104 
105  cout << "\nFilling M(1,Range(1,2)) with junk.\n";
106  fill_with_junk(M(1, Range(1, 2)));
107 
108  cout << "\nM(Range(0,4),Range(0,4)) =\n"
109  << M(Range(0, 4), Range(0, 4)) << "\n";
110 
111  cout << "\nFilling M(Range(4,2,2),Range(6,3)) with junk.\n";
112 
113  MatrixView s = M(Range(4, 2, 2), Range(6, 3));
114  fill_with_junk(s);
115 
116  cout << "\nM =\n" << M << "\n";
117 
118  const Matrix C = M;
119 
120  cout << "\nC(Range(3,4,2),Range(2,3,3)) =\n"
121  << C(Range(3, 4, 2), Range(2, 3, 3)) << "\n";
122 
123  cout << "\nC(Range(3,4,2),Range(2,3,3)).transpose() =\n"
124  << transpose(C(Range(3, 4, 2), Range(2, 3, 3))) << "\n";
125 
126  return 0;
127 }
128 
129 void test2() {
130  Vector v(50000000);
131 
132  cout << "v.nelem() = " << v.nelem() << "\n";
133 
134  cout << "Filling\n";
135  // for (Index i=0; i<v.nelem(); ++i )
136  // v[i] = sqrt(i);
137  v = 1.;
138  cout << "Done\n";
139 }
140 
141 void test4() {
142  Vector a(10);
143  Vector b(a.nelem());
144 
145  for (Index i = 0; i < a.nelem(); ++i) {
146  a[i] = (Numeric)(i + 1);
147  b[i] = (Numeric)(a.nelem() - i);
148  }
149 
150  cout << "a = \n" << a << "\n";
151  cout << "b = \n" << b << "\n";
152  cout << "a*b \n= " << a * b << "\n";
153 
154  Matrix A(11, 6);
155  Matrix B(10, 20);
156  Matrix C(20, 5);
157 
158  B = 2;
159  C = 3;
160  mult(A(Range(1, joker), Range(1, joker)), B, C);
161 
162  // cout << "\nB =\n" << B << "\n";
163  // cout << "\nC =\n" << C << "\n";
164  cout << "\nB*C =\n" << A << "\n";
165 }
166 
167 void test5() {
168  Vector a(10);
169  Vector b(20);
170  Matrix M(10, 20);
171 
172  // Fill b and M with a constant number:
173  b = 1;
174  M = 2;
175 
176  cout << "b = \n" << b << "\n";
177  cout << "M =\n" << M << "\n";
178 
179  mult(a, M, b); // a = M*b
180  cout << "\na = M*b = \n" << a << "\n";
181 
182  mult(transpose((MatrixView)b), transpose((MatrixView)a), M); // b^t = a^t * M
183  cout << "\nb^t = a^t * M = \n" << transpose((MatrixView)b) << "\n";
184 }
185 
186 void test6() {
187  Index n = 5000;
188  Vector x(1, n, 1), y(n);
189  Matrix M(n, n);
190  M = 1;
191  // cout << "x = \n" << x << "\n";
192 
193  cout << "Transforming.\n";
194  // transform(x,sin,x);
195  // transform(transpose(y),sin,transpose(x));
196  // cout << "sin(x) =\n" << y << "\n";
197  for (Index i = 0; i < 1000; ++i) {
198  // mult(y,M,x);
199  transform(y, sin, static_cast<MatrixView>(x));
200  x += 1;
201  }
202  // cout << "y =\n" << y << "\n";
203 
204  cout << "Done.\n";
205 }
206 
207 void test7() {
208  Vector x(1, 20000000, 1);
209  Vector y(x.nelem());
210  transform(y, sin, x);
211  cout << "min(sin(x)), max(sin(x)) = " << min(y) << ", " << max(y) << "\n";
212 }
213 
214 void test8() {
215  Vector x(80000000);
216  for (Index i = 0; i < x.nelem(); ++i) x[i] = (Numeric)i;
217  cout << "Done."
218  << "\n";
219 }
220 
221 void test9() {
222  // Initialization of Matrix with view of other Matrix:
223  Matrix A(4, 8);
224  Matrix B(A(Range(joker), Range(0, 3)));
225  cout << "B = " << B << "\n";
226 }
227 
228 void test10() {
229  // Initialization of Matrix with a vector (giving a 1 column Matrix).
230 
231  // At the moment doing this with a non-const Vector will result in a
232  // warning message.
233  Vector v(1, 8, 1);
234  Matrix M((const Vector)v);
235  cout << "M = " << M << "\n";
236 }
237 
238 void test11() {
239  // Assignment between Vector and Matrix:
240 
241  // At the moment doing this with a non-const Vector will result in a
242  // warning message.
243  Vector v(1, 8, 1);
244  Matrix M(v.nelem(), 1);
245  M = v;
246  cout << "M = " << M << "\n";
247 }
248 
249 void test12() {
250  // Copying of Arrays
251 
252  Array<String> sa(3);
253  sa[0] = "It's ";
254  sa[1] = "a ";
255  sa[2] = "test.";
256 
257  Array<String> sb(sa), sc(sa.nelem());
258 
259  cout << "sb = \n" << sb << "\n";
260 
261  sc = sa;
262 
263  cout << "sc = \n" << sc << "\n";
264 }
265 
266 void test13() {
267  // Mix vector and one-column matrix in += operator.
268  const Vector v(1, 8, 1); // The const is necessary here to
269  // avoid compiler warnings about
270  // different conversion paths.
271  Matrix M(v);
272  M += v;
273  cout << "M = \n" << M << "\n";
274 }
275 
276 void test14() {
277  // Test explicit Array constructors.
278  Array<String> a{"Test"};
279  Array<Index> b{1, 2};
280  Array<Numeric> c{1, 2, 3, 4, 5, 6, 7, 8, 9, 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 0};
281  cout << "a = \n" << a << "\n";
282  cout << "b = \n" << b << "\n";
283  cout << "c = \n" << c << "\n";
284 }
285 
286 void test15() {
287  // Test String.
288  String a = "Nur ein Test.";
289  cout << "a = " << a << "\n";
290  String b(a, 5, -1);
291  cout << "b = " << b << "\n";
292 }
293 
294 /*void test16()
295 {
296  // Test interaction between Array<Numeric> and Vector.
297  Vector a;
298  Array<Numeric> b;
299  b.push_back(1);
300  b.push_back(2);
301  b.push_back(3);
302  a.resize(b.nelem());
303  a = b;
304  cout << "b =\n" << b << "\n";
305  cout << "a =\n" << a << "\n";
306 }*/
307 
308 void test17() {
309  // Test Sum.
310  Vector a(1, 10, 1);
311  cout << "a.sum() = " << a.sum() << "\n";
312 }
313 
314 void test18() {
315  // Test elementvise square of a vector:
316  Vector a(1, 10, 1);
317  a *= a;
318  cout << "a *= a =\n" << a << "\n";
319 }
320 
321 void test19() {
322  // There exists no explicit filling constructor of the form
323  // Vector a(3,1.7).
324  // But you can use the more general filling constructor with 3 arguments.
325 
326  Vector a(1, 10, 1);
327  Vector b(5.3, 10, 0);
328  cout << "a =\n" << a << "\n";
329  cout << "b =\n" << b << "\n";
330 }
331 
332 void test20() {
333  // Test initialization list constructor:
334  Vector a{1, 2, 3, 4, 5, 6, 7, 8, 9, 10};
335  cout << "a =\n" << a << "\n";
336 }
337 
338 void test21() {
339  Numeric s = 0;
340  // Test speed of call by reference:
341  cout << "By reference:\n";
342  for (Index i = 0; i < (Index)1e8; ++i) {
343  s += by_reference(s);
344  s -= by_reference(s);
345  }
346  cout << "s = " << s << "\n";
347 }
348 
349 void test22() {
350  Numeric s = 0;
351  // Test speed of call by value:
352  cout << "By value:\n";
353  for (Index i = 0; i < (Index)1e8; ++i) {
354  s += by_value(s);
355  s -= by_value(s);
356  }
357  cout << "s = " << s << "\n";
358 }
359 
360 void test23() {
361  // Test constructors that fills with constant:
362  Vector a(10, 3.5);
363  cout << "a =\n" << a << "\n";
364  Matrix b(10, 10, 4.5);
365  cout << "b =\n" << b << "\n";
366 }
367 
368 void test24() {
369  // Try element-vise multiplication of Matrix and Vector:
370  Matrix a(5, 1, 2.5);
371  Vector b(1, 5, 1);
372  a *= b;
373  cout << "a*=b =\n" << a << "\n";
374  a /= b;
375  cout << "a/=b =\n" << a << "\n";
376  a += b;
377  cout << "a+=b =\n" << a << "\n";
378  a -= b;
379  cout << "a-=b =\n" << a << "\n";
380 }
381 
382 void test25() {
383  // Test min and max for Array:
384  Array<Index> a{1, 2, 3, 4, 5, 6, 5, 4, 3, 2, 1};
385  cout << "min/max of a = " << min(a) << "/" << max(a) << "\n";
386 }
387 
388 void test26() {
389  cout << "Test filling constructor for Array:\n";
390  Array<String> a(4, "Hello");
391  cout << "a =\n" << a << "\n";
392 }
393 
394 void test27() {
395  cout << "Test Arrays of Vectors:\n";
396  Array<Vector> a;
397  a.push_back({1.0, 2.0});
398  a.push_back(Vector(1.0, 10, 1.0));
399  cout << "a =\n" << a << "\n";
400 }
401 
402 void test28() {
403  cout << "Test default constructor for Matrix:\n";
404  Matrix a;
405  Matrix b(a);
406  cout << "b =\n" << b << "\n";
407 }
408 
409 void test29() {
410  cout << "Test Arrays of Matrix:\n";
411  ArrayOfMatrix a;
412  Matrix b;
413 
414  b.resize(2, 2);
415  b(0, 0) = 1;
416  b(0, 1) = 2;
417  b(1, 0) = 3;
418  b(1, 1) = 4;
419  a.push_back(b);
420  b *= 2;
421  a.push_back(b);
422 
423  a[0].resize(2, 3);
424  a[0] = 4;
425 
426  a.resize(3);
427  a[2].resize(4, 5);
428  a[2] = 5;
429 
430  cout << "a =\n" << a << "\n";
431 }
432 
433 void test30() {
434  cout << "Test Matrices of size 0:\n";
435  Matrix a(0, 0);
436  // cout << "a(0,0) =\n" << a(0,0) << "\n";
437  a.resize(2, 2);
438  a = 1;
439  cout << "a =\n" << a << "\n";
440 
441  Matrix b(3, 0);
442  // cout << "b(0,0) =\n" << b(0,0) << "\n";
443  b.resize(b.nrows(), b.ncols() + 3);
444  b = 2;
445  cout << "b =\n" << b << "\n";
446 
447  Matrix c(0, 3);
448  // cout << "c(0,0) =\n" << c(0,0) << "\n";
449  c.resize(c.nrows() + 3, c.ncols());
450  c = 3;
451  cout << "c =\n" << c << "\n";
452 }
453 
454 void test31() {
455  cout << "Test Tensor3:\n";
456 
457  Tensor3 a(2, 3, 4, 1.0);
458 
459  Index fill = 0;
460 
461  // Fill with some numbers
462  for (Index i = 0; i < a.npages(); ++i)
463  for (Index j = 0; j < a.nrows(); ++j)
464  for (Index k = 0; k < a.ncols(); ++k) a(i, j, k) = (Numeric)(++fill);
465 
466  cout << "a =\n" << a << "\n";
467 
468  cout << "Taking out first row of first page:\n"
469  << a(0, 0, Range(joker)) << "\n";
470 
471  cout << "Taking out last column of second page:\n"
472  << a(1, Range(joker), a.ncols() - 1) << "\n";
473 
474  cout << "Taking out the first letter on every page:\n"
475  << a(Range(joker), 0, 0) << "\n";
476 
477  cout << "Taking out first page:\n"
478  << a(0, Range(joker), Range(joker)) << "\n";
479 
480  cout << "Taking out last row of all pages:\n"
481  << a(Range(joker), a.nrows() - 1, Range(joker)) << "\n";
482 
483  cout << "Taking out second column of all pages:\n"
484  << a(Range(joker), Range(joker), 1) << "\n";
485 
486  a *= 2;
487 
488  cout << "After element-vise multiplication with 2:\n" << a << "\n";
489 
490  transform(a, sqrt, a);
491 
492  cout << "After taking the square-root:\n" << a << "\n";
493 
494  Index s = 200;
495  cout << "Let's allocate a large tensor, "
496  << (Numeric)(s * s * s * 8) / 1024. / 1024. << " MB...\n";
497 
498  a.resize(s, s, s);
499 
500  cout << "Set it to 1...\n";
501 
502  a = 1;
503 
504  cout << "a(900,900,900) = " << a(90, 90, 90) << "\n";
505 
506  fill = 0;
507 
508  cout << "Fill with running numbers, using for loops...\n";
509  for (Index i = 0; i < a.npages(); ++i)
510  for (Index j = 0; j < a.nrows(); ++j)
511  for (Index k = 0; k < a.ncols(); ++k) a(i, j, k) = (Numeric)(++fill);
512 
513  cout << "Max(a) = ...\n";
514 
515  cout << max(a) << "\n";
516 }
517 
518 void test32() {
519  cout << "Test von X = A*X:\n";
520  Matrix X(3, 3), A(3, 3), B(3, 3);
521 
522  for (Index j = 0; j < A.nrows(); ++j)
523  for (Index k = 0; k < A.ncols(); ++k) {
524  X(j, k) = 1;
525  A(j, k) = (Numeric)(j + k);
526  }
527  cout << "A:\n" << A << "\n";
528  cout << "X:\n" << X << "\n";
529 
530  mult(B, A, X);
531  cout << "B = A*X:\n" << B << "\n";
532  mult(X, A, X);
533  cout << "X = A*X:\n" << X << "\n";
534 
535  cout << "This is not the same, and should not be, because you\n"
536  << "are not allowed to use the same matrix as input and output\n"
537  << "for mult!\n";
538 }
539 
540 void test33() {
541  cout << "Making things look bigger than they are...\n";
542 
543  {
544  cout << "1. Make a scalar look like a vector:\n";
545  Numeric a = 3.1415; // Just any number here.
546  VectorView av(a);
547  cout << "a, viewed as a vector: " << av << "\n";
548  cout << "Describe a: " << describe(a) << "\n";
549  av[0] += 1;
550  cout << "a, after the first element\n"
551  << "of the vector has been increased by 1: " << a << "\n";
552  }
553 
554  {
555  cout
556  << "\n2. Make a vector look like a matrix:\n"
557  << "This is an exception, because the new dimension is added at the end.\n";
558  Vector a{1, 2, 3, 4, 5};
559  MatrixView am = a;
560  cout << "a, viewed as a matrix:\n" << am << "\n";
561  cout << "Trasnpose view:\n" << transpose(am) << "\n";
562  cout << "Describe transpose(am): " << describe(transpose(am)) << "\n";
563 
564  cout << "\n3. Make a matrix look like a tensor3:\n";
565  Tensor3View at3 = am;
566  cout << "at3 = \n" << at3 << "\n";
567  cout << "Describe at3: " << describe(at3) << "\n";
568  at3(0, 2, 0) += 1;
569  cout << "a after Increasing element at3(0,2,0) by 1: \n" << a << "\n\n";
570 
571  Tensor4View at4 = at3;
572  cout << "at4 = \n" << at4 << "\n";
573  cout << "Describe at4: " << describe(at4) << "\n";
574 
575  Tensor5View at5 = at4;
576  cout << "at5 = \n" << at5 << "\n";
577  cout << "Describe at5: " << describe(at5) << "\n";
578 
579  Tensor6View at6 = at5;
580  cout << "at6 = \n" << at6 << "\n";
581  cout << "Describe at6: " << describe(at6) << "\n";
582 
583  Tensor7View at7 = at6;
584  cout << "at7 = \n" << at7 << "\n";
585  cout << "Describe at7: " << describe(at7) << "\n";
586 
587  at7(0, 0, 0, 0, 0, 2, 0) -= 1;
588 
589  cout << "After subtracting 1 from at7(0,0,0,0,0,2,0)\n"
590  << "a = " << a << "\n";
591 
592  cout << "\nAll in one go:\n";
593  Numeric b = 3.1415; // Just any number here.
594  Tensor7View bt7 = Tensor6View(
596  cout << "bt7:\n" << bt7 << "\n";
597  cout << "Describe bt7: " << describe(bt7) << "\n";
598  }
599 }
600 
601 void junk4(Tensor4View a) { cout << "Describe a: " << describe(a) << "\n"; }
602 
603 void junk2(ConstVectorView a) { cout << "Describe a: " << describe(a) << "\n"; }
604 
605 void test34() {
606  cout << "Test, if dimension expansion works implicitly.\n";
607 
608  Tensor3 t3(2, 3, 4);
609  junk4(t3);
610 
611  Numeric x;
613 }
614 
615 void test35() {
616  cout << "Test the new copy semantics.\n";
617  Vector a(1, 4, 1);
618  Vector b;
619 
620  b = a;
621  cout << "b = " << b << "\n";
622 
623  Vector aa(1, 5, 1);
624  ConstVectorView c = aa;
625  b = c;
626  cout << "b = " << b << "\n";
627 
628  Vector aaa(1, 6, 1);
629  VectorView d = aaa;
630  b = d;
631  cout << "b = " << b << "\n";
632 }
633 
634 void test36() {
635  cout << "Test using naked joker on Vector.\n";
636  Vector a(1, 4, 1);
637  VectorView b = a[joker];
638  cout << "a = " << a << "\n";
639  cout << "b = " << b << "\n";
640 }
641 
642 void test37(const Index& i) {
643  Vector v1(5e-15, 10, 0.42e-15 / 11);
644  Vector v2 = v1;
645  // Index i = 10000000;
646 
647  v1 /= (Numeric)i;
648  v2 /= (Numeric)i;
649  cout.precision(12);
650  // cout.setf(ios_base::scientific,ios_base::floatfield);
651  v1 *= v1;
652  v2 *= v2;
653  cout << v1 << endl;
654  cout << v2 << endl;
655 }
656 
657 void test38() {
658  Vector v(5, 0.);
659  Numeric* const a = v.get_c_array();
660 
661  a[4] = 5.;
662 
663  cout << v << endl;
664  cout << endl << "========================" << endl << endl;
665 
666  Matrix m(5, 5, 0.);
667  Numeric* const b = m.get_c_array();
668 
669  b[4] = 5.;
670 
671  cout << m << endl;
672  cout << endl << "========================" << endl << endl;
673 
674  Tensor3 t3(5, 6, 7, 0.);
675  Numeric* const c = t3.get_c_array();
676 
677  c[6] = 5.;
678 
679  cout << t3 << endl;
680 }
681 
682 void test39() {
683  Vector v1(1, 5, 1), v2(5);
684 
685  v2 = v1 * 2;
686  // Unfortunately, this thing compiles, but at least it gives an
687  // assertion failure at runtime. I think what happens is that it
688  // implicitly creates a one element vector out of the "2", then
689  // tries to do a scalar product with v1.
690 
691  cout << v2 << endl;
692 }
693 
694 void test40() {
695  Vector v(100);
696 
697  v = 5;
698 
699  cout << v << endl;
700 }
701 
702 void test41() {
703  const Vector v1(10, 1);
704 
705  ConstVectorView vv = v1[Range(0, 5)];
706 
707  cout << "Vector: " << v1 << endl;
708  cout << "VectorView: " << vv << endl;
709 
710  try {
711  vv = 2;
712  } catch (const std::runtime_error& e) {
713  std::cerr << e.what() << endl;
714  exit(EXIT_FAILURE);
715  }
716 
717  cout << "VectorView: " << vv << endl;
718  cout << "Vector: " << v1 << endl;
719 }
720 
721 // Test behaviour of VectorView::operator MatrixView, for which I have fixed
722 // a bug today.
723 // SAB 2013-01-18
724 void test42() {
725  cout << "test42\n";
726  Vector x{1, 2, 3, 4, 5, 6, 7, 8, 9, 10};
727  cout << "x: " << x << endl;
728 
729  VectorView y = x[Range(2, 4, 2)];
730  cout << "y: " << y << endl;
731 
733  cout << "z:\n" << z << endl;
734 
735  cout << "Every other z:\n" << z(Range(1, 2, 2), joker) << endl;
736 
737  // Try write access also:
738  MatrixView zz(y);
739  zz(Range(1, 2, 2), joker) = 0;
740  cout << "zz:\n" << zz << endl;
741  cout << "New x: " << x << endl;
742 }
743 
744 
745 // Test function for internal use, must return 2 as last n2r
746 constexpr Rational test_numeric2rational(const Index i, const Index maxi, const Rational r=0_rat, const Rational n2r=0_rat) {
747  if (i > maxi)
748  return n2r;
749  else {
750  return
751  (r == n2r) ?
752  test_numeric2rational(i+1, maxi, Rational(maxi + i, maxi), numeric2rational(1 + Numeric(i)/Numeric(maxi), 12)) :
753  throw std::logic_error("Fail to initialize");;
754  }
755 }
756 
757 
758 void test43() {
759  // Simple construction compile-time test
760  constexpr Rational r=3_2; // should be 3/2
761  static_assert(r.Nom() == 3, "Rational fail to initialize properly");
762  static_assert(r.Denom() == 2, "Rational fail to initialize properly");
763 
764  // Simple expression compile-time test
765  constexpr Rational r2 = 1 + r; // should be 5/2
766  static_assert(r2.Nom() == 5, "Rational fail to initialize properly");
767  static_assert(r2.Denom() == 2, "Rational fail to initialize properly");
768 
769  // Div-zero by index generates an undefined rational: 0/0
770  constexpr Rational r3 = r / 0; // should be undefined
771  static_assert(r3.Nom() == 0, "Rational fail to initialize properly");
772  static_assert(r3.Denom() == 0, "Rational fail to initialize properly");
773  static_assert(r3 not_eq r3, "Rational fail to initialize properly");
774  static_assert(r3.isUndefined(), "Rational fail to initialize properly");
775 
776  // Mul-zero gives 0/1
777  constexpr Rational r4 = r * 0; // should be undefined
778  static_assert(r4.Nom() == 0, "Rational fail to initialize properly");
779  static_assert(r4.Denom() == 1, "Rational fail to initialize properly");
780 
781  // Complicated expression compile-time test
782  constexpr Rational r5 = (Rational(1, 9) * r) % 3; // should be 1/6
783  static_assert(r5.Nom() == 1, "Rational fail to initialize properly");
784  static_assert(r5.Denom() == 6, "Rational fail to initialize properly");
785 
786  // The simplify operation simplifies expressions
787  constexpr Rational r6 = 10 % r; // should be 1
788  static_assert(r6.Nom() == 1, "Rational fail to initialize properly");
789  static_assert(r6.Denom() == 1, "Rational fail to initialize properly");
790  static_assert(r6.toInt() == 1, "Rational fail to initialize properly");
791  static_assert(r6.toIndex() == 1, "Rational fail to initialize properly");
792  static_assert(r6.toNumeric() == 1e0, "Rational fail to initialize properly");
793 
794  constexpr Index rattest=1<<8;
795  constexpr Rational r7 = test_numeric2rational(0, rattest);
796  static_assert(r7 == 2, "Rational fail to initialize properly");
797 }
798 
799 void test44() {
800 #define docheck(fn, val, expect) \
801  cout << #fn << "(" << val << ") = " << fn(x) << " (expected " << #expect \
802  << ")" << std::endl;
803 
804  Vector x{1, 2, 3};
806  docheck(is_sorted, x, 1)
807 
808  x = {3, 2, 1};
810  docheck(is_sorted, x, 0)
811 
812  x = {1, 2, 2};
814  docheck(is_sorted, x, 1)
815 
816  x = {2, 2, 1};
818  docheck(is_sorted, x, 0)
819 
820  x = {1, NAN, 2};
822  docheck(is_sorted, x, 0)
823 
824  x = {2, NAN, 1};
826  docheck(is_sorted, x, 0)
827 
828  x = {NAN, NAN, NAN};
830  docheck(is_sorted, x, 0)
831 
832 #undef docheck
833 }
834 
835 void test46() {
836  Vector v(5, 0.);
837  nlinspace(v, 1, 10, 10);
838  VectorView v1 = v;
839  auto compare_func = [](Numeric n) { return n != 0; };
840  cout << std::any_of(v1.begin(), v1.end(), compare_func) << endl;
841  v1 = 0.;
842  cout << std::any_of(v1.begin(), v1.end(), compare_func) << endl;
843 }
844 
846 
854 bool test_diagonal(Index ntests) {
855  Rand<Index> rand(1, 100);
856 
857  Index n_diag;
858 
859  bool pass = true;
860 
861  for (Index i = 0; i < ntests; i++) {
862  Matrix A(rand(), rand());
863  n_diag = std::min(A.ncols(), A.nrows());
864 
865  ConstMatrixView AT = transpose(A);
866  ConstMatrixView B = A;
867  if (n_diag > 3) {
868  B = A(Range(1, Joker(), 2), Range(1, Joker(), 2));
869  }
870 
871  ConstVectorView v(A.diagonal());
872  ConstVectorView vt(AT.diagonal());
873  ConstVectorView vb(B.diagonal());
874 
875  random_fill_matrix(A, 10, true);
876 
877  for (Index j = 0; j < n_diag; j++) {
878  pass = pass && (v[j] == A(j, j));
879  pass = pass && (vt[j] == AT(j, j));
880 
881  if (j < vb.nelem()) pass = pass && (vb[j] == B(j, j));
882  }
883  cout << endl;
884  }
885 
886  if (pass)
887  cout << "test diagonal: Passed all tests." << endl;
888  else
889  cout << "test diagonal: Failed." << endl;
890 
891  return pass;
892 }
893 
894 Numeric matrix_vector_mult(Index m, Index n, Index ntests, bool verbose) {
895  Numeric max_err = 0;
896  Matrix A(m, n);
897  Vector x(n), xt, x_ref(n), y(m), y_ref(m);
898 
899  Rand<Index> random_row_stride(1, n / 4);
900  Rand<Index> random_column_stride(1, m / 4);
901 
902  for (Index i = 0; i < ntests; i++) {
903  // A * x
904 
905  random_fill_matrix(A, 1000, false);
906  random_fill_vector(x, 1000, false);
907 
908  mult(y, A, x);
909  mult_general(y_ref, A, x);
910 
911  Numeric err_mul = get_maximum_error(y, y_ref, true);
912  if (err_mul > max_err) max_err = err_mul;
913 
914  if (verbose) {
915  cout << "\t A * x: max. rel. error = " << err_mul << endl;
916  }
917 
918  // A^T * y
919 
920  mult(x, transpose(A), y);
921  mult_general(x_ref, transpose(A), y);
922 
923  err_mul = get_maximum_error(x, x_ref, true);
924  if (err_mul > max_err) max_err = err_mul;
925 
926  if (verbose) {
927  cout << "\t A^T * x: max. rel. error = " << err_mul << endl;
928  }
929 
930  // Random stride
931  Index column_stride(random_column_stride()),
932  row_stride(random_row_stride());
933 
934  Index m_sub, n_sub;
935  m_sub = (m - 1) / row_stride + 1;
936  n_sub = (n - 1) / column_stride + 1;
937 
938  mult(y[Range(0, joker, row_stride)],
939  A(Range(0, m_sub), Range(0, n_sub)),
940  x[Range(0, joker, column_stride)]);
941  mult_general(y_ref[Range(0, joker, row_stride)],
942  A(Range(0, m_sub), Range(0, n_sub)),
943  x[Range(0, joker, column_stride)]);
944 
945  err_mul = get_maximum_error(y[Range(0, joker, row_stride)],
946  y_ref[Range(0, joker, row_stride)],
947  true);
948  if (err_mul > max_err) max_err = err_mul;
949 
950  if (verbose) {
951  cout << "\t Random stride: max. rel. error = " << err_mul << endl << endl;
952  }
953 
954  // Random offset
955  if ((m > 1) && (n > 1)) {
956  Index y_offset = rand() % (m - 1);
957  Index x_offset = rand() % (n - 1);
958 
959  m_sub = m - y_offset - 1;
960  n_sub = n - x_offset - 1;
961 
962  mult(y[Range(y_offset, m_sub)],
963  A(Range(y_offset, m_sub), Range(x_offset, n_sub)),
964  x[Range(x_offset, n_sub)]);
965  mult_general(y_ref[Range(y_offset, m_sub)],
966  A(Range(y_offset, m_sub), Range(x_offset, n_sub)),
967  x[Range(x_offset, n_sub)]);
968 
969  err_mul = get_maximum_error(
970  y[Range(y_offset, m_sub)], y_ref[Range(y_offset, m_sub)], true);
971  if (err_mul > max_err) max_err = err_mul;
972 
973  if (verbose) {
974  cout << "\t Random offset: max. rel. error = " << err_mul << endl
975  << endl;
976  }
977  }
978  }
979 
980  return max_err;
981 }
982 
984  Numeric err, max_err;
985 
986  if (verbose)
987  cout << "Matrix-Vector Multiplication: n = m = 100, ntests = 100" << endl;
988 
989  max_err = matrix_vector_mult(100, 100, 100, verbose);
990  if (verbose) {
991  cout << endl;
992  cout << "Matrix-Vector Multiplication: n = 100, m = 20, ntests = 100"
993  << endl;
994  }
995 
996  err = matrix_vector_mult(100, 20, 100, verbose);
997  if (err > max_err) max_err = err;
998  if (verbose) {
999  cout << endl;
1000  cout << "Matrix-Vector Multiplication: n = 20, m = 100, ntests = 100"
1001  << endl;
1002  }
1003 
1004  err = matrix_vector_mult(20, 100, 100, verbose);
1005  if (err > max_err) max_err = err;
1006  if (verbose) {
1007  if (max_err < 1e-9)
1008  cout << endl << "Matrix Vector Multiplication: PASSED" << endl;
1009  else
1010  cout << endl << "Matrix Vector Multiplication: FAILED" << endl;
1011  }
1012 
1013  err = matrix_vector_mult(100, 1, 100, verbose);
1014  if (err > max_err) max_err = err;
1015  if (verbose) {
1016  cout << endl;
1017  cout << "Matrix-Vector Multiplication: n = 20, m = 100, ntests = 100"
1018  << endl;
1019  }
1020 
1021  err = matrix_vector_mult(1, 100, 100, verbose);
1022  if (err > max_err) max_err = err;
1023  if (verbose) {
1024  cout << endl;
1025  cout << "Matrix-Vector Multiplication: n = 20, m = 100, ntests = 100"
1026  << endl;
1027  }
1028 
1029  if (verbose) {
1030  if (max_err < 1e-9)
1031  cout << endl << "Matrix Vector Multiplication: PASSED" << endl;
1032  else
1033  cout << endl << "Matrix Vector Multiplication: FAILED" << endl;
1034  }
1035  return max_err;
1036 }
1037 
1039 
1059  Index k, Index m, Index n, Index ntests, Index nsubtests, bool verbose) {
1060  Numeric max_err = 0;
1061  Matrix A(m, k), B(k, n), C(m, n), C_ref(m, n);
1062 
1063  for (Index i = 0; i < ntests; i++) {
1064  random_fill_matrix(A, 1000, false);
1065  random_fill_matrix(B, 1000, false);
1066 
1067  if (verbose) {
1068  cout.precision(15);
1069  cout << "MATRIX MULTIPLICATION: m = " << m << ", k = " << k;
1070  cout << ", n = " << n << endl;
1071  }
1072 
1073  // A * B
1074 
1075  random_fill_matrix(C, 10, false);
1076  random_fill_matrix(C_ref, 10, false);
1077 
1078  mult(C, A, B);
1079  mult_general(C_ref, A, B);
1080 
1081  Numeric err_mul = get_maximum_error(C, C_ref, true);
1082  if (err_mul > max_err) max_err = err_mul;
1083 
1084  if (verbose) {
1085  cout << "\t"
1086  << "A * B: max. rel. error = " << err_mul << endl;
1087  }
1088 
1089  // A^T * B
1090 
1091  Matrix AT(k, m);
1092  random_fill_matrix(AT, 100.0, false);
1093 
1094  mult(C, transpose(AT), B);
1095  mult_general(C_ref, transpose(AT), B);
1096  Numeric err_trans1 = get_maximum_error(C, C_ref, true);
1097  if (err_trans1 > max_err) max_err = err_trans1;
1098 
1099  if (verbose) {
1100  cout << "\t"
1101  << "A^T * B: max. rel. err = " << err_trans1 << endl;
1102  }
1103 
1104  // A * B^T
1105 
1106  Matrix BT(n, k);
1107  random_fill_matrix(BT, 100.0, false);
1108 
1109  mult(C, A, transpose(BT));
1110  mult_general(C_ref, A, transpose(BT));
1111  Numeric err_trans2 = get_maximum_error(C, C_ref, true);
1112  if (err_trans2 > max_err) max_err = err_trans2;
1113 
1114  if (verbose) {
1115  cout << "\t"
1116  << "A * B^T: max. rel. err = " << err_trans2 << endl;
1117  }
1118 
1119  // A^T * B^T
1120 
1121  mult(C, transpose(AT), transpose(BT));
1122  mult_general(C_ref, transpose(AT), transpose(BT));
1123  Numeric err_trans3 = get_maximum_error(C, C_ref, true);
1124  if (err_trans3 > max_err) max_err = err_trans3;
1125 
1126  if (verbose) {
1127  cout << "\t"
1128  << "A^T * B^T: max. rel. err = " << err_trans3 << endl;
1129  cout << endl;
1130  }
1131  }
1132 
1133  // Multiplication of submatrices.
1134  Rand<Index> k_rand(1, k - 1), m_rand(1, m - 1), n_rand(1, n - 1);
1135 
1136  for (Index i = 0; i < nsubtests - 1; i++) {
1137  Index k1(k_rand()), m1(m_rand()), n1(n_rand());
1138 
1139  Range r1(random_range(m1));
1140  Range r2(random_range(k1));
1141  Range r3(random_range(n1));
1142 
1143  MatrixView C_sub(C(r1, r3));
1144  MatrixView C_sub_ref(C_ref(r1, r3));
1145 
1146  ConstMatrixView A_sub(A(r1, r2));
1147  ConstMatrixView B_sub(B(r2, r3));
1148 
1149  mult(C_sub, A_sub, B_sub);
1150  mult_general(C_sub_ref, A_sub, B_sub);
1151 
1152  Numeric err = get_maximum_error(C_sub, C_sub_ref, true);
1153  if (err > max_err) max_err = err;
1154 
1155  if (verbose) {
1156  cout << "\t"
1157  << "Submatrix multiplication: max. rel. err = " << err;
1158  cout << endl;
1159  }
1160  }
1161 
1162  if (verbose) {
1163  cout << endl;
1164  }
1165 
1166  return max_err;
1167 }
1168 
1170 
1182  Numeric max_err, err;
1183 
1184  // Trivial case k = m = n = 0.
1185  max_err = matrix_mult(0, 0, 0, 10, 0, verbose);
1186 
1187  // k = 1, m = 1, n = 1.
1188  err = matrix_mult(1, 1, 1, 10, 0, verbose);
1189  if (err > max_err) max_err = err;
1190 
1191  // k = 10, m = 1, n = 10.
1192  err = matrix_mult(10, 1, 10, 20, 20, verbose);
1193  if (err > max_err) max_err = err;
1194 
1195  // k = 10, m = 1, n = 10.
1196  err = matrix_mult(10, 1, 10, 20, 20, verbose);
1197  if (err > max_err) max_err = err;
1198 
1199  // k = 10, m = 1, n = 1.
1200  err = matrix_mult(10, 1, 1, 20, 20, verbose);
1201  if (err > max_err) max_err = err;
1202 
1203  // k = 100, m = 100, n = 100.
1204  err = matrix_mult(200, 200, 200, 20, 20, verbose);
1205  if (err > max_err) max_err = err;
1206 
1207  // k = 10, m = 100, n = 10.
1208  err = matrix_mult(10, 100, 10, 20, 20, verbose);
1209  if (err > max_err) max_err = err;
1210 
1211  // k = 10, m = 100, n = 100.
1212  err = matrix_mult(10, 100, 100, 20, 20, verbose);
1213  if (err > max_err) max_err = err;
1214 
1215  // k = 100, m = 100, n = 100, 100 submatrix multiplications.
1216  err = matrix_mult(100, 100, 100, 0, 100, verbose);
1217  if (err > max_err) max_err = err;
1218 
1219  return max_err;
1220 }
1221 
1223 void test_empty() {
1224  {
1225  cout << "Array" << endl;
1226  ArrayOfIndex v;
1227  v.resize(0);
1228  cout << v.empty() << endl;
1229  v.resize(1);
1230  cout << v.empty() << endl;
1231  }
1232  {
1233  cout << "Vector" << endl;
1234  Vector v;
1235  v.resize(0);
1236  cout << v.empty() << endl;
1237  v.resize(1);
1238  cout << v.empty() << endl;
1239  }
1240  {
1241  cout << "Matrix" << endl;
1242  Matrix v;
1243  v.resize(0, 1);
1244  cout << v.empty() << endl;
1245  v.resize(1, 0);
1246  cout << v.empty() << endl;
1247  v.resize(1, 1);
1248  cout << v.empty() << endl;
1249  }
1250  {
1251  cout << "Tensor3" << endl;
1252  Tensor3 v;
1253  v.resize(0, 1, 1);
1254  cout << v.empty() << endl;
1255  v.resize(1, 0, 1);
1256  cout << v.empty() << endl;
1257  v.resize(1, 1, 0);
1258  cout << v.empty() << endl;
1259  v.resize(1, 1, 1);
1260  cout << v.empty() << endl;
1261  }
1262  {
1263  cout << "Tensor4" << endl;
1264  Tensor4 v;
1265  v.resize(0, 1, 1, 1);
1266  cout << v.empty() << endl;
1267  v.resize(1, 0, 1, 1);
1268  cout << v.empty() << endl;
1269  v.resize(1, 1, 0, 1);
1270  cout << v.empty() << endl;
1271  v.resize(1, 1, 1, 0);
1272  cout << v.empty() << endl;
1273  v.resize(1, 1, 1, 1);
1274  cout << v.empty() << endl;
1275  }
1276  {
1277  cout << "Tensor5" << endl;
1278  Tensor5 v;
1279  v.resize(0, 1, 1, 1, 1);
1280  cout << v.empty() << endl;
1281  v.resize(1, 0, 1, 1, 1);
1282  cout << v.empty() << endl;
1283  v.resize(1, 1, 0, 1, 1);
1284  cout << v.empty() << endl;
1285  v.resize(1, 1, 1, 0, 1);
1286  cout << v.empty() << endl;
1287  v.resize(1, 1, 1, 1, 0);
1288  cout << v.empty() << endl;
1289  v.resize(1, 1, 1, 1, 1);
1290  cout << v.empty() << endl;
1291  }
1292  {
1293  cout << "Tensor6" << endl;
1294  Tensor6 v;
1295  v.resize(0, 1, 1, 1, 1, 1);
1296  cout << v.empty() << endl;
1297  v.resize(1, 0, 1, 1, 1, 1);
1298  cout << v.empty() << endl;
1299  v.resize(1, 1, 0, 1, 1, 1);
1300  cout << v.empty() << endl;
1301  v.resize(1, 1, 1, 0, 1, 1);
1302  cout << v.empty() << endl;
1303  v.resize(1, 1, 1, 1, 0, 1);
1304  cout << v.empty() << endl;
1305  v.resize(1, 1, 1, 1, 1, 0);
1306  cout << v.empty() << endl;
1307  v.resize(1, 1, 1, 1, 1, 1);
1308  cout << v.empty() << endl;
1309  }
1310  {
1311  cout << "Tensor7" << endl;
1312  Tensor7 v;
1313  v.resize(0, 1, 1, 1, 1, 1, 1);
1314  cout << v.empty() << endl;
1315  v.resize(1, 0, 1, 1, 1, 1, 1);
1316  cout << v.empty() << endl;
1317  v.resize(1, 1, 0, 1, 1, 1, 1);
1318  cout << v.empty() << endl;
1319  v.resize(1, 1, 1, 0, 1, 1, 1);
1320  cout << v.empty() << endl;
1321  v.resize(1, 1, 1, 1, 0, 1, 1);
1322  cout << v.empty() << endl;
1323  v.resize(1, 1, 1, 1, 1, 0, 1);
1324  cout << v.empty() << endl;
1325  v.resize(1, 1, 1, 1, 1, 1, 0);
1326  cout << v.empty() << endl;
1327  v.resize(1, 1, 1, 1, 1, 1, 1);
1328  cout << v.empty() << endl;
1329  }
1330 }
1331 
1333  const Numeric start,
1334  const Numeric stop,
1335  const Index n) {
1336  assert(1 < n); // Number of points must be greater 1.
1337  x.resize(n);
1338  Numeric step = (stop - start) / ((double)n - 1);
1339  for (Index i = 0; i < n - 1; i++) x[i] = start + (double)i * step;
1340  x[n - 1] = stop;
1341 }
1342 
1343 void test47() {
1344  // Selecting empty matpack dimensions with Joker shouldn't fail
1345  Vector v;
1346  Matrix m;
1347 
1348  VectorView vv = v[joker];
1349  MatrixView mv = m(joker, joker);
1350 
1351  std::cout << "vv.nelem: " << vv.nelem() << std::endl;
1352  std::cout << "mv.nrows: " << mv.nrows() << " ";
1353  std::cout << "mv.ncols: " << mv.ncols() << std::endl;
1354 
1355  Matrix m2(0, 5);
1356  MatrixView mv2 = m2(joker, 3);
1357  std::cout << "mv2.nrows: " << mv2.nrows() << " ";
1358  std::cout << "mv2.ncols: " << mv2.ncols() << std::endl;
1359 }
1360 
1361 int main() {
1362  // test1();
1363  // test2();
1364  // test3();
1365  // test4();
1366  // test5();
1367  // test6();
1368  // test7();
1369  // test8();
1370  // test9();
1371  // test10();
1372  // test11();
1373  // test12();
1374  // test13();
1375  // test14();
1376  // test15();
1377  // test16();
1378  // test17();
1379  // test18();
1380  // test19();
1381  // test20();
1382  // test21();
1383  // test22();
1384  // test23();
1385  // test24();
1386  // test25();
1387  // test26();
1388  // test27();
1389  // test28();
1390  // test29();
1391  // test30();
1392  // test31();
1393  // test32();
1394  // test33();
1395  // test34();
1396  // test35();
1397  // test36();
1398  // Index i = 10000000;
1399  // test37(i);
1400  // test38();
1401  // test39();
1402  // test40();
1403  // test41();
1404  // test42();
1405  test43();
1406  // test44();
1407  // test45();
1408  // test46();
1409  // test47();
1410 
1411  // const double tolerance = 1e-9;
1412  // double error;
1413  //
1414  // // Matrix Vector Multiplication.
1415  // error = test_matrix_vector_multiplication(false);
1416  // cout << "Matrix Vector Multiplication: ";
1417  // if (error > tolerance)
1418  // cout << "FAILED, maximum error: " << error << endl;
1419  // else
1420  // cout << "PASSED." << endl;
1421  //
1422  // // Matrix Matrix Multiplication.
1423  // error = test_matrix_multiplication(false);
1424  //
1425  // cout << "Matrix Matrix Multiplication: ";
1426  // if (error > tolerance)
1427  // cout << "FAILED, maximum error: " << error << endl;
1428  // else
1429  // cout << "PASSED." << endl;
1430 
1431  //test_diagonal( 100 );
1432  //test_empty();
1433 
1434  return 1;
1435 }
Matrix
The Matrix class.
Definition: matpackI.h:1193
Tensor7View
The Tensor7View class.
Definition: matpackVII.h:1286
transform
void transform(VectorView y, double(&my_func)(double), ConstVectorView x)
A generic transform function for vectors, which can be used to implement mathematical functions opera...
Definition: matpackI.cc:1476
main
int main()
Definition: test_matpack.cc:1361
Tensor4::resize
void resize(Index b, Index p, Index r, Index c)
Resize function.
Definition: matpackIV.cc:1064
test35
void test35()
Definition: test_matpack.cc:615
MatrixView
The MatrixView class.
Definition: matpackI.h:1093
numeric2rational
constexpr Rational numeric2rational(Numeric x, size_t maxdec=4)
Rational from Numeric.
Definition: rational.h:330
exceptions.h
The declarations of all the exception classes.
Tensor6::resize
void resize(Index v, Index s, Index b, Index p, Index r, Index c)
Resize function.
Definition: matpackVI.cc:2175
test13
void test13()
Definition: test_matpack.cc:266
test43
void test43()
Definition: test_matpack.cc:758
rational.h
Contains the rational class definition.
test20
void test20()
Definition: test_matpack.cc:332
test32
void test32()
Definition: test_matpack.cc:518
test_utils.h
Utility functions for testing.
Tensor3
The Tensor3 class.
Definition: matpackIII.h:339
test18
void test18()
Definition: test_matpack.cc:314
wigner_functions.h
Wigner symbol interactions.
test19
void test19()
Definition: test_matpack.cc:321
test38
void test38()
Definition: test_matpack.cc:657
joker
const Joker joker
docheck
#define docheck(fn, val, expect)
describe.h
Header file for describe.cc.
matrix_vector_mult
Numeric matrix_vector_mult(Index m, Index n, Index ntests, bool verbose)
Definition: test_matpack.cc:894
test46
void test46()
Definition: test_matpack.cc:835
mult
void mult(ComplexVectorView y, const ConstComplexMatrixView &M, const ConstComplexVectorView &x)
Matrix-Vector Multiplication.
Definition: complex.cc:1579
test30
void test30()
Definition: test_matpack.cc:433
random_fill_vector
void random_fill_vector(VectorView v, Numeric range, bool positive)
Fill vector with random values.
Definition: test_utils.cc:230
ARTS::Var::y
Vector y(Workspace &ws) noexcept
Definition: autoarts.h:7401
Vector::resize
void resize(Index n)
Resize function.
Definition: matpackI.cc:404
by_reference
Numeric by_reference(const Numeric &x)
Definition: test_matpack.cc:39
test33
void test33()
Definition: test_matpack.cc:540
ConstMatrixView::nrows
Index nrows() const
Returns the number of rows.
Definition: matpackI.cc:429
test25
void test25()
Definition: test_matpack.cc:382
test_matrix_vector_multiplication
Numeric test_matrix_vector_multiplication(bool verbose)
Definition: test_matpack.cc:983
Tensor4
The Tensor4 class.
Definition: matpackIV.h:421
array.h
This file contains the definition of Array.
Tensor3::resize
void resize(Index p, Index r, Index c)
Resize function.
Definition: matpackIII.cc:664
sqrt
Numeric sqrt(const Rational r)
Square root.
Definition: rational.h:620
ConstTensor3View::npages
Index npages() const
Returns the number of pages.
Definition: matpackIII.h:144
test6
void test6()
Definition: test_matpack.cc:186
Tensor3View
The Tensor3View class.
Definition: matpackIII.h:239
Array< String >
test40
void test40()
Definition: test_matpack.cc:694
is_decreasing
bool is_decreasing(ConstVectorView x)
Checks if a vector is sorted in reversed order and is strictly decreasing.
Definition: logic.cc:252
Rand< Index >
Definition: test_utils.h:61
min
#define min(a, b)
Definition: legacy_continua.cc:20628
random_fill_matrix
void random_fill_matrix(MatrixView A, Numeric range, bool positive)
Fill matrix with random values.
Definition: test_utils.cc:59
Tensor7::resize
void resize(Index l, Index v, Index s, Index b, Index p, Index r, Index c)
Resize function.
Definition: matpackVII.cc:5484
ConstTensor5View::empty
bool empty() const
Check if variable is empty.
Definition: matpackV.cc:38
Tensor5::resize
void resize(Index s, Index b, Index p, Index r, Index c)
Resize function.
Definition: matpackV.cc:1743
test2
void test2()
Definition: test_matpack.cc:129
describe
String describe(ConstTensor7View x)
Describe Tensor7.
Definition: describe.cc:38
test42
void test42()
Definition: test_matpack.cc:724
ConstMatrixView::ncols
Index ncols() const
Returns the number of columns.
Definition: matpackI.cc:432
test26
void test26()
Definition: test_matpack.cc:388
my_basic_string< char >
test21
void test21()
Definition: test_matpack.cc:338
ConstTensor4View::empty
bool empty() const
Check if variable is empty.
Definition: matpackIV.cc:52
test41
void test41()
Definition: test_matpack.cc:702
test1
int test1()
Definition: test_matpack.cc:47
test22
void test22()
Definition: test_matpack.cc:349
test34
void test34()
Definition: test_matpack.cc:605
ConstVectorView::nelem
Index nelem() const
Returns the number of elements.
Definition: matpackI.cc:51
VectorView
The VectorView class.
Definition: matpackI.h:610
test_empty
void test_empty()
Check if the empty function is working correctly.
Definition: test_matpack.cc:1223
ARTS::Group::Rational
Rational Rational
Definition: autoarts.h:96
test_numeric2rational
constexpr Rational test_numeric2rational(const Index i, const Index maxi, const Rational r=0_rat, const Rational n2r=0_rat)
Definition: test_matpack.cc:746
ConstTensor7View::empty
bool empty() const
Check if variable is empty.
Definition: matpackVII.cc:36
ConstVectorView::sum
Numeric sum() const
The sum of all elements of a Vector.
Definition: matpackI.cc:53
test10
void test10()
Definition: test_matpack.cc:228
Numeric
NUMERIC Numeric
The type to use for all floating point numbers.
Definition: matpack.h:33
test17
void test17()
Definition: test_matpack.cc:308
Tensor5View
The Tensor5View class.
Definition: matpackV.h:333
junk2
void junk2(ConstVectorView a)
Definition: test_matpack.cc:603
math_funcs.h
max
#define max(a, b)
Definition: legacy_continua.cc:20629
Matrix::resize
void resize(Index r, Index c)
Resize function.
Definition: matpackI.cc:1056
test47
void test47()
Definition: test_matpack.cc:1343
ConstMatrixView::diagonal
ConstVectorView diagonal() const
Matrix diagonal as vector.
Definition: matpackI.cc:489
ConstTensor3View::nrows
Index nrows() const
Returns the number of rows.
Definition: matpackIII.h:147
test14
void test14()
Definition: test_matpack.cc:276
ConstMatrixView::empty
bool empty() const
Returns true if variable size is zero.
Definition: matpackI.cc:426
ConstMatrixView
A constant view of a Matrix.
Definition: matpackI.h:982
test27
void test27()
Definition: test_matpack.cc:394
mult_general
void mult_general(VectorView y, const ConstMatrixView &M, const ConstVectorView &x)
Matrix Vector multiplication.
Definition: matpackI.cc:1181
Tensor5
The Tensor5 class.
Definition: matpackV.h:506
test29
void test29()
Definition: test_matpack.cc:409
oem::Vector
invlib::Vector< ArtsVector > Vector
invlib wrapper type for ARTS vectors.
Definition: oem.h:32
test36
void test36()
Definition: test_matpack.cc:634
Range
The range class.
Definition: matpackI.h:160
junk4
void junk4(Tensor4View a)
Definition: test_matpack.cc:601
nlinspace
void nlinspace(Vector &x, const Numeric start, const Numeric stop, const Index n)
nlinspace
Definition: test_matpack.cc:1332
test39
void test39()
Definition: test_matpack.cc:682
test_diagonal
bool test_diagonal(Index ntests)
Test diagonal vector.
Definition: test_matpack.cc:854
is_increasing
bool is_increasing(ConstVectorView x)
Checks if a vector is sorted and strictly increasing.
Definition: logic.cc:215
logic.h
Header file for logic.cc.
test28
void test28()
Definition: test_matpack.cc:402
Zeeman::start
constexpr Rational start(Rational Ju, Rational Jl, Polarization type) noexcept
Gives the lowest M for a polarization type of this transition.
Definition: zeemandata.h:77
Joker
The Joker class.
Definition: matpackI.h:126
test44
void test44()
Definition: test_matpack.cc:799
ConstTensor3View::ncols
Index ncols() const
Returns the number of columns.
Definition: matpackIII.h:150
test7
void test7()
Definition: test_matpack.cc:207
VectorView::end
Iterator1D end()
Return iterator behind last element.
Definition: matpackI.cc:148
test24
void test24()
Definition: test_matpack.cc:368
ConstVectorView::empty
bool empty() const
Returns true if variable size is zero.
Definition: matpackI.cc:49
random_range
Range random_range(Index n)
Generate random sub-range of the range [0, n-1].
Definition: test_utils.cc:273
test9
void test9()
Definition: test_matpack.cc:221
Tensor4View
The Tensor4View class.
Definition: matpackIV.h:284
matpackII.h
Header file for sparse matrices.
matrix_mult
Numeric matrix_mult(Index k, Index m, Index n, Index ntests, Index nsubtests, bool verbose)
Perform matrix multiplication test.
Definition: test_matpack.cc:1058
get_maximum_error
Numeric get_maximum_error(ConstVectorView v1, ConstVectorView v2, bool relative)
Maximum element-wise error of two vectors.
Definition: test_utils.cc:297
test_matrix_multiplication
Numeric test_matrix_multiplication(bool verbose)
Perform matrix multiplication tests.
Definition: test_matpack.cc:1181
VectorView::begin
Iterator1D begin()
Return iterator to first element.
Definition: matpackI.cc:144
test5
void test5()
Definition: test_matpack.cc:167
ARTS::Var::x
Vector x(Workspace &ws) noexcept
Definition: autoarts.h:7346
fill_with_junk
void fill_with_junk(VectorView x)
Definition: test_matpack.cc:43
Index
INDEX Index
The type to use for all integer numbers and indices.
Definition: matpack.h:39
transpose
ConstComplexMatrixView transpose(ConstComplexMatrixView m)
Const version of transpose.
Definition: complex.cc:1509
ConstTensor3View::empty
bool empty() const
Check if variable is empty.
Definition: matpackIII.cc:38
Tensor6
The Tensor6 class.
Definition: matpackVI.h:1088
MatrixView::get_c_array
const Numeric * get_c_array() const
Conversion to plain C-array.
Definition: matpackI.cc:735
M
#define M
Definition: rng.cc:165
test31
void test31()
Definition: test_matpack.cc:454
test4
void test4()
Definition: test_matpack.cc:141
test37
void test37(const Index &i)
Definition: test_matpack.cc:642
test12
void test12()
Definition: test_matpack.cc:249
by_value
Numeric by_value(Numeric x)
Definition: test_matpack.cc:41
Tensor3View::get_c_array
const Numeric * get_c_array() const
Conversion to plain C-array.
Definition: matpackIII.cc:321
Vector
The Vector class.
Definition: matpackI.h:860
is_sorted
bool is_sorted(ConstVectorView x)
Checks if a vector is sorted in ascending order.
Definition: logic.cc:199
test11
void test11()
Definition: test_matpack.cc:238
ConstVectorView
A constant view of a Vector.
Definition: matpackI.h:476
test15
void test15()
Definition: test_matpack.cc:286
Array::nelem
Index nelem() const
Number of elements.
Definition: array.h:195
Rational
Implements rational numbers to work with other ARTS types.
Definition: rational.h:54
VectorView::get_c_array
const Numeric * get_c_array() const
Conversion to plain C-array, const-version.
Definition: matpackI.cc:272
mystring.h
This file contains the definition of String, the ARTS string class.
Tensor7
The Tensor7 class.
Definition: matpackVII.h:2382
test8
void test8()
Definition: test_matpack.cc:214
matpackVII.h
test23
void test23()
Definition: test_matpack.cc:360
Tensor6View
The Tensor6View class.
Definition: matpackVI.h:621
ConstTensor6View::empty
bool empty() const
Check if variable is empty.
Definition: matpackVI.cc:36