ViSP  2.8.0
vpMatrix_cholesky.cpp
1 /****************************************************************************
2  *
3  * $Id: vpMatrix_lu.cpp 3530 2012-01-03 10:52:12Z fspindle $
4  *
5  * This file is part of the ViSP software.
6  * Copyright (C) 2005 - 2013 by INRIA. All rights reserved.
7  *
8  * This software is free software; you can redistribute it and/or
9  * modify it under the terms of the GNU General Public License
10  * ("GPL") version 2 as published by the Free Software Foundation.
11  * See the file LICENSE.txt at the root directory of this source
12  * distribution for additional information about the GNU GPL.
13  *
14  * For using ViSP with software that can not be combined with the GNU
15  * GPL, please contact INRIA about acquiring a ViSP Professional
16  * Edition License.
17  *
18  * See http://www.irisa.fr/lagadic/visp/visp.html for more information.
19  *
20  * This software was developed at:
21  * INRIA Rennes - Bretagne Atlantique
22  * Campus Universitaire de Beaulieu
23  * 35042 Rennes Cedex
24  * France
25  * http://www.irisa.fr/lagadic
26  *
27  * If you have questions regarding the use of this file, please contact
28  * INRIA at visp@inria.fr
29  *
30  * This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
31  * WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
32  *
33  *
34  * Description:
35  * Matrix Cholesky decomposition.
36  *
37  * Authors:
38  * Filip Novotny
39  *
40  *****************************************************************************/
41 
42 #include <visp/vpConfig.h>
43 
44 #include <visp/vpMatrix.h>
45 #include <visp/vpMath.h>
46 #include <visp/vpColVector.h>
47 
48 // Exception
49 #include <visp/vpException.h>
50 #include <visp/vpMatrixException.h>
51 
52 // Debug trace
53 #include <visp/vpDebug.h>
54 
55 #ifdef VISP_HAVE_LAPACK
56 extern "C" void dpotrf_ (char *uplo, int *n, double *a, int *lda, int *info);
57 extern "C" int dpotri_(char *uplo, int *n, double *a, int *lda, int *info);
58 
59 #endif
60 
61 #ifdef VISP_HAVE_LAPACK
63  int rowNum = (int)this->getRows();
64  int lda = (int)rowNum; //lda is the number of rows because we don't use a submatrix
65  int info;
66 
67  vpMatrix A = *this;
68  dpotrf_((char*)"L",&rowNum,A.data,&lda,&info);
69 
70  if(info!=0)
71  std::cout << "cholesky:dpotrf_:error" << std::endl;
72 
73  dpotri_((char*)"L",&rowNum,A.data,&lda,&info);
74  if(info!=0){
75  std::cout << "cholesky:dpotri_:error" << std::endl;
77  }
78 
79 
80  for(unsigned int i=0;i<A.getRows();i++)
81  for(unsigned int j=0;j<A.getCols();j++)
82  if(i>j) A[i][j] = A[j][i];
83 
84 
85  return A;
86 
87 }
88 #endif
89 
122 #if defined(VISP_HAVE_LAPACK)
123 vpMatrix
125 {
126 
127  if ( rowNum != colNum)
128  {
129  vpERROR_TRACE("\n\t\tCannot invert a non-square vpMatrix") ;
131  "Cannot invert a non-square vpMatrix")) ;
132  }
133 #ifdef VISP_HAVE_LAPACK
134  return inverseByCholeskyLapack();
135 #endif
136 }
137 
138 #endif
Definition of the vpMatrix class.
Definition: vpMatrix.h:96
#define vpERROR_TRACE
Definition: vpDebug.h:379
vpMatrix inverseByCholesky() const
vpMatrix inverseByCholeskyLapack() const
double * data
address of the first element of the data array
Definition: vpMatrix.h:116
unsigned int rowNum
number of rows
Definition: vpMatrix.h:110
unsigned int getCols() const
Return the number of columns of the matrix.
Definition: vpMatrix.h:159
error that can be emited by the vpMatrix class and its derivates
unsigned int colNum
number of columns
Definition: vpMatrix.h:112
unsigned int getRows() const
Return the number of rows of the matrix.
Definition: vpMatrix.h:157