00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018 #ifndef __VNLIterativeSparseSolverTraits_h
00019 #define __VNLIterativeSparseSolverTraits_h
00020
00021 #include <vnl/vnl_vector.h>
00022 #include <vnl/vnl_sparse_matrix.h>
00023 #include <vnl/vnl_sparse_matrix_linear_system.h>
00024 #include <vnl/algo/vnl_lsqr.h>
00025
00026 template< typename T = double >
00027 class VNLIterativeSparseSolverTraits
00028 {
00029 public:
00030 typedef T ValueType;
00031 typedef vnl_sparse_matrix< ValueType > MatrixType;
00032 typedef vnl_vector< ValueType > VectorType;
00033 typedef vnl_lsqr SolverType;
00034
00035 VNLIterativeSparseSolverTraits() {}
00036 ~VNLIterativeSparseSolverTraits() {}
00037
00038 bool IsDirectSolver() const
00039 {
00040 return false;
00041 }
00042
00043 MatrixType InitializeSparseMatrix( const unsigned int& iN ) const
00044 {
00045 return MatrixType( iN, iN );
00046 }
00047
00048 VectorType InitializeVector( const unsigned int& iN ) const
00049 {
00050 return VectorType( iN );
00051 }
00052
00053 void FillMatrix( MatrixType& iA, const unsigned int& iR, const unsigned int& iC, const ValueType& iV ) const
00054 {
00055 iA( iR, iC ) = iV;
00056 }
00057
00058 void AddToMatrix( MatrixType& iA, const unsigned int& iR, const unsigned int& iC, const ValueType& iV ) const
00059 {
00060 iA( iR, iC ) += iV;
00061 }
00062
00063 bool Solve( const MatrixType& iA, const VectorType& iB, VectorType& oX ) const
00064 {
00065 typedef vnl_sparse_matrix_linear_system< ValueType > SparseLinearSystemType;
00066 SparseLinearSystemType system( iA, iB );
00067
00068 SolverType solver( system );
00069 return solver.minimize( oX );
00070 }
00071
00072
00073 bool Solve( const MatrixType& iA,
00074 const VectorType& iBx, const VectorType& iBy,
00075 VectorType& oX, VectorType& oY ) const
00076 {
00077 bool result1 = Solve( iA, iBx, oX );
00078 bool result2 = Solve( iA, iBy, oY );
00079 return( result1 && result2 );
00080 }
00081
00082 bool Solve( const MatrixType& iA,
00083 const VectorType& iB,
00084 const long& iNbIter,
00085 VectorType& oX ) const
00086 {
00087 typedef vnl_sparse_matrix_linear_system< ValueType > SparseLinearSystemType;
00088 SparseLinearSystemType system( iA, iB );
00089
00090 SolverType solver( system );
00091 solver.set_max_iterations( iNbIter );
00092 return solver.minimize( oX );
00093 }
00094
00095 };
00096
00097 #endif
00098