24 #define GRT_DLL_EXPORTS
25 #include "LUDecomposition.h"
29 LUDecomposition::LUDecomposition(
const MatrixFloat &a) : sing(false){
31 debugLog.setProceedingText(
"[DEBUG LUDecomposition]");
32 errorLog.setProceedingText(
"[ERROR LUDecomposition]");
33 warningLog.setProceedingText(
"[WARNING LUDecomposition]");
40 const Float TINY=1.0e-20;
41 unsigned int i,imax,j,k;
49 if ((temp=fabs( lu[i][j] )) > big) big=temp;
52 errorLog <<
"Error in LUDecomposition constructor, big == 0.0" << std::endl;
60 temp=vv[i]*fabs(lu[i][k]);
69 lu[imax][j] = lu[k][j];
76 if (lu[k][k] == 0.0) lu[k][k] = TINY;
77 for (i=k+1; i<N; i++) {
78 temp = lu[i][k] /= lu[k][k];
80 lu[i][j] -= temp * lu[k][j];
86 LUDecomposition::~LUDecomposition(){
92 int i=0,ii=0,ip=0,j=0;
96 if (b.size() != N || x.size() != N){
97 errorLog <<
"solve_vector(const VectorFloat &b,VectorFloat &x) - the size of the two vectors does not match!" << std::endl;
100 for (i=0;i<n;i++) x[i] = b[i];
106 for (j=ii-1;j<i;j++) sum -= lu[i][j] * x[j];
111 for (i=n-1;i>=0;i--) {
113 for (j=i+1;j<n;j++) sum -= lu[i][j] * x[j];
114 x[i] = sum / lu[i][i];
124 errorLog <<
"solve(const MatrixFloat &b,MatrixFloat &x) - the size of the two matrices does not match!" << std::endl;
128 for (
unsigned int j=0; j<m; j++) {
129 for(
unsigned int i=0; i<N; i++) xx[i] = b[i][j];
131 for(
unsigned int i=0; i<N; i++) x[i][j] = xx[i];
141 for (j=0;j<N;j++) ainv[i][j] = 0.0;
144 return solve(ainv,ainv);
147 Float LUDecomposition::det()
150 for (
unsigned int i=0;i<N;i++) dd *= lu[i][i];
162 sdp += (LongFloat) aref[i][j] * (LongFloat)x[j];
165 if( !solve_vector(r,r) ){
168 for (i=0;i<N;i++) x[i] -= r[i];
172 bool LUDecomposition::getIsSingular(){
unsigned int getNumRows() const
unsigned int getNumCols() const
virtual bool resize(const unsigned int r, const unsigned int c)