46#include <visp3/core/vpColVector.h>
47#include <visp3/core/vpMatrix.h>
48#include <visp3/core/vpTime.h>
49#include <visp3/io/vpParseArgv.h>
52#define GETOPTARGS "cdn:i:pf:R:C:vh"
62void usage(
const char *name,
const char *badparam)
65Test matrix pseudo-inverse.\n\
66Outputs a comparison of the results obtained by supported 3rd parties.\n\
69 %s [-n <number of matrices>] [-f <plot filename>]\n\
70 [-R <number of rows>] [-C <number of columns>]\n\
71 [-i <number of iterations>] [-p] [-h]\n",
76 -n <number of matrices> \n\
77 Number of matrices inverted during each test loop.\n\
79 -i <number of iterations> \n\
80 Number of iterations of the test.\n\
82 -f <plot filename> \n\
83 Set output path for plot output.\n\
84 The plot logs the times of \n\
85 the different inversion methods: \n\
86 QR,LU,Cholesky and Pseudo-inverse.\n\
88 -R <number of rows>\n\
89 Number of rows of the automatically generated matrices \n\
92 -C <number of columns>\n\
93 Number of colums of the automatically generated matrices \n\
97 Plot into filename in the gnuplot format. \n\
98 If this option is used, tests results will be logged \n\
99 into a filename specified with -f.\n\
102 Print the help.\n\n");
105 fprintf(stderr,
"ERROR: \n");
106 fprintf(stderr,
"\nBad parameter [%s]\n", badparam);
117bool getOptions(
int argc,
const char **argv,
unsigned int &nb_matrices,
unsigned int &nb_iterations,
118 bool &use_plot_file, std::string &plotfile,
unsigned int &nbrows,
unsigned int &nbcols,
bool &verbose)
126 usage(argv[0], NULL);
130 nb_matrices = (
unsigned int)atoi(optarg_);
133 nb_iterations = (
unsigned int)atoi(optarg_);
137 use_plot_file =
true;
140 use_plot_file =
true;
143 nbrows = (
unsigned int)atoi(optarg_);
146 nbcols = (
unsigned int)atoi(optarg_);
157 usage(argv[0], optarg_);
163 if ((c == 1) || (c == -1)) {
165 usage(argv[0], NULL);
166 std::cerr <<
"ERROR: " << std::endl;
167 std::cerr <<
" Bad argument " << optarg_ << std::endl << std::endl;
174vpMatrix make_random_matrix(
unsigned int nbrows,
unsigned int nbcols)
179 for (
unsigned int i = 0; i < A.
getRows(); i++) {
180 for (
unsigned int j = 0; j < A.
getCols(); j++) {
181 A[i][j] = (double)rand() / (double)RAND_MAX;
188void create_bench_random_matrix(
unsigned int nb_matrices,
unsigned int nb_rows,
unsigned int nb_cols,
bool verbose,
189 std::vector<vpMatrix> &bench)
192 std::cout <<
"Create a bench of " << nb_matrices <<
" " << nb_rows <<
" by " << nb_cols <<
" matrices" << std::endl;
194 for (
unsigned int i = 0; i < nb_matrices; i++) {
195 vpMatrix M = make_random_matrix(nb_rows, nb_cols);
200int test_pseudo_inverse(
const std::vector<vpMatrix> &A,
const std::vector<vpMatrix> &Api)
202 double allowed_error = 1e-3;
204 for (
unsigned int i = 0; i < A.size(); i++) {
205 double error = (A[i] * Api[i] * A[i] - A[i]).frobeniusNorm();
206 if (error > allowed_error) {
207 std::cout <<
"Bad pseudo-inverse [" << i <<
"]: euclidean norm: " << error << std::endl;
214int test_pseudo_inverse(
const std::vector<vpMatrix> &A,
const std::vector<vpMatrix> &Api,
215 const std::vector<vpColVector> &sv,
const std::vector<vpMatrix> &imA,
216 const std::vector<vpMatrix> &imAt,
const std::vector<vpMatrix> &kerAt)
218 double allowed_error = 1e-3;
220 for (
unsigned int i = 0; i < A.size(); i++) {
221 double error = (A[i] * Api[i] * A[i] - A[i]).frobeniusNorm();
222 if (error > allowed_error) {
223 std::cout <<
"Bad pseudo-inverse [" << i <<
"]: euclidean norm: " << error << std::endl;
229 for (
unsigned int i = 0; i < kerAt.size(); i++) {
230 if (kerAt[i].size()) {
231 vpMatrix nullspace = A[i] * kerAt[i].
t();
234 if (error > allowed_error) {
235 std::cout <<
"Bad kernel [" << i <<
"]: euclidean norm: " << error << std::endl;
242 for (
unsigned int i = 0; i < kerAt.size(); i++) {
243 unsigned int rank = imA[i].getCols();
244 vpMatrix U, S(rank, A[i].getCols()), Vt(A[i].getCols(), A[i].getCols());
247 for (
unsigned int j = 0; j < rank; j++)
250 Vt.
insert(imAt[i].t(), 0, 0);
251 Vt.insert(kerAt[i], imAt[i].getCols(), 0);
253 double error = (U * S * Vt - A[i]).frobeniusNorm();
255 if (error > allowed_error) {
256 std::cout <<
"Bad imA, imAt, sv, kerAt [" << i <<
"]: euclidean norm: " << error << std::endl;
264int test_pseudo_inverse_default(
bool verbose,
const std::vector<vpMatrix> &bench, std::vector<double> &time)
267 std::cout <<
"Test pseudo-inverse using default 3rd party" << std::endl;
269 std::cout <<
" Pseudo-inverse on a " << bench[0].getRows() <<
"x" << bench[0].getCols() <<
" matrix" << std::endl;
271 size_t size = bench.size();
272 std::vector<vpMatrix> PI(size), imA(size), imAt(size), kerAt(size);
273 std::vector<vpColVector> sv(size);
274 int ret = EXIT_SUCCESS;
277 unsigned int test = 0;
279 for (
unsigned int i = 0; i < bench.size(); i++) {
280 PI[i] = bench[i].pseudoInverse();
283 for (
unsigned int i = 0; i < time.size(); i++) {
284 ret += test_pseudo_inverse(bench, PI);
290 for (
unsigned int i = 0; i < bench.size(); i++) {
291 bench[i].pseudoInverse(PI[i]);
294 for (
unsigned int i = 0; i < time.size(); i++) {
295 ret += test_pseudo_inverse(bench, PI);
301 for (
unsigned int i = 0; i < bench.size(); i++) {
302 bench[i].pseudoInverse(PI[i], sv[i]);
305 for (
unsigned int i = 0; i < time.size(); i++) {
306 ret += test_pseudo_inverse(bench, PI);
312 for (
unsigned int i = 0; i < bench.size(); i++) {
313 bench[i].pseudoInverse(PI[i], sv[i], 1e-6, imA[i], imAt[i], kerAt[i]);
317 for (
unsigned int i = 0; i < time.size(); i++) {
318 ret += test_pseudo_inverse(bench, PI, sv, imA, imAt, kerAt);
324#if defined(VISP_HAVE_EIGEN3)
325int test_pseudo_inverse_eigen3(
bool verbose,
const std::vector<vpMatrix> &bench, std::vector<double> &time)
328 std::cout <<
"Test pseudo-inverse using Eigen3 3rd party" << std::endl;
330 std::cout <<
" Pseudo-inverse on a " << bench[0].getRows() <<
"x" << bench[0].getCols() <<
" matrix" << std::endl;
332 size_t size = bench.size();
333 std::vector<vpMatrix> PI(size), imA(size), imAt(size), kerAt(size);
334 std::vector<vpColVector> sv(size);
335 int ret = EXIT_SUCCESS;
338 unsigned int test = 0;
340 for (
unsigned int i = 0; i < bench.size(); i++) {
341 PI[i] = bench[i].pseudoInverseEigen3();
344 for (
unsigned int i = 0; i < time.size(); i++) {
345 ret += test_pseudo_inverse(bench, PI);
351 for (
unsigned int i = 0; i < bench.size(); i++) {
352 bench[i].pseudoInverseEigen3(PI[i]);
355 for (
unsigned int i = 0; i < time.size(); i++) {
356 ret += test_pseudo_inverse(bench, PI);
362 for (
unsigned int i = 0; i < bench.size(); i++) {
363 bench[i].pseudoInverseEigen3(PI[i], sv[i]);
366 for (
unsigned int i = 0; i < time.size(); i++) {
367 ret += test_pseudo_inverse(bench, PI);
373 for (
unsigned int i = 0; i < bench.size(); i++) {
374 bench[i].pseudoInverseEigen3(PI[i], sv[i], 1e-6, imA[i], imAt[i], kerAt[i]);
378 for (
unsigned int i = 0; i < time.size(); i++) {
379 ret += test_pseudo_inverse(bench, PI, sv, imA, imAt, kerAt);
386#if defined(VISP_HAVE_LAPACK)
387int test_pseudo_inverse_lapack(
bool verbose,
const std::vector<vpMatrix> &bench, std::vector<double> &time)
390 std::cout <<
"Test pseudo-inverse using Lapack 3rd party" << std::endl;
392 std::cout <<
" Pseudo-inverse on a " << bench[0].getRows() <<
"x" << bench[0].getCols() <<
" matrix" << std::endl;
394 size_t size = bench.size();
395 std::vector<vpMatrix> PI(size), imA(size), imAt(size), kerAt(size);
396 std::vector<vpColVector> sv(size);
397 int ret = EXIT_SUCCESS;
400 unsigned int test = 0;
402 for (
unsigned int i = 0; i < bench.size(); i++) {
403 PI[i] = bench[i].pseudoInverseLapack();
406 for (
unsigned int i = 0; i < time.size(); i++) {
407 ret += test_pseudo_inverse(bench, PI);
413 for (
unsigned int i = 0; i < bench.size(); i++) {
414 bench[i].pseudoInverseLapack(PI[i]);
417 for (
unsigned int i = 0; i < time.size(); i++) {
418 ret += test_pseudo_inverse(bench, PI);
424 for (
unsigned int i = 0; i < bench.size(); i++) {
425 bench[i].pseudoInverseLapack(PI[i], sv[i]);
428 for (
unsigned int i = 0; i < time.size(); i++) {
429 ret += test_pseudo_inverse(bench, PI);
435 for (
unsigned int i = 0; i < bench.size(); i++) {
436 bench[i].pseudoInverseLapack(PI[i], sv[i], 1e-6, imA[i], imAt[i], kerAt[i]);
440 for (
unsigned int i = 0; i < time.size(); i++) {
441 ret += test_pseudo_inverse(bench, PI, sv, imA, imAt, kerAt);
448#if defined(VISP_HAVE_OPENCV)
449int test_pseudo_inverse_opencv(
bool verbose,
const std::vector<vpMatrix> &bench, std::vector<double> &time)
452 std::cout <<
"Test pseudo-inverse using OpenCV 3rd party" << std::endl;
454 std::cout <<
" Pseudo-inverse on a " << bench[0].getRows() <<
"x" << bench[0].getCols() <<
" matrix" << std::endl;
456 size_t size = bench.size();
457 std::vector<vpMatrix> PI(size), imA(size), imAt(size), kerAt(size);
458 std::vector<vpColVector> sv(size);
459 int ret = EXIT_SUCCESS;
462 unsigned int test = 0;
464 for (
unsigned int i = 0; i < bench.size(); i++) {
465 PI[i] = bench[i].pseudoInverseOpenCV();
468 for (
unsigned int i = 0; i < time.size(); i++) {
469 ret += test_pseudo_inverse(bench, PI);
475 for (
unsigned int i = 0; i < bench.size(); i++) {
476 bench[i].pseudoInverseOpenCV(PI[i]);
479 for (
unsigned int i = 0; i < time.size(); i++) {
480 ret += test_pseudo_inverse(bench, PI);
486 for (
unsigned int i = 0; i < bench.size(); i++) {
487 bench[i].pseudoInverseOpenCV(PI[i], sv[i]);
490 for (
unsigned int i = 0; i < time.size(); i++) {
491 ret += test_pseudo_inverse(bench, PI);
497 for (
unsigned int i = 0; i < bench.size(); i++) {
498 bench[i].pseudoInverseOpenCV(PI[i], sv[i], 1e-6, imA[i], imAt[i], kerAt[i]);
502 for (
unsigned int i = 0; i < time.size(); i++) {
503 ret += test_pseudo_inverse(bench, PI, sv, imA, imAt, kerAt);
510void save_time(
const std::string &method,
unsigned int nrows,
unsigned int ncols,
bool verbose,
bool use_plot_file,
511 std::ofstream &of,
const std::vector<double> &time)
513 for (
size_t i = 0; i < time.size(); i++) {
515 of << time[i] <<
"\t";
517 std::cout <<
" " << method <<
" svd(" << nrows <<
"x" << ncols <<
")"
518 <<
" test " << i <<
": " << time[i] << std::endl;
523int main(
int argc,
const char *argv[])
526#if defined(VISP_HAVE_EIGEN3) || defined(VISP_HAVE_LAPACK) || defined(VISP_HAVE_OPENCV)
527 unsigned int nb_matrices = 10;
528 unsigned int nb_iterations = 10;
529 unsigned int nb_rows = 12;
530 unsigned int nb_cols = 6;
531 bool verbose =
false;
532 std::string plotfile(
"plot-pseudo-inv.csv");
533 bool use_plot_file =
false;
536 unsigned int nb_svd_functions = 4;
537 unsigned int nb_test_matrix_size = 3;
538 std::vector<double> time(nb_svd_functions);
539 std::vector<unsigned int> nrows(nb_test_matrix_size), ncols(nb_test_matrix_size);
542 if (getOptions(argc, argv, nb_matrices, nb_iterations, use_plot_file, plotfile, nb_rows, nb_cols, verbose) ==
547 for (
unsigned int s = 0; s < nb_test_matrix_size; s++) {
562 of.open(plotfile.c_str());
566 for (
unsigned int s = 0; s < nb_test_matrix_size; s++) {
567 for (
unsigned int i = 0; i < nb_svd_functions; i++)
568 of <<
"\"default " << nrows[s] <<
"x" << ncols[s] <<
" test " << i <<
"\""
571#if defined(VISP_HAVE_LAPACK)
572 for (
unsigned int i = 0; i < nb_svd_functions; i++)
573 of <<
"\"Lapack " << nrows[s] <<
"x" << ncols[s] <<
" test " << i <<
"\""
576#if defined(VISP_HAVE_EIGEN3)
577 for (
unsigned int i = 0; i < nb_svd_functions; i++)
578 of <<
"\"Eigen3 " << nrows[s] <<
"x" << ncols[s] <<
" test " << i <<
"\""
581#if defined(VISP_HAVE_OPENCV)
582 for (
unsigned int i = 0; i < nb_svd_functions; i++)
583 of <<
"\"OpenCV " << nrows[s] <<
"x" << ncols[s] <<
" test " << i <<
"\""
590 int ret = EXIT_SUCCESS;
591 for (
unsigned int iter = 0; iter < nb_iterations; iter++) {
596 for (
unsigned int s = 0; s < nb_test_matrix_size; s++) {
597 std::vector<vpMatrix> bench_random_matrices;
598 create_bench_random_matrix(nb_matrices, nrows[s], ncols[s], verbose, bench_random_matrices);
600 ret += test_pseudo_inverse_default(verbose, bench_random_matrices, time);
601 save_time(
"default -", nrows[s], ncols[s], verbose, use_plot_file, of, time);
603#if defined(VISP_HAVE_LAPACK)
604 ret += test_pseudo_inverse_lapack(verbose, bench_random_matrices, time);
605 save_time(
"Lapack -", nrows[s], ncols[s], verbose, use_plot_file, of, time);
608#if defined(VISP_HAVE_EIGEN3)
609 ret += test_pseudo_inverse_eigen3(verbose, bench_random_matrices, time);
610 save_time(
"Eigen3 -", nrows[s], ncols[s], verbose, use_plot_file, of, time);
613#if defined(VISP_HAVE_OPENCV)
614 ret += test_pseudo_inverse_opencv(verbose, bench_random_matrices, time);
615 save_time(
"OpenCV -", nrows[s], ncols[s], verbose, use_plot_file, of, time);
623 std::cout <<
"Result saved in " << plotfile << std::endl;
626 if (ret == EXIT_SUCCESS) {
627 std::cout <<
"Test succeed" << std::endl;
629 std::cout <<
"Test failed" << std::endl;
636 std::cout <<
"Test does nothing since you dont't have Lapack, Eigen3 or OpenCV 3rd party" << std::endl;
unsigned int getCols() const
void resize(unsigned int nrows, unsigned int ncols, bool flagNullify=true, bool recopy_=true)
unsigned int getRows() const
error that can be emitted by ViSP classes.
const std::string & getStringMessage() const
Implementation of a matrix and operations on matrices.
double frobeniusNorm() const
void insert(const vpMatrix &A, unsigned int r, unsigned int c)
static bool parse(int *argcPtr, const char **argv, vpArgvInfo *argTable, int flags)
VISP_EXPORT double measureTimeMs()