Main Page   Class Hierarchy   Alphabetical List   Compound List   File List   Compound Members   Related Pages  

AFSymMatrix.hpp

00001 #ifndef __AFSymMatrix_H__
00002 #define __AFSymMatrix_H__
00003 
00004 
00005 
00026 #include <stdio.h>
00027 #include "FloatVector.hpp"
00028 
00029 #include "AFloatMatrix.hpp"
00030 
00031 class FDiagonalMatrix;
00032 class AFProductSeqSMatrix;
00033 
00034 class AFSymMatrix : public AFloatMatrix
00035 {
00036 protected:
00037 
00038 int width;
00039 
00040 public:
00041 AFSymMatrix(int w){ width= w; }
00042 virtual ~AFSymMatrix(){}
00043 
00044 virtual AFSymMatrix* copyAnn()
00045                         {printf("%s() %s: not Implemented\n",__FUNCTION__, __FILE__); return 0;}
00046 virtual AFSymMatrix* copyFull();
00047 
00048 virtual FDiagonalMatrix* copyDiagonal();
00049 // no diagonal
00050 virtual AFSymMatrix* copyL0();                          
00051 virtual AFSymMatrix* copyU0();
00052 // with diagonal
00053 virtual AFSymMatrix* copyL();                           
00054 virtual AFSymMatrix* copyU();
00055 // with unit diagonal
00056 virtual AFSymMatrix* copyL1();                          
00057 virtual AFSymMatrix* copyU1();
00058 
00059 virtual int getWidth() const {return width;}                            // use of abstract is better for dynamic block size,...
00060 virtual int getHeight() const {return width;}
00061 
00062 virtual void setAll(AFloatMatrix* from, bool t=false);
00063 virtual void setAll(float value=0.0);
00064 
00065 // NxN matrix operations  
00066 virtual AFSymMatrix* add_A_B(AFSymMatrix* B, AFSymMatrix* result=0);    // C= A+B
00067 virtual AFSymMatrix* add_B_A(AFSymMatrix* B, AFSymMatrix* result=0)     // C= B+A
00068                 { return B->add_A_B(this, result); }
00069 virtual AFSymMatrix* subst_A_B(AFSymMatrix* B, AFSymMatrix* result=0);  // C= A-B
00070 virtual AFSymMatrix* subst_B_A(AFSymMatrix* B, AFSymMatrix* result=0);  // C= B-A
00071 virtual AFSymMatrix* mult_A_B(AFSymMatrix* B, AFSymMatrix* result=0);   // C= A*B
00072 virtual AFSymMatrix* mult_B_A(AFSymMatrix* B, AFSymMatrix* result=0);   // C= B*A
00073 
00074 // A*x=b <==> 
00075 virtual FloatVector* solve(FloatVector* b, FloatVector* x=0){return 0;} // Ax=b, should be possible if SPD, ... 
00076 virtual AFSymMatrix* invert(){return 0;}                                                // A^-1, is easy if SPD
00077 
00078 
00079 virtual AFSymMatrix* jacobiMatrix();                                            // used to compute SOR omega
00080 virtual FloatVector* jacobiVector(FloatVector* b);                              // xk1= BJ*xk + bJ 
00081 
00082 // direct solvers
00083 virtual FloatVector* linearSolverU0(FloatVector* b, FloatVector* x=0);                  // Ux=b
00084 virtual FloatVector* linearSolverU0(FloatVector* b, int k, FloatVector* x);             // U(k)x(k)=b(k)
00085 virtual FloatVector* linearSolverU1(FloatVector* b, FloatVector* x=0);                  // U1x=b
00086 virtual FloatVector* linearSolverL0(FloatVector* b, FloatVector* x=0);                  // Lx=b
00087 virtual FloatVector* linearSolverL1(FloatVector* b, FloatVector* x=0);                  // L1x=b
00088 
00089 
00090 
00091 
00092 virtual AFSymMatrix* inverseSPD(AFSymMatrix* result=0);                         // SPD inversion
00093 
00094 // Combined solvers
00095 virtual FloatVector* linearSolverL1U0(FloatVector* b, FloatVector* x=0);        // crout LU
00096 virtual FloatVector* linearSolverL0U1(FloatVector* b, FloatVector* x=0);        
00097 
00098 virtual FloatVector* linearSolverLLt(FloatVector* b, FloatVector* x=0);         // cholesky solvers
00099 virtual FloatVector* linearSolverUUt(FloatVector* b, FloatVector* x=0);
00100 
00101 // Cholesky sans racine carree
00102 virtual FloatVector* linearSolverL1DL1t(FloatVector* b, FloatVector* x=0);
00103 virtual FloatVector* linearSolverU1DU1t(FloatVector* b, FloatVector* x=0);
00104 
00105 // factorisation complete
00106 virtual AFSymMatrix* factorisationL1U0();               // L1/U0 Doolittle
00107 virtual AFSymMatrix* factorisationL0U1();               // L0/U1 Crout
00108 
00109 virtual AFSymMatrix* factorisationLLt();                // LLt Cholesky
00110 virtual AFSymMatrix* factorisationUUt();                // UUt Cholesky
00111 virtual AFSymMatrix* factorisationL1DL1t();     // L1DL1t Cholesky sans racine carree
00112 virtual AFSymMatrix* factorisationU1DU1t();     // U1DU1t Cholesky sans racine carree
00113 
00114 virtual AFSymMatrix* choleskyUDUt(AFSymMatrix* result=0);               // in place UDUt Cholesky sans racine carree
00115 virtual AFSymMatrix* inverseUDUt(AFSymMatrix* result=0);                // in place UDUt Cholesky sans racine carree
00116 
00117 // incomplete factorisation & preconditioners
00118 // ------------------------------------------
00119 virtual AFloatMatrix* ilu0()
00120                         {printf("%s() %s: not Implemented\n",__FUNCTION__, __FILE__); return 0; }
00121 virtual AFloatMatrix* milu0()
00122                         {printf("%s() %s: not Implemented\n",__FUNCTION__, __FILE__); return 0; }
00123 virtual AFloatMatrix* rilu0(float omega=0.95)
00124                         {printf("%s() %s: not Implemented\n",__FUNCTION__, __FILE__); return 0; }
00125 virtual AFloatMatrix* iluth(double threshold=0.001)
00126                         {printf("%s() %s: not Implemented\n",__FUNCTION__, __FILE__); return 0; }
00127 virtual AFloatMatrix* ilut(int p, double eps)
00128                         {printf("%s() %s: not Implemented\n",__FUNCTION__, __FILE__); return 0;}                // Threshold is related to magitude in row
00129 virtual AFloatMatrix* icf()
00130                         {printf("%s() %s: not Implemented\n",__FUNCTION__, __FILE__); return 0; }
00131 virtual AFloatMatrix* ildlt()
00132                         {printf("%s() %s: not Implemented\n",__FUNCTION__, __FILE__); return 0; }
00133 
00134 virtual FDiagonalMatrix* jacobiPrecond(FDiagonalMatrix* S=0);
00135 
00136 virtual void sor(double omega, AFSymMatrix** R=0, AFSymMatrix**T=0)
00137                         {printf("%s() %s: not Implemented\n",__FUNCTION__, __FILE__);} 
00138 
00139 // new SSOR
00140 virtual AFProductSeqSMatrix* ssorFactorisation();       // allows to select the typed SSOR
00141 
00142 // old SSOR way
00143 virtual void ssor(double omega, AFSymMatrix** R=0, FDiagonalMatrix** S=0, AFSymMatrix**T=0);
00144 virtual void ssorSPD(double omega, AFSymMatrix** R=0, FDiagonalMatrix** S=0);
00145 
00146 // Scaled matrices require diag(A)=I
00147 virtual void ssorScaled(double omega, float* cste, AFSymMatrix** R=0, AFSymMatrix**T=0);
00148 virtual void ssorScaledSPD(double omega, float* cste, AFSymMatrix** R=0);
00149 
00150 // approximate inversion
00151 virtual AFloatMatrix* apinv(int iterOuter, int iterInner){return 0;}            // column based iterative inversion
00152 
00153 // iterative simple linear solvers: Au=b
00154 // -------------------------------------
00155 virtual FloatVector* gaussSeidelLS(int iterMax, FloatVector* b, FloatVector* x=0);
00156 virtual FloatVector* jacobiLS(int iterMax, FloatVector* b, FloatVector* x=0);
00157 virtual FloatVector* sorLS(int iterMax, FloatVector* b, float omega=1.0, FloatVector* x=0);
00158 
00159 // single step iterative linear solvers: Au=b, single iteration
00160 virtual FloatVector* gaussSeidel(FloatVector* b, FloatVector* x=0);
00161 virtual FloatVector* jacobi(FloatVector* b, FloatVector* x=0, FloatVector* dest=0);
00162 virtual FloatVector* sor(FloatVector* b, float omega=1.0, FloatVector* x=0);
00163 
00164 // finding eigenvectors
00165 virtual FloatVector* power(int nmax, float* firstEigenValue, FloatVector* evMax0=0);
00166 virtual FloatVector* inverseIteration(double eigenvalue, int nmax=3, FloatVector* ev=0); 
00167 
00168 // find eigenvalues
00169 virtual float rayleighQuotient(FloatVector* ev=0);
00170 
00171 virtual void output();
00172 virtual void output(FILE* file);
00173 };
00174 
00175 #endif
00176 
00177          
SourceForge.net Logo
Restoreinpaint sourceforge project `C++/Java Image Processing, Restoration, Inpainting Project'.

Bernard De Cuyper: Open Project Leader: Concept, design and development.
Bernard De Cuyper & Eddy Fraiha 2002, 2003. Bernard De Cuyper 2004. Open and free, for friendly usage only.
Modifications on Belgium ground of this piece of artistic work, by governement institutions or companies, must be notified to Bernard De Cuyper.
bern_bdc@hotmail.com