ARTS 2.5.4 (git: 31ce4f0e)
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 <exception>
22#include <iostream>
23#include "array.h"
24#include "describe.h"
25#include "exceptions.h"
26#include "logic.h"
27#include "math_funcs.h"
28#include "matpackII.h"
29#include "matpackVII.h"
30#include "matpack_eigen.h"
31#include "mystring.h"
32#include "rational.h"
33#include "test_utils.h"
34#include "time.h"
35#include "wigner_functions.h"
36
37using std::cout;
38using std::endl;
39using std::runtime_error;
40
41Numeric by_reference(const Numeric& x) { return x + 1; }
42
43Numeric by_value(Numeric x) { return x + 1; }
44
45void fill_with_junk(VectorView x) { x = 999; }
46
47void fill_with_junk(MatrixView x) { x = 888; }
48
49int test1() {
50 Vector v(20);
51
52 cout << "v.nelem() = " << v.nelem() << "\n";
53
54 for (Index i = 0; i < v.nelem(); ++i) v[i] = (Numeric)i;
55
56 cout << "v.begin() = " << *v.begin() << "\n";
57
58 cout << "v = \n" << v << "\n";
59
60 fill_with_junk(v[Range(1, 8, 2)][Range(2, joker)]);
61 // fill_with_junk(v);
62
63 Vector v2 = v[Range(2, 4)];
64
65 cout << "v2 = \n" << v2 << "\n";
66
67 for (Index i = 0; i < 1000; ++i) {
68 Vector v3(1000);
69 v3 = (Numeric)i;
70 }
71
72 v2[Range(joker)] = 88;
73
74 v2[Range(0, 2)] = 77;
75
76 cout << "v = \n" << v << "\n";
77 cout << "v2 = \n" << v2 << "\n";
78 cout << "v2.nelem() = \n" << v2.nelem() << "\n";
79
80 Vector v3;
81 v3.resize(v2.nelem());
82 v3 = v2;
83
84 cout << "\nv3 = \n" << v3 << "\n";
86 cout << "\nv3 after junking v2 = \n" << v3 << "\n";
87 v3 *= 2;
88 cout << "\nv3 after *2 = \n" << v3 << "\n";
89
90 Matrix M(10, 15);
91 {
92 Numeric n = 0;
93 for (Index i = 0; i < M.nrows(); ++i)
94 for (Index j = 0; j < M.ncols(); ++j) M(i, j) = ++n;
95 }
96
97 cout << "\nM =\n" << M << "\n";
98
99 cout << "\nM(Range(2,4),Range(2,4)) =\n"
100 << M(Range(2, 4), Range(2, 4)) << "\n";
101
102 cout << "\nM(Range(2,4),Range(2,4))(Range(1,2),Range(1,2)) =\n"
103 << M(Range(2, 4), Range(2, 4))(Range(1, 2), Range(1, 2)) << "\n";
104
105 cout << "\nM(1,Range(joker)) =\n" << M(1, Range(joker)) << "\n";
106
107 cout << "\nFilling M(1,Range(1,2)) with junk.\n";
108 fill_with_junk(M(1, Range(1, 2)));
109
110 cout << "\nM(Range(0,4),Range(0,4)) =\n"
111 << M(Range(0, 4), Range(0, 4)) << "\n";
112
113 cout << "\nFilling M(Range(4,2,2),Range(6,3)) with junk.\n";
114
115 MatrixView s = M(Range(4, 2, 2), Range(6, 3));
117
118 cout << "\nM =\n" << M << "\n";
119
120 const Matrix C = M;
121
122 cout << "\nC(Range(3,4,2),Range(2,3,3)) =\n"
123 << C(Range(3, 4, 2), Range(2, 3, 3)) << "\n";
124
125 cout << "\nC(Range(3,4,2),Range(2,3,3)).transpose() =\n"
126 << transpose(C(Range(3, 4, 2), Range(2, 3, 3))) << "\n";
127
128 return 0;
129}
130
131void test2() {
132 Vector v(50000000);
133
134 cout << "v.nelem() = " << v.nelem() << "\n";
135
136 cout << "Filling\n";
137 // for (Index i=0; i<v.nelem(); ++i )
138 // v[i] = sqrt(i);
139 v = 1.;
140 cout << "Done\n";
141}
142
143void test4() {
144 Vector a(10);
145 Vector b(a.nelem());
146
147 for (Index i = 0; i < a.nelem(); ++i) {
148 a[i] = (Numeric)(i + 1);
149 b[i] = (Numeric)(a.nelem() - i);
150 }
151
152 cout << "a = \n" << a << "\n";
153 cout << "b = \n" << b << "\n";
154 cout << "a*b \n= " << a * b << "\n";
155
156 Matrix A(11, 6);
157 Matrix B(10, 20);
158 Matrix C(20, 5);
159
160 B = 2;
161 C = 3;
162 mult(A(Range(1, joker), Range(1, joker)), B, C);
163
164 // cout << "\nB =\n" << B << "\n";
165 // cout << "\nC =\n" << C << "\n";
166 cout << "\nB*C =\n" << A << "\n";
167}
168
169void test5() {
170 Vector a(10);
171 Vector b(20);
172 Matrix M(10, 20);
173
174 // Fill b and M with a constant number:
175 b = 1;
176 M = 2;
177
178 cout << "b = \n" << b << "\n";
179 cout << "M =\n" << M << "\n";
180
181 mult(a, M, b); // a = M*b
182 cout << "\na = M*b = \n" << a << "\n";
183
184 mult(transpose((MatrixView)b), transpose((MatrixView)a), M); // b^t = a^t * M
185 cout << "\nb^t = a^t * M = \n" << transpose((MatrixView)b) << "\n";
186}
187
188void test6() {
189 Index n = 5000;
190 Vector x(1, n, 1), y(n);
191 Matrix M(n, n);
192 M = 1;
193 // cout << "x = \n" << x << "\n";
194
195 cout << "Transforming.\n";
196 // transform(x,sin,x);
197 // transform(transpose(y),sin,transpose(x));
198 // cout << "sin(x) =\n" << y << "\n";
199 for (Index i = 0; i < 1000; ++i) {
200 // mult(y,M,x);
201 transform(y, sin, static_cast<MatrixView>(x));
202 x += 1;
203 }
204 // cout << "y =\n" << y << "\n";
205
206 cout << "Done.\n";
207}
208
209void test7() {
210 Vector x(1, 20000000, 1);
211 Vector y(x.nelem());
212 transform(y, sin, x);
213 cout << "min(sin(x)), max(sin(x)) = " << min(y) << ", " << max(y) << "\n";
214}
215
216void test8() {
217 Vector x(80000000);
218 for (Index i = 0; i < x.nelem(); ++i) x[i] = (Numeric)i;
219 cout << "Done."
220 << "\n";
221}
222
223void test9() {
224 // Initialization of Matrix with view of other Matrix:
225 Matrix A(4, 8);
226 Matrix B(A(Range(joker), Range(0, 3)));
227 cout << "B = " << B << "\n";
228}
229
230void test10() {
231 // Initialization of Matrix with a vector (giving a 1 column Matrix).
232
233 // At the moment doing this with a non-const Vector will result in a
234 // warning message.
235 Vector v(1, 8, 1);
236 Matrix M((const Vector)v);
237 cout << "M = " << M << "\n";
238}
239
240void test11() {
241 // Assignment between Vector and Matrix:
242
243 // At the moment doing this with a non-const Vector will result in a
244 // warning message.
245 Vector v(1, 8, 1);
246 Matrix M(v.nelem(), 1);
247 M = v;
248 cout << "M = " << M << "\n";
249}
250
251void test12() {
252 // Copying of Arrays
253
254 Array<String> sa(3);
255 sa[0] = "It's ";
256 sa[1] = "a ";
257 sa[2] = "test.";
258
259 Array<String> sb(sa), sc(sa.nelem());
260
261 cout << "sb = \n" << sb << "\n";
262
263 sc = sa;
264
265 cout << "sc = \n" << sc << "\n";
266}
267
268void test13() {
269 // Mix vector and one-column matrix in += operator.
270 const Vector v(1, 8, 1); // The const is necessary here to
271 // avoid compiler warnings about
272 // different conversion paths.
273 Matrix M(v);
274 M += v;
275 cout << "M = \n" << M << "\n";
276}
277
278void test14() {
279 // Test explicit Array constructors.
280 Array<String> a{"Test"};
281 Array<Index> b{1, 2};
282 Array<Numeric> c{1, 2, 3, 4, 5, 6, 7, 8, 9, 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 0};
283 cout << "a = \n" << a << "\n";
284 cout << "b = \n" << b << "\n";
285 cout << "c = \n" << c << "\n";
286}
287
288void test15() {
289 // Test String.
290 String a = "Nur ein Test.";
291 cout << "a = " << a << "\n";
292 String b(a, 5, -1);
293 cout << "b = " << b << "\n";
294}
295
296/*void test16()
297{
298 // Test interaction between Array<Numeric> and Vector.
299 Vector a;
300 Array<Numeric> b;
301 b.push_back(1);
302 b.push_back(2);
303 b.push_back(3);
304 a.resize(b.nelem());
305 a = b;
306 cout << "b =\n" << b << "\n";
307 cout << "a =\n" << a << "\n";
308}*/
309
310void test17() {
311 // Test Sum.
312 Vector a(1, 10, 1);
313 cout << "a.sum() = " << a.sum() << "\n";
314}
315
316void test18() {
317 // Test elementvise square of a vector:
318 Vector a(1, 10, 1);
319 a *= a;
320 cout << "a *= a =\n" << a << "\n";
321}
322
323void test19() {
324 // There exists no explicit filling constructor of the form
325 // Vector a(3,1.7).
326 // But you can use the more general filling constructor with 3 arguments.
327
328 Vector a(1, 10, 1);
329 Vector b(5.3, 10, 0);
330 cout << "a =\n" << a << "\n";
331 cout << "b =\n" << b << "\n";
332}
333
334void test20() {
335 // Test initialization list constructor:
336 Vector a{1, 2, 3, 4, 5, 6, 7, 8, 9, 10};
337 cout << "a =\n" << a << "\n";
338}
339
340void test21() {
341 Numeric s = 0;
342 // Test speed of call by reference:
343 cout << "By reference:\n";
344 for (Index i = 0; i < (Index)1e8; ++i) {
345 s += by_reference(s);
346 s -= by_reference(s);
347 }
348 cout << "s = " << s << "\n";
349}
350
351void test22() {
352 Numeric s = 0;
353 // Test speed of call by value:
354 cout << "By value:\n";
355 for (Index i = 0; i < (Index)1e8; ++i) {
356 s += by_value(s);
357 s -= by_value(s);
358 }
359 cout << "s = " << s << "\n";
360}
361
362void test23() {
363 // Test constructors that fills with constant:
364 Vector a(10, 3.5);
365 cout << "a =\n" << a << "\n";
366 Matrix b(10, 10, 4.5);
367 cout << "b =\n" << b << "\n";
368}
369
370void test24() {
371 // Try element-vise multiplication of Matrix and Vector:
372 Matrix a(5, 1, 2.5);
373 Vector b(1, 5, 1);
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 a -= b;
381 cout << "a-=b =\n" << a << "\n";
382}
383
384void test25() {
385 // Test min and max for Array:
386 Array<Index> a{1, 2, 3, 4, 5, 6, 5, 4, 3, 2, 1};
387 cout << "min/max of a = " << min(a) << "/" << max(a) << "\n";
388}
389
390void test26() {
391 cout << "Test filling constructor for Array:\n";
392 Array<String> a(4, "Hello");
393 cout << "a =\n" << a << "\n";
394}
395
396void test27() {
397 cout << "Test Arrays of Vectors:\n";
399 a.push_back({1.0, 2.0});
400 a.push_back(Vector(1.0, 10, 1.0));
401 cout << "a =\n" << a << "\n";
402}
403
404void test28() {
405 cout << "Test default constructor for Matrix:\n";
406 Matrix a;
407 Matrix b(a);
408 cout << "b =\n" << b << "\n";
409}
410
411void test29() {
412 cout << "Test Arrays of Matrix:\n";
414 Matrix b;
415
416 b.resize(2, 2);
417 b(0, 0) = 1;
418 b(0, 1) = 2;
419 b(1, 0) = 3;
420 b(1, 1) = 4;
421 a.push_back(b);
422 b *= 2;
423 a.push_back(b);
424
425 a[0].resize(2, 3);
426 a[0] = 4;
427
428 a.resize(3);
429 a[2].resize(4, 5);
430 a[2] = 5;
431
432 cout << "a =\n" << a << "\n";
433}
434
435void test30() {
436 cout << "Test Matrices of size 0:\n";
437 Matrix a(0, 0);
438 // cout << "a(0,0) =\n" << a(0,0) << "\n";
439 a.resize(2, 2);
440 a = 1;
441 cout << "a =\n" << a << "\n";
442
443 Matrix b(3, 0);
444 // cout << "b(0,0) =\n" << b(0,0) << "\n";
445 b.resize(b.nrows(), b.ncols() + 3);
446 b = 2;
447 cout << "b =\n" << b << "\n";
448
449 Matrix c(0, 3);
450 // cout << "c(0,0) =\n" << c(0,0) << "\n";
451 c.resize(c.nrows() + 3, c.ncols());
452 c = 3;
453 cout << "c =\n" << c << "\n";
454}
455
456void test31() {
457 cout << "Test Tensor3:\n";
458
459 Tensor3 a(2, 3, 4, 1.0);
460
461 Index fill = 0;
462
463 // Fill with some numbers
464 for (Index i = 0; i < a.npages(); ++i)
465 for (Index j = 0; j < a.nrows(); ++j)
466 for (Index k = 0; k < a.ncols(); ++k) a(i, j, k) = (Numeric)(++fill);
467
468 cout << "a =\n" << a << "\n";
469
470 cout << "Taking out first row of first page:\n"
471 << a(0, 0, Range(joker)) << "\n";
472
473 cout << "Taking out last column of second page:\n"
474 << a(1, Range(joker), a.ncols() - 1) << "\n";
475
476 cout << "Taking out the first letter on every page:\n"
477 << a(Range(joker), 0, 0) << "\n";
478
479 cout << "Taking out first page:\n"
480 << a(0, Range(joker), Range(joker)) << "\n";
481
482 cout << "Taking out last row of all pages:\n"
483 << a(Range(joker), a.nrows() - 1, Range(joker)) << "\n";
484
485 cout << "Taking out second column of all pages:\n"
486 << a(Range(joker), Range(joker), 1) << "\n";
487
488 a *= 2;
489
490 cout << "After element-vise multiplication with 2:\n" << a << "\n";
491
492 transform(a, sqrt, a);
493
494 cout << "After taking the square-root:\n" << a << "\n";
495
496 Index s = 200;
497 cout << "Let's allocate a large tensor, "
498 << (Numeric)(s * s * s * 8) / 1024. / 1024. << " MB...\n";
499
500 a.resize(s, s, s);
501
502 cout << "Set it to 1...\n";
503
504 a = 1;
505
506 cout << "a(900,900,900) = " << a(90, 90, 90) << "\n";
507
508 fill = 0;
509
510 cout << "Fill with running numbers, using for loops...\n";
511 for (Index i = 0; i < a.npages(); ++i)
512 for (Index j = 0; j < a.nrows(); ++j)
513 for (Index k = 0; k < a.ncols(); ++k) a(i, j, k) = (Numeric)(++fill);
514
515 cout << "Max(a) = ...\n";
516
517 cout << max(a) << "\n";
518}
519
520void test32() {
521 cout << "Test von X = A*X:\n";
522 Matrix X(3, 3), A(3, 3), B(3, 3);
523
524 for (Index j = 0; j < A.nrows(); ++j)
525 for (Index k = 0; k < A.ncols(); ++k) {
526 X(j, k) = 1;
527 A(j, k) = (Numeric)(j + k);
528 }
529 cout << "A:\n" << A << "\n";
530 cout << "X:\n" << X << "\n";
531
532 mult(B, A, X);
533 cout << "B = A*X:\n" << B << "\n";
534 mult(X, A, X);
535 cout << "X = A*X:\n" << X << "\n";
536
537 cout << "This is not the same, and should not be, because you\n"
538 << "are not allowed to use the same matrix as input and output\n"
539 << "for mult!\n";
540}
541
542void test33() {
543 cout << "Making things look bigger than they are...\n";
544
545 {
546 cout << "1. Make a scalar look like a vector:\n";
547 Numeric a = 3.1415; // Just any number here.
548 VectorView av(a);
549 cout << "a, viewed as a vector: " << av << "\n";
550 cout << "Describe a: " << describe(a) << "\n";
551 av[0] += 1;
552 cout << "a, after the first element\n"
553 << "of the vector has been increased by 1: " << a << "\n";
554 }
555
556 {
557 cout
558 << "\n2. Make a vector look like a matrix:\n"
559 << "This is an exception, because the new dimension is added at the end.\n";
560 Vector a{1, 2, 3, 4, 5};
561 MatrixView am = a;
562 cout << "a, viewed as a matrix:\n" << am << "\n";
563 cout << "Trasnpose view:\n" << transpose(am) << "\n";
564 cout << "Describe transpose(am): " << describe(transpose(am)) << "\n";
565
566 cout << "\n3. Make a matrix look like a tensor3:\n";
567 Tensor3View at3 = am;
568 cout << "at3 = \n" << at3 << "\n";
569 cout << "Describe at3: " << describe(at3) << "\n";
570 at3(0, 2, 0) += 1;
571 cout << "a after Increasing element at3(0,2,0) by 1: \n" << a << "\n\n";
572
573 Tensor4View at4 = at3;
574 cout << "at4 = \n" << at4 << "\n";
575 cout << "Describe at4: " << describe(at4) << "\n";
576
577 Tensor5View at5 = at4;
578 cout << "at5 = \n" << at5 << "\n";
579 cout << "Describe at5: " << describe(at5) << "\n";
580
581 Tensor6View at6 = at5;
582 cout << "at6 = \n" << at6 << "\n";
583 cout << "Describe at6: " << describe(at6) << "\n";
584
585 Tensor7View at7 = at6;
586 cout << "at7 = \n" << at7 << "\n";
587 cout << "Describe at7: " << describe(at7) << "\n";
588
589 at7(0, 0, 0, 0, 0, 2, 0) -= 1;
590
591 cout << "After subtracting 1 from at7(0,0,0,0,0,2,0)\n"
592 << "a = " << a << "\n";
593
594 cout << "\nAll in one go:\n";
595 Numeric b = 3.1415; // Just any number here.
598 cout << "bt7:\n" << bt7 << "\n";
599 cout << "Describe bt7: " << describe(bt7) << "\n";
600 }
601}
602
603void junk4(Tensor4View a) { cout << "Describe a: " << describe(a) << "\n"; }
604
605void junk2(ConstVectorView a) { cout << "Describe a: " << describe(a) << "\n"; }
606
607void test34() {
608 cout << "Test, if dimension expansion works implicitly.\n";
609
610 Tensor3 t3(2, 3, 4);
611 junk4(t3);
612
613 Numeric x;
615}
616
617void test35() {
618 cout << "Test the new copy semantics.\n";
619 Vector a(1, 4, 1);
620 Vector b;
621
622 b = a;
623 cout << "b = " << b << "\n";
624
625 Vector aa(1, 5, 1);
626 ConstVectorView c = aa;
627 b = c;
628 cout << "b = " << b << "\n";
629
630 Vector aaa(1, 6, 1);
631 VectorView d = aaa;
632 b = d;
633 cout << "b = " << b << "\n";
634}
635
636void test36() {
637 cout << "Test using naked joker on Vector.\n";
638 Vector a(1, 4, 1);
639 VectorView b = a[joker];
640 cout << "a = " << a << "\n";
641 cout << "b = " << b << "\n";
642}
643
644void test37(const Index& i) {
645 Vector v1(5e-15, 10, 0.42e-15 / 11);
646 Vector v2 = v1;
647 // Index i = 10000000;
648
649 v1 /= (Numeric)i;
650 v2 /= (Numeric)i;
651 cout.precision(12);
652 // cout.setf(ios_base::scientific,ios_base::floatfield);
653 v1 *= v1;
654 v2 *= v2;
655 cout << v1 << endl;
656 cout << v2 << endl;
657}
658
659void test38() {
660 Vector v(5, 0.);
661 Numeric* const a = v.get_c_array();
662
663 a[4] = 5.;
664
665 cout << v << endl;
666 cout << endl << "========================" << endl << endl;
667
668 Matrix m(5, 5, 0.);
669 Numeric* const b = m.get_c_array();
670
671 b[4] = 5.;
672
673 cout << m << endl;
674 cout << endl << "========================" << endl << endl;
675
676 Tensor3 t3(5, 6, 7, 0.);
677 Numeric* const c = t3.get_c_array();
678
679 c[6] = 5.;
680
681 cout << t3 << endl;
682}
683
684void test39() {
685 Vector v1(1, 5, 1), v2(5);
686
687 v2 = v1 * 2;
688 // Unfortunately, this thing compiles, but at least it gives an
689 // assertion failure at runtime. I think what happens is that it
690 // implicitly creates a one element vector out of the "2", then
691 // tries to do a scalar product with v1.
692
693 cout << v2 << endl;
694}
695
696void test40() {
697 Vector v(100);
698
699 v = 5;
700
701 cout << v << endl;
702}
703
704void test41() {
705 const Vector v1(10, 1);
706
707 ConstVectorView vv = v1[Range(0, 5)];
708
709 cout << "Vector: " << v1 << endl;
710 cout << "VectorView: " << vv << endl;
711
712 try {
713 vv = 2;
714 } catch (const std::runtime_error& e) {
715 std::cerr << e.what() << endl;
716 exit(EXIT_FAILURE);
717 }
718
719 cout << "VectorView: " << vv << endl;
720 cout << "Vector: " << v1 << endl;
721}
722
723// Test behaviour of VectorView::operator MatrixView, for which I have fixed
724// a bug today.
725// SAB 2013-01-18
726void test42() {
727 cout << "test42\n";
728 Vector x{1, 2, 3, 4, 5, 6, 7, 8, 9, 10};
729 cout << "x: " << x << endl;
730
731 VectorView y = x[Range(2, 4, 2)];
732 cout << "y: " << y << endl;
733
735 cout << "z:\n" << z << endl;
736
737 cout << "Every other z:\n" << z(Range(1, 2, 2), joker) << endl;
738
739 // Try write access also:
740 MatrixView zz(y);
741 zz(Range(1, 2, 2), joker) = 0;
742 cout << "zz:\n" << zz << endl;
743 cout << "New x: " << x << endl;
744}
745
746
747// Test function for internal use, must return 2 as last n2r
748constexpr Rational test_numeric2rational(const Index i, const Index maxi, const Rational r=0, const Rational n2r=0) {
749 if (i > maxi)
750 return n2r;
751 else {
752 return
753 (r == n2r) ?
754 test_numeric2rational(i+1, maxi, Rational(maxi + i, maxi), numeric2rational(1 + Numeric(i)/Numeric(maxi), 12)) :
755 throw std::logic_error("Fail to initialize");;
756 }
757}
758
759
760void test43() {
761 // Simple construction compile-time test
762 constexpr Rational r=3_2; // should be 3/2
763 static_assert(r.numer == 3, "Rational fail to initialize properly");
764 static_assert(r.denom == 2, "Rational fail to initialize properly");
765
766 // Simple expression compile-time test
767 constexpr Rational r2 = 1 + r; // should be 5/2
768 static_assert(r2.numer == 5, "Rational fail to initialize properly");
769 static_assert(r2.denom == 2, "Rational fail to initialize properly");
770
771 // Div-zero by index generates an undefined rational: 0/0
772 constexpr Rational r3 = r / 0; // should be undefined
773 static_assert(r3.numer == 0, "Rational fail to initialize properly");
774 static_assert(r3.denom == 0, "Rational fail to initialize properly");
775 static_assert(r3 not_eq r3, "Undefined rational does not equal self");
776 static_assert(r3.isUndefined(), "Rational fail to initialize properly");
777
778 // Mul-zero gives 0/1
779 constexpr Rational r4 = r * 0; // should be undefined
780 static_assert(r4.numer == 0, "Rational fail to initialize properly");
781 static_assert(r4.denom == 1, "Rational fail to initialize properly");
782
783 // Complicated expression compile-time test
784 constexpr Rational r5 = (Rational(1, 9) * r) % 3; // should be 1/6
785 static_assert(r5.numer == 1, "Rational fail to initialize properly");
786 static_assert(r5.denom == 6, "Rational fail to initialize properly");
787
788 // The simplify operation simplifies expressions
789 constexpr Rational r6 = 10 % r; // should be 1
790 static_assert(r6.numer == 1, "Rational fail to initialize properly");
791 static_assert(r6.denom == 1, "Rational fail to initialize properly");
792 static_assert(r6.toInt() == 1, "Rational fail to initialize properly");
793 static_assert(r6.toIndex() == 1, "Rational fail to initialize properly");
794 static_assert(r6.toNumeric() == 1e0, "Rational fail to initialize properly");
795
796 constexpr Index rattest=1<<8;
797 constexpr Rational r7 = test_numeric2rational(0, rattest);
798 static_assert(r7 == 2, "Rational fail to initialize properly");
799}
800
801void test44() {
802#define docheck(fn, val, expect) \
803 cout << #fn << "(" << val << ") = " << fn(x) << " (expected " << #expect \
804 << ")" << std::endl;
805
806 Vector x{1, 2, 3};
808 docheck(is_sorted, x, 1)
809
810 x = {3, 2, 1};
812 docheck(is_sorted, x, 0)
813
814 x = {1, 2, 2};
816 docheck(is_sorted, x, 1)
817
818 x = {2, 2, 1};
820 docheck(is_sorted, x, 0)
821
822 x = {1, NAN, 2};
824 docheck(is_sorted, x, 0)
825
826 x = {2, NAN, 1};
828 docheck(is_sorted, x, 0)
829
830 x = {NAN, NAN, NAN};
832 docheck(is_sorted, x, 0)
833
834#undef docheck
835}
836
837void test46() {
838 Vector v(5, 0.);
839 nlinspace(v, 1, 10, 10);
840 VectorView v1 = v;
841 auto compare_func = [](Numeric n) { return n != 0; };
842 cout << std::any_of(v1.begin(), v1.end(), compare_func) << endl;
843 v1 = 0.;
844 cout << std::any_of(v1.begin(), v1.end(), compare_func) << endl;
845}
846
848
856bool test_diagonal(Index ntests) {
857 Rand<Index> rand(1, 100);
858
859 Index n_diag;
860
861 bool pass = true;
862
863 for (Index i = 0; i < ntests; i++) {
864 Matrix A(rand(), rand());
865 n_diag = std::min(A.ncols(), A.nrows());
866
868 ConstMatrixView B = A;
869 if (n_diag > 3) {
870 B = A(Range(1, Joker(), 2), Range(1, Joker(), 2));
871 }
872
874 ConstVectorView vt(AT.diagonal());
876
877 random_fill_matrix(A, 10, true);
878
879 for (Index j = 0; j < n_diag; j++) {
880 pass = pass && (v[j] == A(j, j));
881 pass = pass && (vt[j] == AT(j, j));
882
883 if (j < vb.nelem()) pass = pass && (vb[j] == B(j, j));
884 }
885 cout << endl;
886 }
887
888 if (pass)
889 cout << "test diagonal: Passed all tests." << endl;
890 else
891 cout << "test diagonal: Failed." << endl;
892
893 return pass;
894}
895
896Numeric matrix_vector_mult(Index m, Index n, Index ntests, bool verbose) {
897 Numeric max_err = 0;
898 Matrix A(m, n);
899 Vector x(n), xt, x_ref(n), y(m), y_ref(m);
900
901 Rand<Index> random_row_stride(1, n / 4);
902 Rand<Index> random_column_stride(1, m / 4);
903
904 for (Index i = 0; i < ntests; i++) {
905 // A * x
906
907 random_fill_matrix(A, 1000, false);
908 random_fill_vector(x, 1000, false);
909
910 mult(y, A, x);
911 mult_general(y_ref, A, x);
912
913 Numeric err_mul = get_maximum_error(y, y_ref, true);
914 if (err_mul > max_err) max_err = err_mul;
915
916 if (verbose) {
917 cout << "\t A * x: max. rel. error = " << err_mul << endl;
918 }
919
920 // A^T * y
921
922 mult(x, transpose(A), y);
923 mult_general(x_ref, transpose(A), y);
924
925 err_mul = get_maximum_error(x, x_ref, true);
926 if (err_mul > max_err) max_err = err_mul;
927
928 if (verbose) {
929 cout << "\t A^T * x: max. rel. error = " << err_mul << endl;
930 }
931
932 // Random stride
933 Index column_stride(random_column_stride()),
934 row_stride(random_row_stride());
935
936 Index m_sub, n_sub;
937 m_sub = (m - 1) / row_stride + 1;
938 n_sub = (n - 1) / column_stride + 1;
939
940 mult(y[Range(0, joker, row_stride)],
941 A(Range(0, m_sub), Range(0, n_sub)),
942 x[Range(0, joker, column_stride)]);
943 mult_general(y_ref[Range(0, joker, row_stride)],
944 A(Range(0, m_sub), Range(0, n_sub)),
945 x[Range(0, joker, column_stride)]);
946
947 err_mul = get_maximum_error(y[Range(0, joker, row_stride)],
948 y_ref[Range(0, joker, row_stride)],
949 true);
950 if (err_mul > max_err) max_err = err_mul;
951
952 if (verbose) {
953 cout << "\t Random stride: max. rel. error = " << err_mul << endl << endl;
954 }
955
956 // Random offset
957 if ((m > 1) && (n > 1)) {
958 Index y_offset = rand() % (m - 1);
959 Index x_offset = rand() % (n - 1);
960
961 m_sub = m - y_offset - 1;
962 n_sub = n - x_offset - 1;
963
964 mult(y[Range(y_offset, m_sub)],
965 A(Range(y_offset, m_sub), Range(x_offset, n_sub)),
966 x[Range(x_offset, n_sub)]);
967 mult_general(y_ref[Range(y_offset, m_sub)],
968 A(Range(y_offset, m_sub), Range(x_offset, n_sub)),
969 x[Range(x_offset, n_sub)]);
970
971 err_mul = get_maximum_error(
972 y[Range(y_offset, m_sub)], y_ref[Range(y_offset, m_sub)], true);
973 if (err_mul > max_err) max_err = err_mul;
974
975 if (verbose) {
976 cout << "\t Random offset: max. rel. error = " << err_mul << endl
977 << endl;
978 }
979 }
980 }
981
982 return max_err;
983}
984
986 Numeric err, max_err;
987
988 if (verbose)
989 cout << "Matrix-Vector Multiplication: n = m = 100, ntests = 100" << endl;
990
991 max_err = matrix_vector_mult(100, 100, 100, verbose);
992 if (verbose) {
993 cout << endl;
994 cout << "Matrix-Vector Multiplication: n = 100, m = 20, ntests = 100"
995 << endl;
996 }
997
998 err = matrix_vector_mult(100, 20, 100, verbose);
999 if (err > max_err) max_err = err;
1000 if (verbose) {
1001 cout << endl;
1002 cout << "Matrix-Vector Multiplication: n = 20, m = 100, ntests = 100"
1003 << endl;
1004 }
1005
1006 err = matrix_vector_mult(20, 100, 100, verbose);
1007 if (err > max_err) max_err = err;
1008 if (verbose) {
1009 if (max_err < 1e-9)
1010 cout << endl << "Matrix Vector Multiplication: PASSED" << endl;
1011 else
1012 cout << endl << "Matrix Vector Multiplication: FAILED" << endl;
1013 }
1014
1015 err = matrix_vector_mult(100, 1, 100, verbose);
1016 if (err > max_err) max_err = err;
1017 if (verbose) {
1018 cout << endl;
1019 cout << "Matrix-Vector Multiplication: n = 20, m = 100, ntests = 100"
1020 << endl;
1021 }
1022
1023 err = matrix_vector_mult(1, 100, 100, verbose);
1024 if (err > max_err) max_err = err;
1025 if (verbose) {
1026 cout << endl;
1027 cout << "Matrix-Vector Multiplication: n = 20, m = 100, ntests = 100"
1028 << endl;
1029 }
1030
1031 if (verbose) {
1032 if (max_err < 1e-9)
1033 cout << endl << "Matrix Vector Multiplication: PASSED" << endl;
1034 else
1035 cout << endl << "Matrix Vector Multiplication: FAILED" << endl;
1036 }
1037 return max_err;
1038}
1039
1041
1061 Index k, Index m, Index n, Index ntests, Index nsubtests, bool verbose) {
1062 Numeric max_err = 0;
1063 Matrix A(m, k), B(k, n), C(m, n), C_ref(m, n);
1064
1065 for (Index i = 0; i < ntests; i++) {
1066 random_fill_matrix(A, 1000, false);
1067 random_fill_matrix(B, 1000, false);
1068
1069 if (verbose) {
1070 cout.precision(15);
1071 cout << "MATRIX MULTIPLICATION: m = " << m << ", k = " << k;
1072 cout << ", n = " << n << endl;
1073 }
1074
1075 // A * B
1076
1077 random_fill_matrix(C, 10, false);
1078 random_fill_matrix(C_ref, 10, false);
1079
1080 mult(C, A, B);
1081 mult_general(C_ref, A, B);
1082
1083 Numeric err_mul = get_maximum_error(C, C_ref, true);
1084 if (err_mul > max_err) max_err = err_mul;
1085
1086 if (verbose) {
1087 cout << "\t"
1088 << "A * B: max. rel. error = " << err_mul << endl;
1089 }
1090
1091 // A^T * B
1092
1093 Matrix AT(k, m);
1094 random_fill_matrix(AT, 100.0, false);
1095
1096 mult(C, transpose(AT), B);
1097 mult_general(C_ref, transpose(AT), B);
1098 Numeric err_trans1 = get_maximum_error(C, C_ref, true);
1099 if (err_trans1 > max_err) max_err = err_trans1;
1100
1101 if (verbose) {
1102 cout << "\t"
1103 << "A^T * B: max. rel. err = " << err_trans1 << endl;
1104 }
1105
1106 // A * B^T
1107
1108 Matrix BT(n, k);
1109 random_fill_matrix(BT, 100.0, false);
1110
1111 mult(C, A, transpose(BT));
1112 mult_general(C_ref, A, transpose(BT));
1113 Numeric err_trans2 = get_maximum_error(C, C_ref, true);
1114 if (err_trans2 > max_err) max_err = err_trans2;
1115
1116 if (verbose) {
1117 cout << "\t"
1118 << "A * B^T: max. rel. err = " << err_trans2 << endl;
1119 }
1120
1121 // A^T * B^T
1122
1123 mult(C, transpose(AT), transpose(BT));
1124 mult_general(C_ref, transpose(AT), transpose(BT));
1125 Numeric err_trans3 = get_maximum_error(C, C_ref, true);
1126 if (err_trans3 > max_err) max_err = err_trans3;
1127
1128 if (verbose) {
1129 cout << "\t"
1130 << "A^T * B^T: max. rel. err = " << err_trans3 << endl;
1131 cout << endl;
1132 }
1133 }
1134
1135 // Multiplication of submatrices.
1136 Rand<Index> k_rand(1, k - 1), m_rand(1, m - 1), n_rand(1, n - 1);
1137
1138 for (Index i = 0; i < nsubtests - 1; i++) {
1139 Index k1(k_rand()), m1(m_rand()), n1(n_rand());
1140
1141 Range r1(random_range(m1));
1142 Range r2(random_range(k1));
1143 Range r3(random_range(n1));
1144
1145 MatrixView C_sub(C(r1, r3));
1146 MatrixView C_sub_ref(C_ref(r1, r3));
1147
1148 ConstMatrixView A_sub(A(r1, r2));
1149 ConstMatrixView B_sub(B(r2, r3));
1150
1151 mult(C_sub, A_sub, B_sub);
1152 mult_general(C_sub_ref, A_sub, B_sub);
1153
1154 Numeric err = get_maximum_error(C_sub, C_sub_ref, true);
1155 if (err > max_err) max_err = err;
1156
1157 if (verbose) {
1158 cout << "\t"
1159 << "Submatrix multiplication: max. rel. err = " << err;
1160 cout << endl;
1161 }
1162 }
1163
1164 if (verbose) {
1165 cout << endl;
1166 }
1167
1168 return max_err;
1169}
1170
1172
1184 Numeric max_err, err;
1185
1186 // Trivial case k = m = n = 0.
1187 max_err = matrix_mult(0, 0, 0, 10, 0, verbose);
1188
1189 // k = 1, m = 1, n = 1.
1190 err = matrix_mult(1, 1, 1, 10, 0, verbose);
1191 if (err > max_err) max_err = err;
1192
1193 // k = 10, m = 1, n = 10.
1194 err = matrix_mult(10, 1, 10, 20, 20, verbose);
1195 if (err > max_err) max_err = err;
1196
1197 // k = 10, m = 1, n = 10.
1198 err = matrix_mult(10, 1, 10, 20, 20, verbose);
1199 if (err > max_err) max_err = err;
1200
1201 // k = 10, m = 1, n = 1.
1202 err = matrix_mult(10, 1, 1, 20, 20, verbose);
1203 if (err > max_err) max_err = err;
1204
1205 // k = 100, m = 100, n = 100.
1206 err = matrix_mult(200, 200, 200, 20, 20, verbose);
1207 if (err > max_err) max_err = err;
1208
1209 // k = 10, m = 100, n = 10.
1210 err = matrix_mult(10, 100, 10, 20, 20, verbose);
1211 if (err > max_err) max_err = err;
1212
1213 // k = 10, m = 100, n = 100.
1214 err = matrix_mult(10, 100, 100, 20, 20, verbose);
1215 if (err > max_err) max_err = err;
1216
1217 // k = 100, m = 100, n = 100, 100 submatrix multiplications.
1218 err = matrix_mult(100, 100, 100, 0, 100, verbose);
1219 if (err > max_err) max_err = err;
1220
1221 return max_err;
1222}
1223
1226 {
1227 cout << "Array" << endl;
1229 v.resize(0);
1230 cout << v.empty() << endl;
1231 v.resize(1);
1232 cout << v.empty() << endl;
1233 }
1234 {
1235 cout << "Vector" << endl;
1236 Vector v;
1237 v.resize(0);
1238 cout << v.empty() << endl;
1239 v.resize(1);
1240 cout << v.empty() << endl;
1241 }
1242 {
1243 cout << "Matrix" << endl;
1244 Matrix v;
1245 v.resize(0, 1);
1246 cout << v.empty() << endl;
1247 v.resize(1, 0);
1248 cout << v.empty() << endl;
1249 v.resize(1, 1);
1250 cout << v.empty() << endl;
1251 }
1252 {
1253 cout << "Tensor3" << endl;
1254 Tensor3 v;
1255 v.resize(0, 1, 1);
1256 cout << v.empty() << endl;
1257 v.resize(1, 0, 1);
1258 cout << v.empty() << endl;
1259 v.resize(1, 1, 0);
1260 cout << v.empty() << endl;
1261 v.resize(1, 1, 1);
1262 cout << v.empty() << endl;
1263 }
1264 {
1265 cout << "Tensor4" << endl;
1266 Tensor4 v;
1267 v.resize(0, 1, 1, 1);
1268 cout << v.empty() << endl;
1269 v.resize(1, 0, 1, 1);
1270 cout << v.empty() << endl;
1271 v.resize(1, 1, 0, 1);
1272 cout << v.empty() << endl;
1273 v.resize(1, 1, 1, 0);
1274 cout << v.empty() << endl;
1275 v.resize(1, 1, 1, 1);
1276 cout << v.empty() << endl;
1277 }
1278 {
1279 cout << "Tensor5" << endl;
1280 Tensor5 v;
1281 v.resize(0, 1, 1, 1, 1);
1282 cout << v.empty() << endl;
1283 v.resize(1, 0, 1, 1, 1);
1284 cout << v.empty() << endl;
1285 v.resize(1, 1, 0, 1, 1);
1286 cout << v.empty() << endl;
1287 v.resize(1, 1, 1, 0, 1);
1288 cout << v.empty() << endl;
1289 v.resize(1, 1, 1, 1, 0);
1290 cout << v.empty() << endl;
1291 v.resize(1, 1, 1, 1, 1);
1292 cout << v.empty() << endl;
1293 }
1294 {
1295 cout << "Tensor6" << endl;
1296 Tensor6 v;
1297 v.resize(0, 1, 1, 1, 1, 1);
1298 cout << v.empty() << endl;
1299 v.resize(1, 0, 1, 1, 1, 1);
1300 cout << v.empty() << endl;
1301 v.resize(1, 1, 0, 1, 1, 1);
1302 cout << v.empty() << endl;
1303 v.resize(1, 1, 1, 0, 1, 1);
1304 cout << v.empty() << endl;
1305 v.resize(1, 1, 1, 1, 0, 1);
1306 cout << v.empty() << endl;
1307 v.resize(1, 1, 1, 1, 1, 0);
1308 cout << v.empty() << endl;
1309 v.resize(1, 1, 1, 1, 1, 1);
1310 cout << v.empty() << endl;
1311 }
1312 {
1313 cout << "Tensor7" << endl;
1314 Tensor7 v;
1315 v.resize(0, 1, 1, 1, 1, 1, 1);
1316 cout << v.empty() << endl;
1317 v.resize(1, 0, 1, 1, 1, 1, 1);
1318 cout << v.empty() << endl;
1319 v.resize(1, 1, 0, 1, 1, 1, 1);
1320 cout << v.empty() << endl;
1321 v.resize(1, 1, 1, 0, 1, 1, 1);
1322 cout << v.empty() << endl;
1323 v.resize(1, 1, 1, 1, 0, 1, 1);
1324 cout << v.empty() << endl;
1325 v.resize(1, 1, 1, 1, 1, 0, 1);
1326 cout << v.empty() << endl;
1327 v.resize(1, 1, 1, 1, 1, 1, 0);
1328 cout << v.empty() << endl;
1329 v.resize(1, 1, 1, 1, 1, 1, 1);
1330 cout << v.empty() << endl;
1331 }
1332}
1333
1335 const Numeric start,
1336 const Numeric stop,
1337 const Index n) {
1338 ARTS_ASSERT(1 < n); // Number of points must be greater 1.
1339 x.resize(n);
1340 Numeric step = (stop - start) / ((double)n - 1);
1341 for (Index i = 0; i < n - 1; i++) x[i] = start + (double)i * step;
1342 x[n - 1] = stop;
1343}
1344
1345void test47() {
1346 // Selecting empty matpack dimensions with Joker shouldn't fail
1347 Vector v;
1348 Matrix m;
1349
1350 VectorView vv = v[joker];
1351 MatrixView mv = m(joker, joker);
1352
1353 std::cout << "vv.nelem: " << vv.nelem() << std::endl;
1354 std::cout << "mv.nrows: " << mv.nrows() << " ";
1355 std::cout << "mv.ncols: " << mv.ncols() << std::endl;
1356
1357 Matrix m2(0, 5);
1358 MatrixView mv2 = m2(joker, 3);
1359 std::cout << "mv2.nrows: " << mv2.nrows() << " ";
1360 std::cout << "mv2.ncols: " << mv2.ncols() << std::endl;
1361}
1362
1363void test48() {
1364
1365 // Test each simple reduction
1366 std::cout << Tensor7(2, 2, 2, 2, 2, 2, 1, 8).reduce_rank<0, 1, 2, 3, 4, 5>() << '\n';
1367 std::cout << Tensor7(2, 2, 2, 2, 2, 1, 1, 9).reduce_rank<0, 1, 2, 3, 4>() << '\n';
1368 std::cout << Tensor7(2, 2, 2, 2, 1, 1, 1, 10).reduce_rank<0, 1, 2, 3>() << '\n';
1369 std::cout << Tensor7(2, 2, 2, 1, 1, 1, 1, 11).reduce_rank<0, 1, 2>() << '\n';
1370 std::cout << Tensor7(2, 2, 1, 1, 1, 1, 1, 12).reduce_rank<0, 1>() << '\n';
1371 std::cout << Tensor7(2, 1, 1, 1, 1, 1, 1, 13).reduce_rank<0>() << '\n';
1372 std::cout << Tensor7(1, 1, 1, 1, 1, 1, 1, 14) << '\n';
1373 std::cout << Tensor6(2, 2, 2, 2, 2, 1, 15).reduce_rank<0, 1, 2, 3, 4>() << '\n';
1374 std::cout << Tensor6(2, 2, 2, 2, 1, 1, 16).reduce_rank<0, 1, 2, 3>() << '\n';
1375 std::cout << Tensor6(2, 2, 2, 1, 1, 1, 17).reduce_rank<0, 1, 2>() << '\n';
1376 std::cout << Tensor6(2, 2, 1, 1, 1, 1, 18).reduce_rank<0, 1>() << '\n';
1377 std::cout << Tensor6(2, 1, 1, 1, 1, 1, 19).reduce_rank<0>() << '\n';
1378 std::cout << Tensor6(1, 1, 1, 1, 1, 1, 20) << '\n';
1379 std::cout << Tensor5(2, 2, 2, 2, 1, 21).reduce_rank<0, 1, 2, 3>() << '\n';
1380 std::cout << Tensor5(2, 2, 2, 1, 1, 22).reduce_rank<0, 1, 2>() << '\n';
1381 std::cout << Tensor5(2, 2, 1, 1, 1, 23).reduce_rank<0, 1>() << '\n';
1382 std::cout << Tensor5(2, 1, 1, 1, 1, 24).reduce_rank<0>() << '\n';
1383 std::cout << Tensor5(1, 1, 1, 1, 1, 25) << '\n';
1384 std::cout << Tensor4(2, 2, 2, 1, 26).reduce_rank<0, 1, 2>() << '\n';
1385 std::cout << Tensor4(2, 2, 1, 1, 27).reduce_rank<0, 1>() << '\n';
1386 std::cout << Tensor4(2, 1, 1, 1, 28).reduce_rank<0>() << '\n';
1387 std::cout << Tensor4(1, 1, 1, 1, 29) << '\n';
1388 std::cout << Tensor3(2, 2, 1, 30).reduce_rank<0, 1>() << '\n';
1389 std::cout << Tensor3(2, 1, 1, 31).reduce_rank<0>() << '\n';
1390 std::cout << Tensor3(1, 1, 1, 32) << '\n';
1391 std::cout << Matrix(2, 1, 33).reduce_rank<0>() << '\n';
1392 std::cout << Matrix(1, 1, 34) << '\n';
1393
1394 // Test that the reductions work along different axis for Vector
1395 std::cout << matpack::eigen::mat(Tensor7(2, 1, 1, 1, 1, 1, 1, 35).reduce_rank<0, 1>()).transpose() << '\n'
1396 << matpack::eigen::mat(Tensor7(2, 1, 1, 1, 1, 1, 1, 35).reduce_rank<0, 2>()).transpose() << '\n'
1397 << matpack::eigen::mat(Tensor7(2, 1, 1, 1, 1, 1, 1, 35).reduce_rank<0, 3>()).transpose() << '\n'
1398 << matpack::eigen::mat(Tensor7(2, 1, 1, 1, 1, 1, 1, 35).reduce_rank<0, 4>()).transpose() << '\n'
1399 << matpack::eigen::mat(Tensor7(2, 1, 1, 1, 1, 1, 1, 35).reduce_rank<0, 5>()).transpose() << '\n'
1400 << matpack::eigen::mat(Tensor7(2, 1, 1, 1, 1, 1, 1, 35).reduce_rank<0, 6>()).transpose() << '\n';
1401
1402 // Test that the reductions work along different axis for Matrix
1403 std::cout << matpack::eigen::mat(Tensor7(2, 2, 1, 1, 1, 1, 1, 36).reduce_rank<0, 1>()) << '\n'
1404 << matpack::eigen::mat(Tensor7(2, 1, 2, 1, 1, 1, 1, 36).reduce_rank<0, 2>()) << '\n'
1405 << matpack::eigen::mat(Tensor7(2, 1, 1, 2, 1, 1, 1, 36).reduce_rank<0, 3>()) << '\n'
1406 << matpack::eigen::mat(Tensor7(2, 1, 1, 1, 2, 1, 1, 36).reduce_rank<0, 4>()) << '\n'
1407 << matpack::eigen::mat(Tensor7(2, 1, 1, 1, 1, 2, 1, 36).reduce_rank<0, 5>()) << '\n'
1408 << matpack::eigen::mat(Tensor7(2, 1, 1, 1, 1, 1, 2, 36).reduce_rank<0, 6>()) << '\n';
1409
1410 // Test that the right elements are accessed
1411 {
1412 Index val=1;
1413 Tensor7 testvar(9, 2, 4, 3, 5, 7, 11);for (Index i=0; i<9; i++)
1414 for (Index j=0; j<2; j++) {
1415 for (Index k=0; k<4; k++) {
1416 for (Index l=0; l<3; l++) {
1417 for (Index m=0; m<5; m++) {
1418 for (Index n=0; n<7; n++) {
1419 for (Index o=0; o<11; o++) {
1420 testvar(i, j, k, l, m, n, o) = Numeric(val++);
1421 }
1422 }
1423 }
1424 }
1425 }
1426 }
1427
1428 // Vector test
1429 for (Index j=0; j<2; j++) {
1430 for (Index k=0; k<4; k++) {
1431 for (Index l=0; l<3; l++) {
1432 for (Index m=0; m<5; m++) {
1433 for (Index n=0; n<7; n++) {
1434 for (Index o=0; o<11; o++) {
1435 std::cout << Tensor7(testvar(joker, Range(j, 1), Range(k, 1), Range(l, 1), Range(m, 1), Range(n, 1), Range(o, 1))).reduce_rank<0>() << '\n';
1436 std::cout << testvar(joker, j, k, l, m, n, o) << '\n';
1437 std::cout << '\n';
1438 }
1439 }
1440 }
1441 }
1442 }
1443 }
1444 for (Index i=0; i<9; i++) {
1445 for (Index k=0; k<4; k++) {
1446 for (Index l=0; l<3; l++) {
1447 for (Index m=0; m<5; m++) {
1448 for (Index n=0; n<7; n++) {
1449 for (Index o=0; o<11; o++) {
1450 std::cout << Tensor7(testvar(Range(i, 1), joker, Range(k, 1), Range(l, 1), Range(m, 1), Range(n, 1), Range(o, 1))).reduce_rank<1>() << '\n';
1451 std::cout << testvar(i, joker, k, l, m, n, o) << '\n';
1452 std::cout << '\n';
1453 }
1454 }
1455 }
1456 }
1457 }
1458 }
1459 for (Index i=0; i<9; i++) {
1460 for (Index j=0; j<2; j++) {
1461 for (Index l=0; l<3; l++) {
1462 for (Index m=0; m<5; m++) {
1463 for (Index n=0; n<7; n++) {
1464 for (Index o=0; o<11; o++) {
1465 std::cout << Tensor7(testvar(Range(i, 1), Range(j, 1), joker, Range(l, 1), Range(m, 1), Range(n, 1), Range(o, 1))).reduce_rank<2>() << '\n';
1466 std::cout << testvar(i, j, joker, l, m, n, o) << '\n';
1467 std::cout << '\n';
1468 }
1469 }
1470 }
1471 }
1472 }
1473 }
1474 for (Index i=0; i<9; i++) {
1475 for (Index j=0; j<2; j++) {
1476 for (Index k=0; k<4; k++) {
1477 for (Index m=0; m<5; m++) {
1478 for (Index n=0; n<7; n++) {
1479 for (Index o=0; o<11; o++) {
1480 std::cout << Tensor7(testvar(Range(i, 1), Range(j, 1), Range(k, 1), joker, Range(m, 1), Range(n, 1), Range(o, 1))).reduce_rank<3>() << '\n';
1481 std::cout << testvar(i, j, k, joker, m, n, o) << '\n';
1482 std::cout << '\n';
1483 }
1484 }
1485 }
1486 }
1487 }
1488 }
1489 for (Index i=0; i<9; i++) {
1490 for (Index j=0; j<2; j++) {
1491 for (Index k=0; k<4; k++) {
1492 for (Index l=0; l<3; l++) {
1493 for (Index n=0; n<7; n++) {
1494 for (Index o=0; o<11; o++) {
1495 std::cout << Tensor7(testvar(Range(i, 1), Range(j, 1), Range(k, 1), Range(l, 1), joker, Range(n, 1), Range(o, 1))).reduce_rank<4>() << '\n';
1496 std::cout << testvar(i, j, k, l, joker, n, o) << '\n';
1497 std::cout << '\n';
1498 }
1499 }
1500 }
1501 }
1502 }
1503 }
1504 for (Index i=0; i<9; i++) {
1505 for (Index j=0; j<2; j++) {
1506 for (Index k=0; k<4; k++) {
1507 for (Index l=0; l<3; l++) {
1508 for (Index m=0; m<5; m++) {
1509 for (Index o=0; o<11; o++) {
1510 std::cout << Tensor7(testvar(Range(i, 1), Range(j, 1), Range(k, 1), Range(l, 1), Range(m, 1), joker, Range(o, 1))).reduce_rank<5>() << '\n';
1511 std::cout << testvar(i, j, k, l, m, joker, o) << '\n';
1512 std::cout << '\n';
1513 }
1514 }
1515 }
1516 }
1517 }
1518 }
1519 for (Index i=0; i<9; i++) {
1520 for (Index j=0; j<2; j++) {
1521 for (Index k=0; k<4; k++) {
1522 for (Index l=0; l<3; l++) {
1523 for (Index m=0; m<5; m++) {
1524 for (Index n=0; n<7; n++) {
1525 std::cout << Tensor7(testvar(Range(i, 1), Range(j, 1), Range(k, 1), Range(l, 1), Range(m, 1), Range(n, 1), joker)).reduce_rank<6>() << '\n';
1526 std::cout << testvar(i, j, k, l, m, n, joker) << '\n';
1527 std::cout << '\n';
1528 }
1529 }
1530 }
1531 }
1532 }
1533 }
1534
1535 // Matrix test
1536 for (Index k=0; k<4; k++) {
1537 for (Index l=0; l<3; l++) {
1538 for (Index m=0; m<5; m++) {
1539 for (Index n=0; n<7; n++) {
1540 for (Index o=0; o<11; o++) {
1541 std::cout << Tensor7(testvar(joker, joker, Range(k, 1), Range(l, 1), Range(m, 1), Range(n, 1), Range(o, 1))).reduce_rank<0, 1>() << '\n';
1542 std::cout << testvar(joker, joker, k, l, m, n, o) << '\n';
1543 std::cout << '\n';
1544 }
1545 }
1546 }
1547 }
1548 }
1549 for (Index j=0; j<2; j++) {
1550 for (Index l=0; l<3; l++) {
1551 for (Index m=0; m<5; m++) {
1552 for (Index n=0; n<7; n++) {
1553 for (Index o=0; o<11; o++) {
1554 std::cout << Tensor7(testvar(joker, Range(j, 1), joker, Range(l, 1), Range(m, 1), Range(n, 1), Range(o, 1))).reduce_rank<0, 2>() << '\n';
1555 std::cout << testvar(joker, j, joker, l, m, n, o) << '\n';
1556 std::cout << '\n';
1557 }
1558 }
1559 }
1560 }
1561 }
1562 for (Index j=0; j<2; j++) {
1563 for (Index k=0; k<4; k++) {
1564 for (Index m=0; m<5; m++) {
1565 for (Index n=0; n<7; n++) {
1566 for (Index o=0; o<11; o++) {
1567 std::cout << Tensor7(testvar(joker, Range(j, 1), Range(k, 1), joker, Range(m, 1), Range(n, 1), Range(o, 1))).reduce_rank<0, 3>() << '\n';
1568 std::cout << testvar(joker, j, k, joker, m, n, o) << '\n';
1569 std::cout << '\n';
1570 }
1571 }
1572 }
1573 }
1574 }
1575 for (Index j=0; j<2; j++) {
1576 for (Index k=0; k<4; k++) {
1577 for (Index l=0; l<3; l++) {
1578 for (Index n=0; n<7; n++) {
1579 for (Index o=0; o<11; o++) {
1580 std::cout << Tensor7(testvar(joker, Range(j, 1), Range(k, 1), Range(l, 1), joker, Range(n, 1), Range(o, 1))).reduce_rank<0, 4>() << '\n';
1581 std::cout << testvar(joker, j, k, l, joker, n, o) << '\n';
1582 std::cout << '\n';
1583 }
1584 }
1585 }
1586 }
1587 }
1588 for (Index j=0; j<2; j++) {
1589 for (Index k=0; k<4; k++) {
1590 for (Index l=0; l<3; l++) {
1591 for (Index m=0; m<5; m++) {
1592 for (Index o=0; o<11; o++) {
1593 std::cout << Tensor7(testvar(joker, Range(j, 1), Range(k, 1), Range(l, 1), Range(m, 1), joker, Range(o, 1))).reduce_rank<0, 5>() << '\n';
1594 std::cout << testvar(joker, j, k, l, m, joker, o) << '\n';
1595 std::cout << '\n';
1596 }
1597 }
1598 }
1599 }
1600 }
1601 for (Index j=0; j<2; j++) {
1602 for (Index k=0; k<4; k++) {
1603 for (Index l=0; l<3; l++) {
1604 for (Index m=0; m<5; m++) {
1605 for (Index n=0; n<7; n++) {
1606 std::cout << Tensor7(testvar(joker, Range(j, 1), Range(k, 1), Range(l, 1), Range(m, 1), Range(n, 1), joker)).reduce_rank<0, 6>() << '\n';
1607 std::cout << testvar(joker, j, k, l, m, n, joker) << '\n';
1608 std::cout << '\n';
1609 }
1610 }
1611 }
1612 }
1613 }
1614 for (Index i=0; i<9; i++) {
1615 for (Index l=0; l<3; l++) {
1616 for (Index m=0; m<5; m++) {
1617 for (Index n=0; n<7; n++) {
1618 for (Index o=0; o<11; o++) {
1619 std::cout << Tensor7(testvar(Range(i, 1), joker, joker, Range(l, 1), Range(m, 1), Range(n, 1), Range(o, 1))).reduce_rank<1, 2>() << '\n';
1620 std::cout << testvar(i, joker, joker, l, m, n, o) << '\n';
1621 std::cout << '\n';
1622 }
1623 }
1624 }
1625 }
1626 }
1627 for (Index i=0; i<9; i++) {
1628 for (Index k=0; k<4; k++) {
1629 for (Index m=0; m<5; m++) {
1630 for (Index n=0; n<7; n++) {
1631 for (Index o=0; o<11; o++) {
1632 std::cout << Tensor7(testvar(Range(i, 1), joker, Range(k, 1), joker, Range(m, 1), Range(n, 1), Range(o, 1))).reduce_rank<1, 3>() << '\n';
1633 std::cout << testvar(i, joker, k, joker, m, n, o) << '\n';
1634 std::cout << '\n';
1635 }
1636 }
1637 }
1638 }
1639 }
1640 for (Index i=0; i<9; i++) {
1641 for (Index k=0; k<4; k++) {
1642 for (Index l=0; l<3; l++) {
1643 for (Index n=0; n<7; n++) {
1644 for (Index o=0; o<11; o++) {
1645 std::cout << Tensor7(testvar(Range(i, 1), joker, Range(k, 1), Range(l, 1), joker, Range(n, 1), Range(o, 1))).reduce_rank<1, 4>() << '\n';
1646 std::cout << testvar(i, joker, k, l, joker, n, o) << '\n';
1647 std::cout << '\n';
1648 }
1649 }
1650 }
1651 }
1652 }
1653 for (Index i=0; i<9; i++) {
1654 for (Index k=0; k<4; k++) {
1655 for (Index l=0; l<3; l++) {
1656 for (Index m=0; m<5; m++) {
1657 for (Index o=0; o<11; o++) {
1658 std::cout << Tensor7(testvar(Range(i, 1), joker, Range(k, 1), Range(l, 1), Range(m, 1), joker, Range(o, 1))).reduce_rank<1, 5>() << '\n';
1659 std::cout << testvar(i, joker, k, l, m, joker, o) << '\n';
1660 std::cout << '\n';
1661 }
1662 }
1663 }
1664 }
1665 }
1666 for (Index i=0; i<9; i++) {
1667 for (Index k=0; k<4; k++) {
1668 for (Index l=0; l<3; l++) {
1669 for (Index m=0; m<5; m++) {
1670 for (Index n=0; n<7; n++) {
1671 std::cout << Tensor7(testvar(Range(i, 1), joker, Range(k, 1), Range(l, 1), Range(m, 1), Range(n, 1), joker)).reduce_rank<1, 6>() << '\n';
1672 std::cout << testvar(i, joker, k, l, m, n, joker) << '\n';
1673 std::cout << '\n';
1674 }
1675 }
1676 }
1677 }
1678 }
1679 for (Index i=0; i<9; i++) {
1680 for (Index j=0; j<2; j++) {
1681 for (Index m=0; m<5; m++) {
1682 for (Index n=0; n<7; n++) {
1683 for (Index o=0; o<11; o++) {
1684 std::cout << Tensor7(testvar(Range(i, 1), Range(j, 1), joker, joker, Range(m, 1), Range(n, 1), Range(o, 1))).reduce_rank<2, 3>() << '\n';
1685 std::cout << testvar(i, j, joker, joker, m, n, o) << '\n';
1686 std::cout << '\n';
1687 }
1688 }
1689 }
1690 }
1691 }
1692 for (Index i=0; i<9; i++) {
1693 for (Index j=0; j<2; j++) {
1694 for (Index l=0; l<3; l++) {
1695 for (Index n=0; n<7; n++) {
1696 for (Index o=0; o<11; o++) {
1697 std::cout << Tensor7(testvar(Range(i, 1), Range(j, 1), joker, Range(l, 1), joker, Range(n, 1), Range(o, 1))).reduce_rank<2, 4>() << '\n';
1698 std::cout << testvar(i, j, joker, l, joker, n, o) << '\n';
1699 std::cout << '\n';
1700 }
1701 }
1702 }
1703 }
1704 }
1705 for (Index i=0; i<9; i++) {
1706 for (Index j=0; j<2; j++) {
1707 for (Index l=0; l<3; l++) {
1708 for (Index m=0; m<5; m++) {
1709 for (Index o=0; o<11; o++) {
1710 std::cout << Tensor7(testvar(Range(i, 1), Range(j, 1), joker, Range(l, 1), Range(m, 1), joker, Range(o, 1))).reduce_rank<2, 5>() << '\n';
1711 std::cout << testvar(i, j, joker, l, m, joker, o) << '\n';
1712 std::cout << '\n';
1713 }
1714 }
1715 }
1716 }
1717 }
1718 for (Index i=0; i<9; i++) {
1719 for (Index j=0; j<2; j++) {
1720 for (Index l=0; l<3; l++) {
1721 for (Index m=0; m<5; m++) {
1722 for (Index n=0; n<7; n++) {
1723 std::cout << Tensor7(testvar(Range(i, 1), Range(j, 1), joker, Range(l, 1), Range(m, 1), Range(n, 1), joker)).reduce_rank<2, 6>() << '\n';
1724 std::cout << testvar(i, j, joker, l, m, n, joker) << '\n';
1725 std::cout << '\n';
1726 }
1727 }
1728 }
1729 }
1730 }
1731 for (Index i=0; i<9; i++) {
1732 for (Index j=0; j<2; j++) {
1733 for (Index k=0; k<4; k++) {
1734 for (Index n=0; n<7; n++) {
1735 for (Index o=0; o<11; o++) {
1736 std::cout << Tensor7(testvar(Range(i, 1), Range(j, 1), Range(k, 1), joker, joker, Range(n, 1), Range(o, 1))).reduce_rank<3, 4>() << '\n';
1737 std::cout << testvar(i, j, k, joker, joker, n, o) << '\n';
1738 std::cout << '\n';
1739 }
1740 }
1741 }
1742 }
1743 }
1744 for (Index i=0; i<9; i++) {
1745 for (Index j=0; j<2; j++) {
1746 for (Index k=0; k<4; k++) {
1747 for (Index m=0; m<5; m++) {
1748 for (Index o=0; o<11; o++) {
1749 std::cout << Tensor7(testvar(Range(i, 1), Range(j, 1), Range(k, 1), joker, Range(m, 1), joker, Range(o, 1))).reduce_rank<3, 5>() << '\n';
1750 std::cout << testvar(i, j, k, joker, m, joker, o) << '\n';
1751 std::cout << '\n';
1752 }
1753 }
1754 }
1755 }
1756 }
1757 for (Index i=0; i<9; i++) {
1758 for (Index j=0; j<2; j++) {
1759 for (Index k=0; k<4; k++) {
1760 for (Index m=0; m<5; m++) {
1761 for (Index n=0; n<7; n++) {
1762 std::cout << Tensor7(testvar(Range(i, 1), Range(j, 1), Range(k, 1), joker, Range(m, 1), Range(n, 1), joker)).reduce_rank<3, 6>() << '\n';
1763 std::cout << testvar(i, j, k, joker, m, n, joker) << '\n';
1764 std::cout << '\n';
1765 }
1766 }
1767 }
1768 }
1769 }
1770 for (Index i=0; i<9; i++) {
1771 for (Index j=0; j<2; j++) {
1772 for (Index k=0; k<4; k++) {
1773 for (Index l=0; l<3; l++) {
1774 for (Index o=0; o<11; o++) {
1775 std::cout << Tensor7(testvar(Range(i, 1), Range(j, 1), Range(k, 1), Range(l, 1), joker, joker, Range(o, 1))).reduce_rank<4, 5>() << '\n';
1776 std::cout << testvar(i, j, k, l, joker, joker, o) << '\n';
1777 std::cout << '\n';
1778 }
1779 }
1780 }
1781 }
1782 }
1783 for (Index i=0; i<9; i++) {
1784 for (Index j=0; j<2; j++) {
1785 for (Index k=0; k<4; k++) {
1786 for (Index l=0; l<3; l++) {
1787 for (Index n=0; n<7; n++) {
1788 std::cout << Tensor7(testvar(Range(i, 1), Range(j, 1), Range(k, 1), Range(l, 1), joker, Range(n, 1), joker)).reduce_rank<4, 6>() << '\n';
1789 std::cout << testvar(i, j, k, l, joker, n, joker) << '\n';
1790 std::cout << '\n';
1791 }
1792 }
1793 }
1794 }
1795 }
1796 for (Index i=0; i<9; i++) {
1797 for (Index j=0; j<2; j++) {
1798 for (Index k=0; k<4; k++) {
1799 for (Index l=0; l<3; l++) {
1800 for (Index m=0; m<5; m++) {
1801 std::cout << Tensor7(testvar(Range(i, 1), Range(j, 1), Range(k, 1), Range(l, 1), Range(m, 1), joker, joker)).reduce_rank<5, 6>() << '\n';
1802 std::cout << testvar(i, j, k, l, m, joker, joker) << '\n';
1803 std::cout << '\n';
1804 }
1805 }
1806 }
1807 }
1808 }
1809 }
1810}
1811
1813 try {
1814 wigner3j(1, 0, 1, 0, 0, 0);
1815 } catch(std::exception& e) {
1816 std::cerr << e.what() << '\n';
1817 }
1818
1819 make_wigner_ready(250, 20000000, 3);
1820 wigner3j(1, 0, 1, 0, 0, 0);
1821}
1822
1824 std::vector<Index> x(30);
1825 std::iota(x.begin(), x.end(), -15);
1826 for (auto& i: x) std::cout << "-1^" << i << '=' << pow_negative_one(i) << '\n';
1827}
1828
1829int main() {
1830 // test1();
1831 // test2();
1832 // test3();
1833 // test4();
1834 // test5();
1835 // test6();
1836 // test7();
1837 // test8();
1838 // test9();
1839 // test10();
1840 // test11();
1841 // test12();
1842 // test13();
1843 // test14();
1844 // test15();
1845 // test16();
1846 // test17();
1847 // test18();
1848 // test19();
1849 // test20();
1850 // test21();
1851 // test22();
1852 // test23();
1853 // test24();
1854 // test25();
1855 // test26();
1856 // test27();
1857 // test28();
1858 // test29();
1859 // test30();
1860 // test31();
1861 // test32();
1862 // test33();
1863 // test34();
1864 // test35();
1865 // test36();
1866 // Index i = 10000000;
1867 // test37(i);
1868 // test38();
1869 // test39();
1870 // test40();
1871 // test41();
1872 // test42();
1873 // test43();
1874 // test44();
1875 // test45();
1876 // test46();
1877 // test47();
1878 //test48();
1879 //test_wigner_error();
1881
1882 // const double tolerance = 1e-9;
1883 // double error;
1884 //
1885 // // Matrix Vector Multiplication.
1886 // error = test_matrix_vector_multiplication(false);
1887 // cout << "Matrix Vector Multiplication: ";
1888 // if (error > tolerance)
1889 // cout << "FAILED, maximum error: " << error << endl;
1890 // else
1891 // cout << "PASSED." << endl;
1892 //
1893 // // Matrix Matrix Multiplication.
1894 // error = test_matrix_multiplication(false);
1895 //
1896 // cout << "Matrix Matrix Multiplication: ";
1897 // if (error > tolerance)
1898 // cout << "FAILED, maximum error: " << error << endl;
1899 // else
1900 // cout << "PASSED." << endl;
1901
1902 //test_diagonal( 100 );
1903 //test_empty();
1904
1905 return 1;
1906}
This file contains the definition of Array.
Index nelem() const ARTS_NOEXCEPT
Definition: array.h:92
A constant view of a Matrix.
Definition: matpackI.h:1043
Numeric * get_c_array() const noexcept
Definition: matpackI.h:1133
Index nrows() const noexcept
Definition: matpackI.h:1055
Index ncols() const noexcept
Definition: matpackI.h:1056
ConstVectorView diagonal() const ARTS_NOEXCEPT
Matrix diagonal as vector.
Definition: matpackI.cc:468
A constant view of a Vector.
Definition: matpackI.h:512
Index nelem() const noexcept
Returns the number of elements.
Definition: matpackI.h:536
The MatrixView class.
Definition: matpackI.h:1164
The Matrix class.
Definition: matpackI.h:1261
The range class.
Definition: matpackI.h:159
The Tensor3View class.
Definition: matpackIII.h:246
const Numeric * get_c_array() const ARTS_NOEXCEPT
Conversion to plain C-array.
Definition: matpackIII.cc:315
The Tensor3 class.
Definition: matpackIII.h:346
Vector reduce_rank() &&ARTS_NOEXCEPT
Definition: matpackIII.h:391
The Tensor4View class.
Definition: matpackIV.h:292
The Tensor4 class.
Definition: matpackIV.h:429
Vector reduce_rank() &&ARTS_NOEXCEPT
Definition: matpackIV.h:482
The Tensor5View class.
Definition: matpackV.h:343
The Tensor5 class.
Definition: matpackV.h:516
Vector reduce_rank() &&ARTS_NOEXCEPT
Definition: matpackV.h:569
The Tensor6View class.
Definition: matpackVI.h:632
The Tensor6 class.
Definition: matpackVI.h:1099
Vector reduce_rank() &&ARTS_NOEXCEPT
Definition: matpackVI.h:1156
The Tensor7View class.
Definition: matpackVII.h:1303
The Tensor7 class.
Definition: matpackVII.h:2399
Vector reduce_rank() &&ARTS_NOEXCEPT
Definition: matpackVII.h:2434
The VectorView class.
Definition: matpackI.h:663
Iterator1D begin() ARTS_NOEXCEPT
Return iterator to first element.
Definition: matpackI.cc:144
Iterator1D end() ARTS_NOEXCEPT
Return iterator behind last element.
Definition: matpackI.cc:148
The Vector class.
Definition: matpackI.h:899
void resize(Index n)
Resize function.
Definition: matpackI.cc:390
void mult(MatrixView C, ConstMatrixView A, const Block &B)
#define ARTS_ASSERT(condition,...)
Definition: debug.h:83
String describe(ConstTensor7View x)
Describe Tensor7.
Definition: describe.cc:38
Header file for describe.cc.
The declarations of all the exception classes.
#define min(a, b)
#define max(a, b)
bool is_increasing(ConstVectorView x)
Checks if a vector is sorted and strictly increasing.
Definition: logic.cc:218
bool is_sorted(ConstVectorView x)
Checks if a vector is sorted in ascending order.
Definition: logic.cc:202
bool is_decreasing(ConstVectorView x)
Checks if a vector is sorted in reversed order and is strictly decreasing.
Definition: logic.cc:255
Header file for logic.cc.
constexpr Index pow_negative_one(Index x) noexcept
Computes std::pow(-1, x) without std::pow.
Definition: math_funcs.h:158
Header file for sparse matrices.
void mult_general(VectorView y, const ConstMatrixView &M, const ConstVectorView &x) ARTS_NOEXCEPT
Matrix Vector multiplication.
Definition: matpackI.cc:1136
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:1431
NUMERIC Numeric
The type to use for all floating point numbers.
Definition: matpack.h:33
INDEX Index
The type to use for all integer numbers and indices.
Definition: matpack.h:39
ConstComplexMatrixView transpose(ConstComplexMatrixView m)
Const version of transpose.
const Joker joker
This file contains the definition of String, the ARTS string class.
constexpr Numeric k
Boltzmann constant convenience name [J/K].
constexpr Numeric e
Elementary charge convenience name [C].
constexpr Numeric l(const Index p0, const Index n, const Numeric x, const SortedVectorType &xi, const Index j, const std::pair< Numeric, Numeric > cycle={ -180, 180}) noexcept
constexpr std::string_view Joker
Definition: isotopologues.h:13
constexpr Rational start(Rational Ju, Rational Jl, Polarization type) noexcept
Gives the lowest M for a polarization type of this transition.
Definition: zeemandata.h:80
auto mat(matpack::matrix auto &&x)
Map the input to a non-owning const-correct Eigen Map representing a row matrix.
Definition: matpack_eigen.h:76
invlib::Matrix< ArtsMatrix > Matrix
invlib wrapper type for ARTS matrices.
Definition: oem.h:33
invlib::Vector< ArtsVector > Vector
invlib wrapper type for ARTS vectors.
Definition: oem.h:31
Contains the rational class definition.
Numeric sqrt(const Rational r)
Square root.
Definition: rational.h:705
constexpr Rational numeric2rational(Numeric x, size_t maxdec=4) noexcept
Rational from Numeric.
Definition: rational.h:329
#define M
Definition: rng.cc:165
Implements rational numbers to work with other ARTS types.
Definition: rational.h:43
constexpr bool isUndefined() const noexcept
Is the object not defined.
Definition: rational.h:77
constexpr Numeric toNumeric() const noexcept
Converts this to a Numeric.
Definition: rational.h:111
constexpr Index toIndex(int n=1) const noexcept
Converts the value to index by n-scaled division.
Definition: rational.h:103
Index numer
Definition: rational.h:44
Index denom
Definition: rational.h:45
constexpr int toInt(int n=1) const noexcept
Converts the value to int by n-scaled division in Index form.
Definition: rational.h:122
void test26()
void test2()
void test29()
void nlinspace(Vector &x, const Numeric start, const Numeric stop, const Index n)
nlinspace
void test30()
void test46()
void junk2(ConstVectorView a)
void test18()
void test_wigner_error()
#define docheck(fn, val, expect)
void test8()
void test4()
Numeric matrix_vector_mult(Index m, Index n, Index ntests, bool verbose)
void junk4(Tensor4View a)
void test9()
void test42()
void test20()
void test27()
Numeric test_matrix_vector_multiplication(bool verbose)
void test15()
void test48()
int test1()
Definition: test_matpack.cc:49
void test_pow_negative_one()
Numeric matrix_mult(Index k, Index m, Index n, Index ntests, Index nsubtests, bool verbose)
Perform matrix multiplication test.
void test25()
Numeric test_matrix_multiplication(bool verbose)
Perform matrix multiplication tests.
void test31()
Numeric by_value(Numeric x)
Definition: test_matpack.cc:43
void test44()
void test39()
void test23()
void fill_with_junk(VectorView x)
Definition: test_matpack.cc:45
void test5()
void test7()
void test6()
Numeric by_reference(const Numeric &x)
Definition: test_matpack.cc:41
void test33()
void test10()
void test24()
void test47()
void test12()
void test19()
constexpr Rational test_numeric2rational(const Index i, const Index maxi, const Rational r=0, const Rational n2r=0)
void test17()
void test13()
void test38()
void test35()
void test40()
void test14()
bool test_diagonal(Index ntests)
Test diagonal vector.
void test37(const Index &i)
void test11()
int main()
void test28()
void test43()
void test41()
void test22()
void test34()
void test21()
void test_empty()
Check if the empty function is working correctly.
void test32()
void test36()
void random_fill_matrix(MatrixView A, Numeric range, bool positive)
Fill matrix with random values.
Definition: test_utils.cc:59
Numeric get_maximum_error(ConstVectorView v1, ConstVectorView v2, bool relative)
Maximum element-wise error of two vectors.
Definition: test_utils.cc:297
void random_fill_vector(VectorView v, Numeric range, bool positive)
Fill vector with random values.
Definition: test_utils.cc:230
Range random_range(Index n)
Generate random sub-range of the range [0, n-1].
Definition: test_utils.cc:273
Utility functions for testing.
#define d
#define v
#define a
#define c
#define b
Index make_wigner_ready(int largest, int fastest, int size)
Ready Wigner.
Numeric wigner3j(const Rational j1, const Rational j2, const Rational j3, const Rational m1, const Rational m2, const Rational m3)
Wigner 3J symbol.
Wigner symbol interactions.