46 #include <visp3/core/vpTime.h>
48 #include <visp3/core/vpMatrix.h>
49 #include <visp3/core/vpColVector.h>
50 #include <visp3/io/vpParseArgv.h>
58 #define GETOPTARGS "cdn:i:pf:R:C:vh"
60 void usage(
const char *name,
const char *badparam);
61 bool getOptions(
int argc,
const char **argv,
62 unsigned int& nb_matrices,
unsigned int& nb_iterations,
63 bool& use_plot_file, std::string& plotfile,
64 unsigned int& nbrows,
unsigned int& nbcols,
bool& verbose);
65 void writeTime(
const char *name,
double time);
66 vpMatrix makeRandomMatrix(
unsigned int nbrows,
unsigned int nbcols);
76 void usage(
const char *name,
const char *badparam)
79 Test matrix inversions\n\
80 using LU, QR and Cholesky methods as well as Pseudo-inverse.\n\
81 Outputs a comparison of these methods.\n\
84 %s [-n <number of matrices>] [-f <plot filename>]\n\
85 [-R <number of rows>] [-C <number of columns>]\n\
86 [-i <number of iterations>] [-p] [-h]\n", name);
90 -n <number of matrices> \n\
91 Number of matrices inverted during each test loop.\n\
93 -i <number of iterations> \n\
94 Number of iterations of the test.\n\
96 -f <plot filename> \n\
97 Set output path for plot output.\n\
98 The plot logs the times of \n\
99 the different inversion methods: \n\
100 QR,LU,Cholesky and Pseudo-inverse.\n\
102 -R <number of rows>\n\
103 Number of rows of the automatically generated matrices \n\
106 -C <number of columns>\n\
107 Number of colums of the automatically generated matrices \n\
111 Plot into filename in the gnuplot format. \n\
112 If this option is used, tests results will be logged \n\
113 into a filename specified with -f.\n\
116 Print the help.\n\n");
119 fprintf(stderr,
"ERROR: \n" );
120 fprintf(stderr,
"\nBad parameter [%s]\n", badparam);
131 bool getOptions(
int argc,
const char **argv,
132 unsigned int& nb_matrices,
unsigned int& nb_iterations,
133 bool& use_plot_file, std::string& plotfile,
134 unsigned int& nbrows,
unsigned int& nbcols,
bool& verbose)
141 case 'h': usage(argv[0], NULL);
return false;
break;
143 nb_matrices = (
unsigned int)atoi(optarg_);
146 nb_iterations = (
unsigned int)atoi(optarg_);
150 use_plot_file =
true;
153 use_plot_file =
true;
156 nbrows = (
unsigned int)atoi(optarg_);
159 nbcols = (
unsigned int)atoi(optarg_);
170 usage(argv[0], optarg_);
175 if ((c == 1) || (c == -1)) {
177 usage(argv[0], NULL);
178 std::cerr <<
"ERROR: " << std::endl;
179 std::cerr <<
" Bad argument " << optarg_ << std::endl << std::endl;
186 void writeTime(
const char *name,
double time)
188 std::cout << name <<
" time=" << time <<
" ms." << std::endl;
191 vpMatrix makeRandomMatrix(
unsigned int nbrows,
unsigned int nbcols)
196 for (
unsigned int i=0 ; i < A.
getRows() ; i++)
197 for (
unsigned int j=0 ; j < A.
getCols() ; j++)
198 A[i][j] = (
double)rand()/(double)RAND_MAX;
204 main(
int argc,
const char ** argv)
206 #ifdef VISP_HAVE_LAPACK_C
208 unsigned int nb_matrices=1000;
209 unsigned int nb_iterations=10;
210 unsigned int nb_rows = 6;
211 unsigned int nb_cols = 6;
212 bool verbose =
false;
213 std::string plotfile(
"plot.txt");
214 bool use_plot_file=
false;
217 double t, qr_time, lu_time,pi_time,chol_time;
219 if (getOptions(argc, argv, nb_matrices,nb_iterations,use_plot_file,plotfile,nb_rows,nb_cols,verbose) ==
false) {
224 of.open(plotfile.c_str());
227 for(
unsigned int iter=0;iter<nb_iterations;iter++){
228 std::vector<vpMatrix> benchQR;
229 std::vector<vpMatrix> benchLU;
230 std::vector<vpMatrix> benchCholesky;
231 std::vector<vpMatrix> benchPseudoInverse;
233 std::cout <<
"********* generating matrices for iteration " << iter <<
"." << std::endl;
234 for(
unsigned int i=0;i<nb_matrices;i++){
238 for(cur=makeRandomMatrix(nb_rows,nb_cols);std::abs(det=cur.
AtA().
det())<.01;cur = makeRandomMatrix(nb_rows,nb_cols))
240 std::cout <<
"Generated random matrix A*tA=" << std::endl << cur.
AtA() << std::endl;
241 std::cout <<
"generated random matrix not invertibleL: det="<<det<<
". Retrying..." << std::endl;
243 benchCholesky.push_back(cur);
244 benchQR.push_back(cur);
245 benchLU.push_back(cur);
246 benchPseudoInverse.push_back(cur);
250 std::cout <<
"\t Inverting " << benchCholesky[0].
AtA().
getRows() <<
"x" << benchCholesky[0].AtA().getCols() <<
" matrix using cholesky decomposition." << std::endl;
252 for(
unsigned int i=0;i<nb_matrices;i++){
253 benchCholesky[i]=benchCholesky[i].AtA().inverseByCholesky()*benchCholesky[i].transpose();
258 std::cout <<
"\t Inverting " << benchLU[0].AtA().getRows() <<
"x" << benchLU[0].AtA().getCols() <<
" matrix using LU decomposition." << std::endl;
260 for(
unsigned int i=0;i<nb_matrices;i++)
261 benchLU[i] = benchLU[i].AtA().inverseByLU()*benchLU[i].transpose();
265 std::cout <<
"\t Inverting " << benchQR[0].AtA().getRows() <<
"x" << benchQR[0].AtA().getCols() <<
" matrix using QR decomposition." << std::endl;
267 for(
unsigned int i=0;i<nb_matrices;i++){
268 benchQR[i]=benchQR[i].AtA().inverseByQR()*benchQR[i].transpose();
273 std::cout <<
"\t Inverting " << benchPseudoInverse[0].AtA().getRows() <<
"x" << benchPseudoInverse[0].AtA().getCols() <<
" matrix while computing pseudo-inverse." << std::endl;
275 for(
unsigned int i=0;i<nb_matrices;i++){
276 benchPseudoInverse[i]=benchPseudoInverse[i].pseudoInverse();
280 double avg_err_lu_qr=0.;
281 double avg_err_lu_pi=0.;
282 double avg_err_lu_chol=0.;
283 double avg_err_qr_pi=0.;
284 double avg_err_qr_chol=0.;
285 double avg_err_pi_chol=0.;
287 for(
unsigned int i=0;i<nb_matrices;i++){
288 avg_err_lu_qr+= (benchQR[i]-benchLU[i]).euclideanNorm();
289 avg_err_lu_pi+= (benchPseudoInverse[i]-benchLU[i]).euclideanNorm();
290 avg_err_qr_pi+= (benchPseudoInverse[i]-benchQR[i]).euclideanNorm();
291 avg_err_qr_chol+= (benchCholesky[i]-benchQR[i]).euclideanNorm();
292 avg_err_lu_chol+= (benchCholesky[i]-benchLU[i]).euclideanNorm();
293 avg_err_pi_chol+= (benchCholesky[i]-benchPseudoInverse[i]).euclideanNorm();
296 avg_err_lu_qr/=nb_matrices;
297 avg_err_lu_pi/=nb_matrices;
298 avg_err_qr_pi/=nb_matrices;
301 of << iter <<
"\t" << lu_time <<
"\t" << qr_time <<
"\t" << pi_time <<
"\t" << chol_time <<
"\t" << avg_err_lu_qr <<
"\t" << avg_err_qr_pi <<
"\t" << avg_err_lu_pi <<
"\t" << avg_err_qr_chol <<
"\t" << avg_err_lu_chol <<
"\t" << avg_err_pi_chol << std::endl;
303 if(verbose || !use_plot_file){
304 writeTime(
"LU",lu_time);
305 writeTime(
"QR",qr_time);
306 writeTime(
"Pseudo-inverse",pi_time);
307 writeTime(
"Cholesky",chol_time);
313 std::cout <<
"Catch an exception: " << e << std::endl;
320 std::cout <<
"You don't have lapack installed" << std::endl;
Implementation of a matrix and operations on matrices.
void resize(const unsigned int nrows, const unsigned int ncols, const bool flagNullify=true)
error that can be emited by ViSP classes.
unsigned int getCols() const
Return the number of columns of the 2D array.
static bool parse(int *argcPtr, const char **argv, vpArgvInfo *argTable, int flags)
unsigned int getRows() const
Return the number of rows of the 2D array.
double det(vpDetMethod method=LU_DECOMPOSITION) const
VISP_EXPORT double measureTimeMs()