Visual Servoing Platform  version 3.6.1 under development (2024-10-15)
testMatrixInverse.cpp
1 /*
2  * ViSP, open source Visual Servoing Platform software.
3  * Copyright (C) 2005 - 2024 by Inria. All rights reserved.
4  *
5  * This software is free software; you can redistribute it and/or modify
6  * it under the terms of the GNU General Public License as published by
7  * the Free Software Foundation; either version 2 of the License, or
8  * (at your option) any later version.
9  * See the file LICENSE.txt at the root directory of this source
10  * distribution for additional information about the GNU GPL.
11  *
12  * For using ViSP with software that can not be combined with the GNU
13  * GPL, please contact Inria about acquiring a ViSP Professional
14  * Edition License.
15  *
16  * See https://visp.inria.fr for more information.
17  *
18  * This software was developed at:
19  * Inria Rennes - Bretagne Atlantique
20  * Campus Universitaire de Beaulieu
21  * 35042 Rennes Cedex
22  * France
23  *
24  * If you have questions regarding the use of this file, please contact
25  * Inria at visp@inria.fr
26  *
27  * This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
28  * WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
29  *
30  * Description:
31  * Test various inversions.
32  */
33 
39 #include <cmath>
40 #include <fstream>
41 #include <iostream>
42 #include <stdio.h>
43 #include <stdlib.h>
44 #include <vector>
45 
46 #include <visp3/core/vpColVector.h>
47 #include <visp3/core/vpMatrix.h>
48 #include <visp3/core/vpTime.h>
49 #include <visp3/io/vpParseArgv.h>
50 
51 // List of allowed command line options
52 #define GETOPTARGS "cdn:i:pf:R:C:vh"
53 
54 #ifdef ENABLE_VISP_NAMESPACE
55 using namespace VISP_NAMESPACE_NAME;
56 #endif
57 
66 void usage(const char *name, const char *badparam)
67 {
68  fprintf(stdout, "\n\
69 Test matrix inversions\n\
70 using LU, QR and Cholesky methods as well as Pseudo-inverse.\n\
71 Outputs a comparison of these methods.\n\
72 \n\
73 SYNOPSIS\n\
74  %s [-n <number of matrices>] [-f <plot filename>]\n\
75  [-R <number of rows>] [-C <number of columns>]\n\
76  [-i <number of iterations>] [-p] [-h]\n",
77  name);
78 
79  fprintf(stdout, "\n\
80 OPTIONS: Default\n\
81  -n <number of matrices> \n\
82  Number of matrices inverted during each test loop.\n\
83 \n\
84  -i <number of iterations> \n\
85  Number of iterations of the test.\n\
86 \n\
87  -f <plot filename> \n\
88  Set output path for plot output.\n\
89  The plot logs the times of \n\
90  the different inversion methods: \n\
91  QR,LU,Cholesky and Pseudo-inverse.\n\
92 \n\
93  -R <number of rows>\n\
94  Number of rows of the automatically generated matrices \n\
95  we test on.\n\
96 \n\
97  -C <number of columns>\n\
98  Number of colums of the automatically generated matrices \n\
99  we test on.\n\
100 \n\
101  -p \n\
102  Plot into filename in the gnuplot format. \n\
103  If this option is used, tests results will be logged \n\
104  into a filename specified with -f.\n\
105 \n\
106  -h\n\
107  Print the help.\n\n");
108 
109  if (badparam) {
110  fprintf(stderr, "ERROR: \n");
111  fprintf(stderr, "\nBad parameter [%s]\n", badparam);
112  }
113 }
114 
122 bool getOptions(int argc, const char **argv, unsigned int &nb_matrices, unsigned int &nb_iterations,
123  bool &use_plot_file, std::string &plotfile, unsigned int &nbrows, unsigned int &nbcols, bool &verbose)
124 {
125  const char *optarg_;
126  int c;
127  while ((c = vpParseArgv::parse(argc, argv, GETOPTARGS, &optarg_)) > 1) {
128 
129  switch (c) {
130  case 'h':
131  usage(argv[0], nullptr);
132  return false;
133  break;
134  case 'n':
135  nb_matrices = (unsigned int)atoi(optarg_);
136  break;
137  case 'i':
138  nb_iterations = (unsigned int)atoi(optarg_);
139  break;
140  case 'f':
141  plotfile = optarg_;
142  use_plot_file = true;
143  break;
144  case 'p':
145  use_plot_file = true;
146  break;
147  case 'R':
148  nbrows = (unsigned int)atoi(optarg_);
149  break;
150  case 'C':
151  nbcols = (unsigned int)atoi(optarg_);
152  break;
153  case 'v':
154  verbose = true;
155  break;
156  // add default options -c -d
157  case 'c':
158  break;
159  case 'd':
160  break;
161  default:
162  usage(argv[0], optarg_);
163  return false;
164  break;
165  }
166  }
167 
168  if ((c == 1) || (c == -1)) {
169  // standalone param or error
170  usage(argv[0], nullptr);
171  std::cerr << "ERROR: " << std::endl;
172  std::cerr << " Bad argument " << optarg_ << std::endl << std::endl;
173  return false;
174  }
175 
176  return true;
177 }
178 
179 vpMatrix make_random_matrix(unsigned int nbrows, unsigned int nbcols)
180 {
181  vpMatrix A;
182  A.resize(nbrows, nbcols);
183 
184  for (unsigned int i = 0; i < A.getRows(); i++)
185  for (unsigned int j = 0; j < A.getCols(); j++)
186  A[i][j] = (double)rand() / (double)RAND_MAX;
187  return A;
188 }
189 
190 vpMatrix make_random_symmetric_positive_matrix(unsigned int n)
191 {
192  vpMatrix A;
193  A.resize(n, n);
194  vpMatrix I;
195  I.eye(n);
196 
197  for (unsigned int i = 0; i < A.getRows(); i++)
198  for (unsigned int j = 0; j < A.getCols(); j++)
199  A[i][j] = (double)rand() / (double)RAND_MAX;
200 
201  A = 0.5 * (A + A.t());
202  A = A + n * I;
203  return A;
204 }
205 
206 vpMatrix make_random_triangular_matrix(unsigned int nbrows)
207 {
208  vpMatrix A;
209  A.resize(nbrows, nbrows);
210 
211  for (unsigned int i = 0; i < A.getRows(); i++) {
212  for (unsigned int j = i; j < A.getCols(); j++) {
213  A[i][j] = static_cast<double>(rand()) / static_cast<double>(RAND_MAX);
214  if (i != j) {
215  A[j][i] = 0;
216  }
217  }
218  }
219 
220  return A;
221 }
222 
223 void create_bench_random_matrix(unsigned int nb_matrices, unsigned int nb_rows, unsigned int nb_cols, bool verbose,
224  std::vector<vpMatrix> &bench)
225 {
226  if (verbose)
227  std::cout << "Create a bench of " << nb_matrices << " " << nb_rows << " by " << nb_cols << " matrices" << std::endl;
228  bench.clear();
229  for (unsigned int i = 0; i < nb_matrices; i++) {
230  vpMatrix M;
231 #if defined(VISP_HAVE_EIGEN3) || defined(VISP_HAVE_LAPACK) || defined(VISP_HAVE_OPENCV)
232  double det = 0.;
233  // don't put singular matrices in the benchmark
234  for (M = make_random_matrix(nb_rows, nb_cols); std::fabs(det = M.AtA().det()) < .01;
235  M = make_random_matrix(nb_rows, nb_cols)) {
236  if (verbose) {
237  std::cout << " Generated random matrix AtA=" << std::endl << M.AtA() << std::endl;
238  std::cout << " Generated random matrix not invertible: det=" << det << ". Retrying..." << std::endl;
239  }
240  }
241 #else
242  M = make_random_matrix(nb_rows, nb_cols);
243 #endif
244  bench.push_back(M);
245  }
246 }
247 
248 void create_bench_symmetric_positive_matrix(unsigned int nb_matrices, unsigned int n, bool verbose,
249  std::vector<vpMatrix> &bench)
250 {
251  if (verbose)
252  std::cout << "Create a bench of " << nb_matrices << " " << n << " by " << n << " symmetric positive matrices"
253  << std::endl;
254  bench.clear();
255  for (unsigned int i = 0; i < nb_matrices; i++) {
256  vpMatrix M;
257 #if defined(VISP_HAVE_EIGEN3) || defined(VISP_HAVE_LAPACK) || defined(VISP_HAVE_OPENCV)
258  double det = 0.;
259  // don't put singular matrices in the benchmark
260  for (M = make_random_symmetric_positive_matrix(n); std::fabs(det = M.det()) < .01;
261  M = make_random_symmetric_positive_matrix(n)) {
262  if (verbose) {
263  std::cout << " Generated random symmetric positive matrix A=" << std::endl << M << std::endl;
264  std::cout << " Generated random symmetric positive matrix not "
265  "invertibleL: det="
266  << det << ". Retrying..." << std::endl;
267  }
268  }
269 #else
270  M = make_random_symmetric_positive_matrix(n);
271 #endif
272  bench.push_back(M);
273  }
274 }
275 
276 void create_bench_random_triangular_matrix(unsigned int nb_matrices, unsigned int n, bool verbose,
277  std::vector<vpMatrix> &bench)
278 {
279  if (verbose)
280  std::cout << "Create a bench of " << nb_matrices << " " << n << " by " << n << " triangular matrices" << std::endl;
281  bench.clear();
282  for (unsigned int i = 0; i < nb_matrices; i++) {
283  vpMatrix M;
284 #if defined(VISP_HAVE_EIGEN3) || defined(VISP_HAVE_LAPACK) || defined(VISP_HAVE_OPENCV)
285  double det = 0.;
286  // don't put singular matrices in the benchmark
287  for (M = make_random_triangular_matrix(n); std::fabs(det = M.det()) < .01; M = make_random_triangular_matrix(n)) {
288  if (verbose) {
289  std::cout << " Generated random symmetric positive matrix A=" << std::endl << M << std::endl;
290  std::cout << " Generated random symmetric positive matrix not "
291  "invertibleL: det="
292  << det << ". Retrying..." << std::endl;
293  }
294  }
295 #else
296  M = make_random_triangular_matrix(n);
297 #endif
298  bench.push_back(M);
299  }
300 }
301 
302 int test_inverse(const std::vector<vpMatrix> &bench, const std::vector<vpMatrix> &result)
303 {
304  double epsilon = 1e-10;
305  for (unsigned int i = 0; i < bench.size(); i++) {
306  vpMatrix I = bench[i] * result[i];
307  if (std::fabs(I.frobeniusNorm() - sqrt(static_cast<double>(bench[0].AtA().getRows()))) > epsilon) {
308  std::cout << "Bad inverse[" << i << "]: " << I.frobeniusNorm() << " " << sqrt((double)bench[0].AtA().getRows())
309  << std::endl;
310  return EXIT_FAILURE;
311  }
312  }
313  return EXIT_SUCCESS;
314 }
315 
316 int test_inverse_lu_small(bool verbose, const std::vector<vpMatrix> &bench, double &time)
317 {
318  if (verbose)
319  std::cout << "Test inverse by LU on small matrices" << std::endl;
320  // Compute inverse
321  if (verbose)
322  std::cout << " Inverting " << bench[0].getRows() << "x" << bench[0].getCols() << " small matrix." << std::endl;
323  std::vector<vpMatrix> result(bench.size());
324  double t = vpTime::measureTimeMs();
325  for (unsigned int i = 0; i < bench.size(); i++) {
326  result[i] = bench[i].inverseByLU();
327  }
328  time = vpTime::measureTimeMs() - t;
329 
330  // Test inverse
331  return test_inverse(bench, result);
332 }
333 
334 #if defined(VISP_HAVE_EIGEN3)
335 int test_inverse_lu_eigen3(bool verbose, const std::vector<vpMatrix> &bench, double &time)
336 {
337  if (verbose)
338  std::cout << "Test inverse by LU using Eigen3 3rd party" << std::endl;
339  // Compute inverse
340  if (verbose)
341  std::cout << " Inverting " << bench[0].AtA().getRows() << "x" << bench[0].AtA().getCols()
342  << " matrix using LU decomposition (Eigen3)." << std::endl;
343  std::vector<vpMatrix> result(bench.size());
344  double t = vpTime::measureTimeMs();
345  for (unsigned int i = 0; i < bench.size(); i++) {
346  result[i] = bench[i].AtA().inverseByLUEigen3() * bench[i].transpose();
347  }
348  time = vpTime::measureTimeMs() - t;
349 
350  // Test inverse
351  return test_inverse(bench, result);
352 }
353 #endif
354 
355 #if defined(VISP_HAVE_LAPACK)
356 int test_inverse_lu_lapack(bool verbose, const std::vector<vpMatrix> &bench, double &time)
357 {
358  if (verbose)
359  std::cout << "Test inverse by LU using Lapack 3rd party" << std::endl;
360  // Compute inverse
361  if (verbose)
362  std::cout << " Inverting " << bench[0].AtA().getRows() << "x" << bench[0].AtA().getCols()
363  << " matrix using LU decomposition (Lapack)." << std::endl;
364  std::vector<vpMatrix> result(bench.size());
365  double t = vpTime::measureTimeMs();
366  for (unsigned int i = 0; i < bench.size(); i++) {
367  result[i] = bench[i].AtA().inverseByLULapack() * bench[i].transpose();
368  }
369  time = vpTime::measureTimeMs() - t;
370 
371  // Test inverse
372  return test_inverse(bench, result);
373 }
374 
375 int test_inverse_cholesky_lapack(bool verbose, const std::vector<vpMatrix> &bench, double &time)
376 {
377  if (verbose)
378  std::cout << "Test inverse by Cholesky using Lapack 3rd party" << std::endl;
379  // Compute inverse
380  if (verbose)
381  std::cout << " Inverting " << bench[0].AtA().getRows() << "x" << bench[0].AtA().getCols()
382  << " matrix using cholesky decomposition (Lapack)." << std::endl;
383  std::vector<vpMatrix> result(bench.size());
384  double t = vpTime::measureTimeMs();
385  for (unsigned int i = 0; i < bench.size(); i++) {
386  result[i] = bench[i].AtA().inverseByCholeskyLapack() * bench[i].transpose();
387  }
388  time = vpTime::measureTimeMs() - t;
389 
390  // Test inverse
391  return test_inverse(bench, result);
392 }
393 
394 int test_inverse_qr_lapack(bool verbose, const std::vector<vpMatrix> &bench, double &time)
395 {
396  if (verbose)
397  std::cout << "Test inverse by QR using Lapack 3rd party" << std::endl;
398  // Compute inverse
399  if (verbose)
400  std::cout << " Inverting " << bench[0].AtA().getRows() << "x" << bench[0].AtA().getCols()
401  << " matrix using QR decomposition (Lapack)" << std::endl;
402  std::vector<vpMatrix> result(bench.size());
403  double t = vpTime::measureTimeMs();
404  for (unsigned int i = 0; i < bench.size(); i++) {
405  result[i] = bench[i].AtA().inverseByQRLapack() * bench[i].transpose();
406  }
407  time = vpTime::measureTimeMs() - t;
408 
409  // Test inverse
410  return test_inverse(bench, result);
411 }
412 #endif
413 
414 #if defined(VISP_HAVE_OPENCV)
415 int test_inverse_lu_opencv(bool verbose, const std::vector<vpMatrix> &bench, double &time)
416 {
417  if (verbose)
418  std::cout << "Test inverse by LU using OpenCV 3rd party" << std::endl;
419  // Compute inverse
420  if (verbose)
421  std::cout << " Inverting " << bench[0].AtA().getRows() << "x" << bench[0].AtA().getCols()
422  << " matrix using LU decomposition (OpenCV)" << std::endl;
423  std::vector<vpMatrix> result(bench.size());
424  double t = vpTime::measureTimeMs();
425  for (unsigned int i = 0; i < bench.size(); i++) {
426  result[i] = bench[i].AtA().inverseByLUOpenCV() * bench[i].transpose();
427  }
428  time = vpTime::measureTimeMs() - t;
429 
430  // Test inverse
431  return test_inverse(bench, result);
432 }
433 
434 int test_inverse_cholesky_opencv(bool verbose, const std::vector<vpMatrix> &bench, double &time)
435 {
436  if (verbose)
437  std::cout << "Test inverse by Cholesky using OpenCV 3rd party" << std::endl;
438  // Compute inverse
439  if (verbose)
440  std::cout << " Inverting " << bench[0].AtA().getRows() << "x" << bench[0].AtA().getCols()
441  << " matrix using Cholesky decomposition (OpenCV)" << std::endl;
442  std::vector<vpMatrix> result(bench.size());
443  double t = vpTime::measureTimeMs();
444  for (unsigned int i = 0; i < bench.size(); i++) {
445  result[i] = bench[i].AtA().inverseByCholeskyOpenCV() * bench[i].transpose();
446  }
447  time = vpTime::measureTimeMs() - t;
448 
449  // Test inverse
450  return test_inverse(bench, result);
451 }
452 #endif
453 
454 #if defined(VISP_HAVE_LAPACK) || defined(VISP_HAVE_EIGEN3) || defined(VISP_HAVE_OPENCV)
455 // SVD is only available for these 3rd parties
456 int test_pseudo_inverse(bool verbose, const std::vector<vpMatrix> &bench, double &time)
457 {
458  if (verbose)
459  std::cout << "Test pseudo inverse using either Lapack, Eigen3 or OpenCV 3rd party" << std::endl;
460  // Compute inverse
461  if (verbose)
462  std::cout << " Pseudo inverting " << bench[0].AtA().getRows() << "x" << bench[0].AtA().getCols() << " matrix"
463  << std::endl;
464  std::vector<vpMatrix> result(bench.size());
465  double t = vpTime::measureTimeMs();
466  for (unsigned int i = 0; i < bench.size(); i++) {
467  result[i] = bench[i].AtA().pseudoInverse() * bench[i].transpose();
468  }
469  time = vpTime::measureTimeMs() - t;
470 
471  // Test inverse
472  return test_inverse(bench, result);
473 }
474 
475 int test_inverse_triangular(bool verbose, const std::vector<vpMatrix> &bench, double &time)
476 {
477  if (verbose)
478  std::cout << "Test inverse triangular using Lapack" << std::endl;
479  // Compute inverse
480  if (verbose)
481  std::cout << " Triangular inverse " << bench[0].getRows() << "x" << bench[0].getCols() << " matrix" << std::endl;
482  std::vector<vpMatrix> result(bench.size());
483  double t = vpTime::measureTimeMs();
484  for (unsigned int i = 0; i < bench.size(); i++) {
485  result[i] = bench[i].inverseTriangular(true);
486  }
487  time = vpTime::measureTimeMs() - t;
488 
489  // Test inverse
490  return test_inverse(bench, result);
491 }
492 #endif
493 
494 void save_time(const std::string &method, bool verbose, bool use_plot_file, std::ofstream &of, double time)
495 {
496  if (use_plot_file)
497  of << time << "\t";
498  if (verbose || !use_plot_file) {
499  std::cout << method << time << std::endl;
500  }
501 }
502 
503 int main(int argc, const char *argv[])
504 {
505  try {
506  unsigned int nb_matrices = 1000;
507  unsigned int nb_iterations = 10;
508  unsigned int nb_rows = 6;
509  unsigned int nb_cols = 6;
510  bool verbose = false;
511  std::string plotfile("plot-inv.csv");
512  bool use_plot_file = false;
513  std::ofstream of;
514 
515  // Read the command line options
516  if (getOptions(argc, argv, nb_matrices, nb_iterations, use_plot_file, plotfile, nb_rows, nb_cols, verbose) ==
517  false) {
518  return EXIT_FAILURE;
519  }
520 
521  if (use_plot_file) {
522  of.open(plotfile.c_str());
523  of << "iter"
524  << "\t";
525 
526 #if defined(VISP_HAVE_LAPACK)
527  of << "\"LU Lapack\""
528  << "\t";
529 #endif
530 #if defined(VISP_HAVE_EIGEN3)
531  of << "\"LU Eigen3\""
532  << "\t";
533 #endif
534 #if defined(VISP_HAVE_OPENCV)
535  of << "\"LU OpenCV\""
536  << "\t";
537 #endif
538 
539 #if defined(VISP_HAVE_LAPACK)
540  of << "\"Cholesky Lapack\""
541  << "\t";
542 #endif
543 
544 #if defined(VISP_HAVE_OPENCV)
545  of << "\"Cholesky OpenCV\""
546  << "\t";
547 #endif
548 
549 #if defined(VISP_HAVE_LAPACK)
550  of << "\"QR Lapack\""
551  << "\t";
552 #endif
553 
554 #if defined(VISP_HAVE_LAPACK) || defined(VISP_HAVE_EIGEN3) || defined(VISP_HAVE_OPENCV)
555  of << "\"Pseudo inverse (Lapack, Eigen3 or OpenCV)\""
556  << "\t";
557 #endif
558  of << std::endl;
559  }
560 
561  int ret = EXIT_SUCCESS;
562  for (unsigned int iter = 0; iter < nb_iterations; iter++) {
563  std::vector<vpMatrix> bench_random_matrices_11;
564  create_bench_random_matrix(nb_matrices, 1, 1, verbose, bench_random_matrices_11);
565  std::vector<vpMatrix> bench_random_matrices_22;
566  create_bench_random_matrix(nb_matrices, 2, 2, verbose, bench_random_matrices_22);
567  std::vector<vpMatrix> bench_random_matrices_33;
568  create_bench_random_matrix(nb_matrices, 3, 3, verbose, bench_random_matrices_33);
569 #if defined(VISP_HAVE_EIGEN3) || defined(VISP_HAVE_LAPACK) || defined(VISP_HAVE_OPENCV)
570  std::vector<vpMatrix> bench_random_matrices;
571  create_bench_random_matrix(nb_matrices, nb_rows, nb_cols, verbose, bench_random_matrices);
572  std::vector<vpMatrix> bench_symmetric_positive_matrices;
573  create_bench_symmetric_positive_matrix(nb_matrices, nb_rows, verbose, bench_symmetric_positive_matrices);
574  std::vector<vpMatrix> bench_triangular_matrices;
575  create_bench_random_triangular_matrix(nb_matrices, nb_rows, verbose, bench_triangular_matrices);
576 #endif
577 
578  if (use_plot_file)
579  of << iter << "\t";
580 
581  double time;
582 
583  // LU inverse on 1 by 1 matrices
584  ret += test_inverse_lu_small(verbose, bench_random_matrices_11, time);
585  save_time("Inverse by LU 1x1: ", verbose, use_plot_file, of, time);
586  // LU inverse on 2 by 2 matrices
587  ret += test_inverse_lu_small(verbose, bench_random_matrices_22, time);
588  save_time("Inverse by LU 2x2: ", verbose, use_plot_file, of, time);
589  // LU inverse on 3 by 3 matrices
590  ret += test_inverse_lu_small(verbose, bench_random_matrices_33, time);
591  save_time("Inverse by LU 3x3: ", verbose, use_plot_file, of, time);
592 
593  // LU decomposition
594 #if defined(VISP_HAVE_LAPACK)
595  ret += test_inverse_lu_lapack(verbose, bench_random_matrices, time);
596  save_time("Inverse by LU (Lapack): ", verbose, use_plot_file, of, time);
597 #endif
598 
599 #if defined(VISP_HAVE_EIGEN3)
600  ret += test_inverse_lu_eigen3(verbose, bench_random_matrices, time);
601  save_time("Inverse by LU (Eigen3): ", verbose, use_plot_file, of, time);
602 #endif
603 
604 #if defined(VISP_HAVE_OPENCV)
605  ret += test_inverse_lu_opencv(verbose, bench_random_matrices, time);
606  save_time("Inverse by LU (OpenCV): ", verbose, use_plot_file, of, time);
607 #endif
608 
609  // Cholesky for symmetric positive matrices
610 #if defined(VISP_HAVE_LAPACK)
611  ret += test_inverse_cholesky_lapack(verbose, bench_symmetric_positive_matrices, time);
612  save_time("Inverse by Cholesly (Lapack): ", verbose, use_plot_file, of, time);
613 #endif
614 
615 #if defined(VISP_HAVE_OPENCV)
616  ret += test_inverse_cholesky_opencv(verbose, bench_symmetric_positive_matrices, time);
617  save_time("Inverse by Cholesky (OpenCV): ", verbose, use_plot_file, of, time);
618 #endif
619 
620  // QR decomposition
621 #if defined(VISP_HAVE_LAPACK)
622  ret += test_inverse_qr_lapack(verbose, bench_random_matrices, time);
623  save_time("Inverse by QR (Lapack): ", verbose, use_plot_file, of, time);
624 #endif
625 
626  // Pseudo-inverse with SVD
627 #if defined(VISP_HAVE_LAPACK) || defined(VISP_HAVE_EIGEN3) || defined(VISP_HAVE_OPENCV)
628  ret += test_pseudo_inverse(verbose, bench_random_matrices, time);
629  save_time("Pseudo inverse (Lapack, Eigen3, OpenCV): ", verbose, use_plot_file, of, time);
630 #endif
631 
632  // Test inverse triangular
633 #if defined(VISP_HAVE_LAPACK)
634  ret += test_inverse_triangular(verbose, bench_triangular_matrices, time);
635  save_time("Triangular inverse (Lapack): ", verbose, use_plot_file, of, time);
636 #endif
637 
638  if (use_plot_file)
639  of << std::endl;
640  }
641  if (use_plot_file) {
642  of.close();
643  std::cout << "Result saved in " << plotfile << std::endl;
644  }
645 
646  if (ret == EXIT_SUCCESS) {
647  std::cout << "Test succeed" << std::endl;
648  }
649  else {
650  std::cout << "Test failed" << std::endl;
651  }
652 
653  return ret;
654  }
655  catch (const vpException &e) {
656  std::cout << "Catch an exception: " << e.getStringMessage() << std::endl;
657  return EXIT_FAILURE;
658  }
659 }
unsigned int getCols() const
Definition: vpArray2D.h:337
void resize(unsigned int nrows, unsigned int ncols, bool flagNullify=true, bool recopy_=true)
Definition: vpArray2D.h:362
unsigned int getRows() const
Definition: vpArray2D.h:347
error that can be emitted by ViSP classes.
Definition: vpException.h:60
const std::string & getStringMessage() const
Definition: vpException.cpp:67
Implementation of a matrix and operations on matrices.
Definition: vpMatrix.h:169
vpMatrix AtA() const
double det(vpDetMethod method=LU_DECOMPOSITION) const
Definition: vpMatrix.cpp:1641
vpMatrix t() const
static bool parse(int *argcPtr, const char **argv, vpArgvInfo *argTable, int flags)
Definition: vpParseArgv.cpp:70
VISP_EXPORT double measureTimeMs()