ARTS  2.4.0(git:4fb77825)
test_covariance_matrix.cc
Go to the documentation of this file.
1 /* Copyright (C) 2017
2  Simon Pfreundschuh <simon.pfreundschuh@chalmers.se>
3 
4  This program is free software; you can redistribute it and/or modify it
5  under the terms of the GNU General Public License as published by the
6  Free Software Foundation; either version 2, or (at your option) any
7  later version.
8 
9  This program is distributed in the hope that it will be useful,
10  but WITHOUT ANY WARRANTY; without even the implied warranty of
11  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12  GNU General Public License for more details.
13 
14  You should have received a copy of the GNU General Public License
15  along with this program; if not, write to the Free Software
16  Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307,
17  USA. */
18 
27 #include <cstdlib>
28 #include <limits>
29 #include <random>
30 #include <utility>
31 
32 #include "array.h"
33 #include "covariance_matrix.h"
34 #include "jacobian.h"
35 #include "lin_alg.h"
36 #include "test_utils.h"
37 #include "xml_io.h"
38 
39 // Type alias describing a retrieval setup for testing.
40 using RetrievalData = std::tuple<ArrayOfRetrievalQuantity, ArrayOfArrayOfIndex>;
41 
54 template <typename F>
55 std::shared_ptr<Matrix> create_covariance_matrix_1D(
56  Index i, Index j, ArrayOfRetrievalQuantity jqs, F f) {
57  if (i > j) {
58  std::swap(i, j);
59  }
60 
61  const RetrievalQuantity& q1 = jqs[i];
62  const RetrievalQuantity& q2 = jqs[j];
63 
64  const Vector& gv1 = q1.Grids()[0];
65  const Vector& gv2 = q2.Grids()[0];
66 
67  std::shared_ptr<Matrix> c =
68  std::make_shared<Matrix>(gv1.nelem(), gv2.nelem());
69 
70  for (Index k = 0; k < gv1.nelem(); k++) {
71  for (Index l = 0; l < gv2.nelem(); l++) {
72  (*c)(k, l) = f(gv1[k], gv2[l]);
73  }
74  }
75  return c;
76 }
77 
81 template <typename F>
82 std::shared_ptr<Sparse> create_sparse_covariance_matrix_1D(
83  Index i, Index j, ArrayOfRetrievalQuantity jqs, F f) {
84  if (i > j) {
85  std::swap(i, j);
86  }
87 
88  const RetrievalQuantity& q1 = jqs[i];
89  const RetrievalQuantity& q2 = jqs[j];
90 
91  const Vector& gv1 = q1.Grids()[0];
92  const Vector& gv2 = q2.Grids()[0];
93 
94  std::shared_ptr<Sparse> s =
95  std::make_shared<Sparse>(gv1.nelem(), gv2.nelem());
96 
97  Index nelem = gv1.nelem() * gv2.nelem();
98  ArrayOfIndex row_indices{}, col_indices{};
99  ArrayOfNumeric elements{};
100  row_indices.reserve(nelem);
101  col_indices.reserve(nelem);
102  elements.reserve(nelem);
103  for (Index k = 0; k < gv1.nelem(); k++) {
104  for (Index l = 0; l < gv2.nelem(); l++) {
105  row_indices.push_back(k);
106  col_indices.push_back(l);
107  elements.push_back(f(gv1[k], gv2[l]));
108  }
109  }
110  s->insert_elements(nelem, row_indices, col_indices, Vector(elements));
111  return s;
112 }
113 
121  std::random_device rd;
122  std::mt19937 gen(rd());
123  std::uniform_int_distribution<> n_rqs_dist(2, 10), n_gs_dist(10, 100);
124 
125  Index n_rqs = n_rqs_dist(gen);
126  Numeric g_begin = 0.0;
127  Numeric g_end = 100.0;
129 
130  for (Index i = 0; i < n_rqs; i++) {
131  Index n_gs = n_gs_dist(gen);
132  ArrayOfVector gvs(1);
133  gvs[0] = Vector(n_gs);
134  for (Index j = 0; j < n_gs; j++) {
135  Numeric frac = static_cast<Numeric>(j) / static_cast<Numeric>(n_gs);
136  gvs[0][j] = g_begin + (g_end - g_begin) * frac;
137  }
138  rqs.emplace_back("maintag", "subtag", "subsubtag", "mode", 1, 0.0, gvs);
139  }
140 
141  ArrayOfArrayOfIndex jis(n_rqs);
142  Index ji = 0;
143  for (Index i = 0; i < n_rqs; i++) {
144  Index n_grid_points = 1;
145  const ArrayOfVector& gs = rqs[i].Grids();
146  for (auto& g : gs) {
147  n_grid_points *= g.nelem();
148  }
149  jis[i] = ArrayOfIndex{ji, ji + n_grid_points - 1};
150  ji += n_grid_points;
151  }
152 
153  return std::make_tuple(rqs, jis);
154 }
155 
164  const ArrayOfArrayOfIndex& jis) {
165  std::random_device rd;
166  std::mt19937 gen(rd());
167  std::uniform_real_distribution<> dis(0.0, 1.0);
168  CovarianceMatrix covmat{};
169 
170  auto g = [](Numeric x, Numeric y) -> Numeric {
171  if (x == y)
172  return 1.0;
173  else
174  return 0.0;
175  };
176  auto g2 = [](Numeric x, Numeric y) -> Numeric {
177  if (x == y)
178  return 0.1;
179  else
180  return 0.0;
181  };
182 
183  for (size_t i = 0; i < rqs.size(); i++) {
184  Range range(jis[i][0], jis[i][1] - jis[i][0] + 1);
185  auto inds = std::make_pair(i, i);
186 
187  Numeric p = dis(gen);
188  if (p < 0.0) {
189  std::shared_ptr<Matrix> m = create_covariance_matrix_1D(i, i, rqs, g);
190  covmat.add_correlation(Block(range, range, inds, m));
191  } else {
192  std::shared_ptr<Sparse> m =
194  covmat.add_correlation(Block(range, range, inds, m));
195  }
196  }
197 
198  for (size_t i = 0; i < rqs.size(); i++) {
199  for (size_t j = i + 1; j < rqs.size(); j++) {
200  Range row_range(jis[i][0], jis[i][1] - jis[i][0] + 1);
201  Range col_range(jis[j][0], jis[j][1] - jis[j][0] + 1);
202  auto inds = std::make_pair(i, j);
203  Numeric p = dis(gen);
204  if (p > 0.8) {
205  p = dis(gen);
206  if (p < 0.0) {
207  std::shared_ptr<Matrix> m = create_covariance_matrix_1D(i, j, rqs, g);
208  covmat.add_correlation(Block(row_range, col_range, inds, m));
209  } else {
210  std::shared_ptr<Sparse> m =
211  create_sparse_covariance_matrix_1D(i, j, rqs, g2);
212  covmat.add_correlation(Block(row_range, col_range, inds, m));
213  }
214  }
215  }
216  }
217  return covmat;
218 }
219 
228  Numeric e = 0.0;
229  for (Index i = 0; i < n_tests; i++) {
232  std::tie(rqs, jis) = setup_retrieval_1D();
233  CovarianceMatrix covmat(random_covariance_matrix(rqs, jis));
234 
235  Index n = covmat.ncols();
236  Matrix A(covmat);
237  Vector v(n), w(n), w_ref(n);
238  random_fill_vector(v, 10.0, false);
239 
240  mult(w, covmat, v);
241  mult(w_ref, A, v);
242 
243  e = std::max(e, get_maximum_error(w, w_ref, true));
244  }
245  return e;
246 }
247 
256  Numeric e = 0.0;
257  for (Index i = 0; i < n_tests; i++) {
260  std::tie(rqs, jis) = setup_retrieval_1D();
261  CovarianceMatrix covmat(random_covariance_matrix(rqs, jis));
262 
263  Index n = covmat.ncols();
264  Matrix A(covmat), B(n, n), C(n, n), C_ref(n, n);
265  random_fill_matrix(B, 10.0, false);
266 
267  mult(C, covmat, B);
268  mult(C_ref, A, B);
269  e = std::max(e, get_maximum_error(C, C_ref, true));
270 
271  mult(C, B, covmat);
272  mult(C_ref, B, A);
273  e = std::max(e, get_maximum_error(C, C_ref, true));
274  }
275  return e;
276 }
277 
289  Numeric e = 0.0;
290  for (Index i = 0; i < n_tests; i++) {
293  std::tie(rqs, jis) = setup_retrieval_1D();
294  CovarianceMatrix covmat(random_covariance_matrix(rqs, jis));
295 
296  Index n = covmat.ncols();
297  Matrix A(covmat), B(n, n), B_ref(n, n), C(n, n);
298  random_fill_matrix(B, 10.0, false);
299 
300  covmat.compute_inverse();
301  mult_inv(B, covmat, A);
302  id_mat(B_ref);
303  mult(C, covmat, B_ref);
304  e = std::max(e, get_maximum_error(B, B_ref, true));
305  id_mat(B_ref);
306  mult_inv(B, A, covmat);
307  e = std::max(e, get_maximum_error(B, B_ref, true));
308  }
309  return e;
310 }
311 
320  Numeric e = 0.0;
321  for (Index i = 0; i < n_tests; i++) {
324  std::tie(rqs, jis) = setup_retrieval_1D();
325  CovarianceMatrix covmat(random_covariance_matrix(rqs, jis));
326 
327  Index n = covmat.ncols();
328  Matrix A(covmat), B(n, n), B_ref(n, n), C(n, n);
329  B = 0.0;
330  B_ref = B;
331 
332  covmat.compute_inverse();
333 
334  B += covmat;
335  B_ref += A;
336  e = std::max(e, get_maximum_error(B, B_ref, true));
337 
338  add_inv(B, covmat);
339  inv(C, A);
340  B_ref += C;
341  e = std::max(e, get_maximum_error(B, B_ref, true));
342  }
343  return e;
344 }
345 
355  Numeric e = 0.0;
356  for (Index i = 0; i < n_tests; i++) {
359  std::tie(rqs, jis) = setup_retrieval_1D();
360  CovarianceMatrix covmat(random_covariance_matrix(rqs, jis));
361 
362  Vector diag_1 = covmat.diagonal();
363  Vector diag_2 = Matrix(covmat).diagonal();
364  e = std::max(e, get_maximum_error(diag_1, diag_2, true));
365  }
366  return e;
367 }
368 
377  Numeric e(0.0);
378  for (Index i = 0; i < n_tests; i++) {
381  std::tie(rqs, jis) = setup_retrieval_1D();
382  CovarianceMatrix covmat_1(random_covariance_matrix(rqs, jis)), covmat_2{};
383  covmat_1.compute_inverse();
384 
385  xml_write_to_file("test.xml", covmat_1, FILE_TYPE_ASCII, 0, Verbosity());
386  xml_read_from_file("test.xml", covmat_2, Verbosity());
387 
388  Matrix A(covmat_1), B(covmat_2), C(covmat_1);
389  e = std::max(e, get_maximum_error(A, B, true));
390 
391  id_mat(C);
392  mult_inv(A, covmat_1, C);
393  mult_inv(B, covmat_2, C);
394  e = std::max(e, get_maximum_error(A, B, true));
395  }
396  return 0.0;
397 }
398 
399 template <typename MatrixType>
400 void covmat_seSet(CovarianceMatrix& covmat,
401  const MatrixType& block,
402  const Verbosity& /*v*/);
403 
404 template <typename MatrixType>
405 void covmat_sxSet(CovarianceMatrix& covmat,
406  const MatrixType& block,
407  const Verbosity& /*v*/);
408 
409 template <typename MatrixType>
411  const MatrixType& block,
412  const Index& i,
413  const Index& j,
414  const Verbosity&);
415 
416 template <typename MatrixType>
418  const MatrixType& block,
419  const Index& i,
420  const Index& j,
421  const Verbosity&);
422 
423 template <typename MatrixType>
425  const ArrayOfRetrievalQuantity& jq,
426  const MatrixType& block,
427  const Index& i,
428  const Index& j,
429  const Verbosity&);
430 
431 template <typename MatrixType>
433  const ArrayOfRetrievalQuantity& jq,
434  const MatrixType& block,
435  const Index& i,
436  const Index& j,
437  const Verbosity&);
438 
443  CovarianceMatrix covmat{};
444  Matrix A(10, 10);
445  Sparse A_sparse(10, 10);
446  Matrix B(20, 20);
447  Matrix C(10, 20);
448 
449  // covmat_seSet
450 
451  covmat_seSet(covmat, A, Verbosity());
452  assert(covmat.ncols() == 10);
453 
454  // covmat_seAddBlock
455 
456  covmat_seAddBlock(covmat, A_sparse, -1, -1, Verbosity());
457  assert(covmat.ncols() == 20);
458 
459  try {
460  covmat_seAddBlock(covmat, A, 3, 3, Verbosity());
461  // This should fail.
462  assert(false);
463  } catch (const std::runtime_error&) {
464  }
465 
466  covmat_seAddBlock(covmat, A, 0, 1, Verbosity());
467 
468  covmat_seAddBlock(covmat, A, 2, 2, Verbosity());
469 
470  try {
471  covmat_seAddBlock(covmat, B, 1, 2, Verbosity());
472  // This should fail.
473  assert(false);
474  } catch (const std::runtime_error&) {
475  }
476 
477  covmat_seAddInverseBlock(covmat, A, 0, 1, Verbosity());
478 
479  try {
480  covmat_seAddInverseBlock(covmat, B, 3, 3, Verbosity());
481  // This should fail.
482  assert(false);
483  } catch (const std::runtime_error&) {
484  }
485 
486  // covmat_sxSet
487  covmat_sxSet(covmat, A, Verbosity());
488  assert(covmat.ncols() == 10);
489 
490  // covmat_sxAddBlock
491 
492  A.resize(1000, 1000);
493  A_sparse.resize(8000, 8000);
494  C.resize(1000, 8000);
495 
496  ArrayOfVector grids_1{Vector(10), Vector(10), Vector(10)};
497  ArrayOfVector grids_2{Vector(20), Vector(20), Vector(20)};
498 
499  RetrievalQuantity rq_1("mt", "st", "sst", "m", 1, 0.1, grids_1);
500  RetrievalQuantity rq_2("mt", "st", "sst", "m", 1, 0.1, grids_2);
501 
503  rqs.push_back(rq_1);
504  rqs.push_back(rq_2);
505 
506  covmat_sxAddBlock(covmat, rqs, A_sparse, -1, -1, Verbosity());
507 
508  try {
509  covmat_sxAddBlock(covmat, rqs, A_sparse, 0, 1, Verbosity());
510  // This should fail.
511  assert(false);
512  } catch (const std::runtime_error&) {
513  }
514 
515  covmat_sxAddBlock(covmat, rqs, C, 0, 1, Verbosity());
516  covmat_sxAddInverseBlock(covmat, rqs, A, 0, 0, Verbosity());
517 
518  try {
519  covmat_sxAddInverseBlock(covmat, rqs, C, 1, 1, Verbosity());
520  // This should fail.
521  assert(false);
522  } catch (const std::runtime_error&) {
523  }
524 }
525 
533 namespace {
534 #pragma GCC diagnostic push
535 #pragma GCC diagnostic ignored "-Wunused-function"
536 #pragma GCC diagnostic ignored "-Wconversion"
537 
538 #include "invlib/algebra.h"
539 #include "invlib/interfaces/arts_wrapper.h"
540 
541 Numeric test_invlib_wrapper(Index n_tests) {
542  Numeric e = 0.0;
543  for (Index i = 0; i < n_tests; i++) {
546  std::tie(rqs, jis) = setup_retrieval_1D();
547  CovarianceMatrix covmat(random_covariance_matrix(rqs, jis));
548  invlib::Matrix<ArtsCovarianceMatrixWrapper> wrapper(covmat);
549  Index n = covmat.ncols();
550 
551  // Multiplication by Vector
552  invlib::Vector<ArtsVector> v{}, w{}, w_ref{};
553  v.resize(n);
554  w.resize(n);
555  w_ref.resize(n);
556  random_fill_vector(v, 10.0, false);
557 
558  mult(w, covmat, v);
559  w_ref = wrapper * v;
560  e = std::max(e, get_maximum_error(w, w_ref, true));
561 
562  solve(w, covmat, v);
563  w_ref = inv(wrapper) * v;
564  e = std::max(e, get_maximum_error(w, w_ref, true));
565  w_ref = inv(inv(wrapper)) * v;
566 
567  // Multiplication by Matrix
568  invlib::Matrix<ArtsMatrix> A{static_cast<Matrix>(covmat)}, B{}, C{},
569  C_ref{};
570  B.resize(n, n);
571  C.resize(n, n);
572  C_ref.resize(n, n);
573 
574  random_fill_matrix(B, 10.0, false);
575  mult(C, covmat, B);
576  C_ref = wrapper * B;
577  e = std::max(e, get_maximum_error(C, C_ref, true));
578 
579  // Inverse
580  mult_inv(C, covmat, B);
581  C_ref = invlib::inv(wrapper) * B;
582  e = std::max(e, get_maximum_error(C, C_ref, true));
583 
584  // Addition to Matrix
585  C += covmat;
586  C_ref = C_ref + wrapper;
587  e = std::max(e, get_maximum_error(C, C_ref, true));
588 
589  // Addition of Inverse to Matrix
590  add_inv(C, covmat);
591  C_ref = C_ref + inv(wrapper);
592  e = std::max(e, get_maximum_error(C, C_ref, true));
593  }
594  return e;
595 }
596 
597 #pragma GCC diagnostic pop
598 } // namespace
599 
600 int main() {
601  Numeric e = 0.0, e_max = 0.0;
602 
603  std::cout << "Testing covariance matrix: " << std::endl;
604  e = test_addition(10);
605  std::cout << "\tAddition: " << e << std::endl;
606  e_max = std::max(e, e_max);
607  if (e_max > 1e-5) {
608  return -1;
609  }
610 
612  std::cout << "\tVector Multiplication: " << e << std::endl;
613  e_max = std::max(e, e_max);
614  if (e_max > 1e-5) {
615  return -1;
616  }
617 
619  std::cout << "\tMatrix Multiplication: " << e << std::endl;
620  e_max = std::max(e, e_max);
621  if (e_max > 1e-5) {
622  return -1;
623  }
624 
625  e = test_inverse(10);
626  std::cout << "\tInverse: " << e << std::endl;
627  e_max = std::max(e, e_max);
628  if (e_max > 1e-5) {
629  return -1;
630  }
631 
632  e = test_io(10);
633  std::cout << "\tXML IO: " << e << std::endl;
634  e_max = std::max(e, e_max);
635  if (e_max > 1e-5) {
636  return -1;
637  }
638 
639  e = test_invlib_wrapper(10);
640  std::cout << "\tinvlib Wrapper: " << e << std::endl;
641  e_max = std::max(e, e_max);
642  if (e_max > 1e-5) {
643  return -1;
644  }
645 
646  e = test_diagonal(10);
647  std::cout << "\tdiagonal " << e << std::endl;
648  e_max = std::max(e, e_max);
649  if (e_max > 1e-5) {
650  return -1;
651  }
652 
653  std::cout << std::endl << "\tTesting workspace functions ... ";
655  std::cout << " DONE." << std::endl;
656 
657  return 0;
658 }
Matrix
The Matrix class.
Definition: matpackI.h:1193
CovarianceMatrix::ncols
Index ncols() const
Definition: covariance_matrix.cc:222
FILE_TYPE_ASCII
@ FILE_TYPE_ASCII
Definition: xml_io.h:38
test_multiplication_by_matrix
Numeric test_multiplication_by_matrix(Index n_tests)
Tests the multiplication of covariance matrices with matrices.
Definition: test_covariance_matrix.cc:255
test_multiplication_by_vector
Numeric test_multiplication_by_vector(Index n_tests)
Tests the multiplication of covariance matrices with vectors.
Definition: test_covariance_matrix.cc:227
w
Complex w(Complex z) noexcept
The Faddeeva function.
Definition: linefunctions.cc:42
RetrievalQuantity::Grids
const ArrayOfVector & Grids() const
Returns the grids of the retrieval.
Definition: jacobian.h:255
test_utils.h
Utility functions for testing.
test_diagonal
Numeric test_diagonal(Index n_tests)
Test extraction of diagonal.
Definition: test_covariance_matrix.cc:354
id_mat
void id_mat(MatrixView I)
Identity Matrix.
Definition: lin_alg.cc:2240
test_workspace_methods
void test_workspace_methods()
Test general functionality of CovarianceMatrix class.
Definition: test_covariance_matrix.cc:442
inv
void inv(MatrixView Ainv, ConstMatrixView A)
Matrix Inverse.
Definition: lin_alg.cc:167
oem::Matrix
invlib::Matrix< ArtsMatrix > Matrix
invlib wrapper type for ARTS matrices.
Definition: oem.h:34
covmat_seSet
void covmat_seSet(CovarianceMatrix &covmat, const MatrixType &block, const Verbosity &)
Definition: m_retrieval.cc:561
mult
void mult(ComplexVectorView y, const ConstComplexMatrixView &M, const ConstComplexVectorView &x)
Matrix-Vector Multiplication.
Definition: complex.cc:1579
main
int main()
Definition: test_covariance_matrix.cc:600
random_covariance_matrix
CovarianceMatrix random_covariance_matrix(const ArrayOfRetrievalQuantity &rqs, const ArrayOfArrayOfIndex &jis)
Create a covariance matrix with correlations between retrieval quantitites given in rqs and jis.
Definition: test_covariance_matrix.cc:163
swap
void swap(ComplexVector &v1, ComplexVector &v2)
Swaps two objects.
Definition: complex.cc:731
Sparse
The Sparse class.
Definition: matpackII.h:60
random_fill_vector
void random_fill_vector(VectorView v, Numeric range, bool positive)
Fill vector with random values.
Definition: test_utils.cc:230
CovarianceMatrix
Definition: covariance_matrix.h:226
ARTS::Var::y
Vector y(Workspace &ws) noexcept
Definition: autoarts.h:7401
ArrayOfRetrievalQuantity
Array< RetrievalQuantity > ArrayOfRetrievalQuantity
Definition: jacobian.h:402
CovarianceMatrix::diagonal
Vector diagonal() const
Diagonal elements as vector.
Definition: covariance_matrix.cc:516
covmat_sxAddInverseBlock
void covmat_sxAddInverseBlock(CovarianceMatrix &covmat_se, const ArrayOfRetrievalQuantity &jq, const MatrixType &block, const Index &i, const Index &j, const Verbosity &)
Definition: m_retrieval.cc:678
covariance_matrix.h
Header files of CovarianceMatrix class.
array.h
This file contains the definition of Array.
CovarianceMatrix::compute_inverse
void compute_inverse() const
Compute the inverse of this correlation matrix.
Definition: covariance_matrix.cc:391
covmat_sxAddBlock
void covmat_sxAddBlock(CovarianceMatrix &covmat_se, const ArrayOfRetrievalQuantity &jq, const MatrixType &block, const Index &i, const Index &j, const Verbosity &)
Definition: m_retrieval.cc:599
Array
This can be used to make arrays out of anything.
Definition: array.h:108
test_io
Numeric test_io(Index n_tests)
Test input and output of covariance matrices.
Definition: test_covariance_matrix.cc:376
Absorption::nelem
Index nelem(const Lines &l)
Number of lines.
Definition: absorptionlines.h:1820
random_fill_matrix
void random_fill_matrix(MatrixView A, Numeric range, bool positive)
Fill matrix with random values.
Definition: test_utils.cc:59
xml_read_from_file
void xml_read_from_file(const String &filename, T &type, const Verbosity &verbosity)
Reads data from XML file.
Definition: xml_io.cc:901
test_addition
Numeric test_addition(Index n_tests)
Test addition of covariance matrices and inverse covariance matrices.
Definition: test_covariance_matrix.cc:319
ConstVectorView::nelem
Index nelem() const
Returns the number of elements.
Definition: matpackI.cc:51
q1
#define q1
Definition: legacy_continua.cc:21438
create_covariance_matrix_1D
std::shared_ptr< Matrix > create_covariance_matrix_1D(Index i, Index j, ArrayOfRetrievalQuantity jqs, F f)
Create a block of a covariance matrix representating a correlation between the retrieval quantities i...
Definition: test_covariance_matrix.cc:55
jacobian.h
Routines for setting up the jacobian.
Numeric
NUMERIC Numeric
The type to use for all floating point numbers.
Definition: matpack.h:33
lin_alg.h
Linear algebra functions.
Verbosity
Definition: messages.h:49
ARTS::Var::covmat_se
CovarianceMatrix covmat_se(Workspace &ws) noexcept
Definition: autoarts.h:2879
RetrievalData
std::tuple< ArrayOfRetrievalQuantity, ArrayOfArrayOfIndex > RetrievalData
Definition: test_covariance_matrix.cc:40
max
#define max(a, b)
Definition: legacy_continua.cc:20629
MatrixType
Eigen::Matrix< Numeric, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > MatrixType
Definition: matpackI.h:109
Matrix::resize
void resize(Index r, Index c)
Resize function.
Definition: matpackI.cc:1056
mult_inv
void mult_inv(MatrixView C, ConstMatrixView A, const CovarianceMatrix &B)
Definition: covariance_matrix.cc:574
Sparse::resize
void resize(Index r, Index c)
Resize function.
Definition: matpackII.cc:361
oem::Vector
invlib::Vector< ArtsVector > Vector
invlib wrapper type for ARTS vectors.
Definition: oem.h:32
covmat_sxSet
void covmat_sxSet(CovarianceMatrix &covmat, const MatrixType &block, const Verbosity &)
Definition: m_retrieval.cc:580
Range
The range class.
Definition: matpackI.h:160
xml_write_to_file
void xml_write_to_file(const String &filename, const T &type, const FileType ftype, const Index no_clobber, const Verbosity &verbosity)
Write data to XML file.
Definition: xml_io.cc:972
covmat_seAddBlock
void covmat_seAddBlock(CovarianceMatrix &covmat_se, const MatrixType &block, const Index &i, const Index &j, const Verbosity &)
Definition: m_retrieval.cc:397
setup_retrieval_1D
RetrievalData setup_retrieval_1D()
Setup a random retrieval case for testing.
Definition: test_covariance_matrix.cc:120
add_inv
void add_inv(MatrixView A, const CovarianceMatrix &B)
Definition: covariance_matrix.cc:611
Block
Definition: covariance_matrix.h:62
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_inverse
Numeric test_inverse(Index n_tests)
Tests the inversion of covariance matrices by computing the products of inverses of covariance matric...
Definition: test_covariance_matrix.cc:288
ARTS::Var::x
Vector x(Workspace &ws) noexcept
Definition: autoarts.h:7346
Index
INDEX Index
The type to use for all integer numbers and indices.
Definition: matpack.h:39
RetrievalQuantity
Deals with internal derivatives, Jacobian definition, and OEM calculations.
Definition: jacobian.h:120
Vector
The Vector class.
Definition: matpackI.h:860
Array::nelem
Index nelem() const
Number of elements.
Definition: array.h:195
ARTS::Group::Verbosity
Verbosity Verbosity
Definition: autoarts.h:114
covmat_seAddInverseBlock
void covmat_seAddInverseBlock(CovarianceMatrix &covmat_se, const MatrixType &block, const Index &i, const Index &j, const Verbosity &)
Definition: m_retrieval.cc:490
create_sparse_covariance_matrix_1D
std::shared_ptr< Sparse > create_sparse_covariance_matrix_1D(Index i, Index j, ArrayOfRetrievalQuantity jqs, F f)
Same as create_covariance_matrix_1D but creates a matrix of type sparse.
Definition: test_covariance_matrix.cc:82
xml_io.h
This file contains basic functions to handle XML data files.
solve
void solve(VectorView w, const CovarianceMatrix &A, ConstVectorView v)
Definition: covariance_matrix.cc:594