ITK
4.1.0
Insight Segmentation and Registration Toolkit
|
00001 /*========================================================================= 00002 * 00003 * Copyright Insight Software Consortium 00004 * 00005 * Licensed under the Apache License, Version 2.0 (the "License"); 00006 * you may not use this file except in compliance with the License. 00007 * You may obtain a copy of the License at 00008 * 00009 * http://www.apache.org/licenses/LICENSE-2.0.txt 00010 * 00011 * Unless required by applicable law or agreed to in writing, software 00012 * distributed under the License is distributed on an "AS IS" BASIS, 00013 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. 00014 * See the License for the specific language governing permissions and 00015 * limitations under the License. 00016 * 00017 *=========================================================================*/ 00018 00019 #ifndef __itkFEMLinearSystemWrapperDenseVNL_h 00020 #define __itkFEMLinearSystemWrapperDenseVNL_h 00021 #include "itkFEMLinearSystemWrapper.h" 00022 #include "vnl/vnl_matrix.h" 00023 #include "vnl/vnl_vector.h" 00024 #include "vnl/algo/vnl_svd.h" 00025 #include <vector> 00026 00027 namespace itk 00028 { 00029 namespace fem 00030 { 00038 class LinearSystemWrapperDenseVNL : public LinearSystemWrapper 00039 { 00040 public: 00041 00042 /* values stored in matrices & vectors */ 00043 typedef LinearSystemWrapper::Float Float; 00044 00045 /* superclass */ 00046 typedef LinearSystemWrapper SuperClass; 00047 00048 /* matrix typedef */ 00049 typedef vnl_matrix<Float> MatrixRepresentation; 00050 00051 /* matrix holder typedef */ 00052 typedef std::vector<MatrixRepresentation *> MatrixHolder; 00053 00054 /* constructor & destructor */ 00055 LinearSystemWrapperDenseVNL() : LinearSystemWrapper(), m_Matrices(0), m_Vectors(0), m_Solutions(0) 00056 { 00057 } 00058 virtual ~LinearSystemWrapperDenseVNL(); 00059 00060 /* memory management routines */ 00061 virtual void InitializeMatrix(unsigned int matrixIndex); 00062 00063 virtual bool IsMatrixInitialized(unsigned int matrixIndex); 00064 00065 virtual void DestroyMatrix(unsigned int matrixIndex); 00066 00067 virtual void InitializeVector(unsigned int vectorIndex); 00068 00069 virtual bool IsVectorInitialized(unsigned int vectorIndex); 00070 00071 virtual void DestroyVector(unsigned int vectorIndex); 00072 00073 virtual void InitializeSolution(unsigned int solutionIndex); 00074 00075 virtual bool IsSolutionInitialized(unsigned int solutionIndex); 00076 00077 virtual void DestroySolution(unsigned int solutionIndex); 00078 00079 virtual void SetMaximumNonZeroValuesInMatrix(unsigned int, unsigned int) 00080 { 00081 } 00082 00083 /* assembly & solving routines */ 00084 virtual Float GetMatrixValue(unsigned int i, unsigned int j, 00085 unsigned int matrixIndex) const 00086 { 00087 return ( *( ( *m_Matrices )[matrixIndex] ) )(i, j); 00088 } 00089 virtual void SetMatrixValue(unsigned int i, unsigned int j, Float value, 00090 unsigned int matrixIndex) 00091 { 00092 ( *( ( *m_Matrices )[matrixIndex] ) )(i, j) = value; 00093 } 00094 virtual void AddMatrixValue(unsigned int i, unsigned int j, Float value, 00095 unsigned int matrixIndex) 00096 { 00097 ( *( ( *m_Matrices )[matrixIndex] ) )(i, j) += value; 00098 } 00099 virtual Float GetVectorValue(unsigned int i, 00100 unsigned int vectorIndex) const 00101 { 00102 return ( *( ( *m_Vectors )[vectorIndex] ) )[i]; 00103 } 00104 virtual void SetVectorValue(unsigned int i, Float value, 00105 unsigned int vectorIndex) 00106 { 00107 ( *( ( *m_Vectors )[vectorIndex] ) )(i) = value; 00108 } 00109 virtual void AddVectorValue(unsigned int i, Float value, 00110 unsigned int vectorIndex) 00111 { 00112 ( *( ( *m_Vectors )[vectorIndex] ) )(i) += value; 00113 } 00114 virtual Float GetSolutionValue(unsigned int i, unsigned int solutionIndex) const; 00115 00116 virtual void SetSolutionValue(unsigned int i, Float value, 00117 unsigned int solutionIndex) 00118 { 00119 ( *( ( *m_Solutions )[solutionIndex] ) )(i) = value; 00120 } 00121 virtual void AddSolutionValue(unsigned int i, Float value, 00122 unsigned int solutionIndex) 00123 { 00124 ( *( ( *m_Solutions )[solutionIndex] ) )(i) += value; 00125 } 00126 virtual void Solve(void); 00127 00128 /* matrix & vector manipulation routines */ 00129 virtual void ScaleMatrix(Float scale, unsigned int matrixIndex); 00130 00131 virtual void ScaleVector(Float scale, unsigned int vectorIndex); 00132 00133 virtual void ScaleSolution(Float scale, unsigned int solutionIndex); 00134 00135 virtual void SwapMatrices(unsigned int matrixIndex1, unsigned int matrixIndex2); 00136 00137 virtual void SwapVectors(unsigned int vectorIndex1, unsigned int vectorIndex2); 00138 00139 virtual void SwapSolutions(unsigned int solutionIndex1, unsigned int solutionIndex2); 00140 00141 virtual void CopySolution2Vector(unsigned solutionIndex, unsigned int vectorIndex); 00142 00143 virtual void CopyVector2Solution(unsigned int vectorIndex, unsigned int solutionIndex); 00144 00145 virtual void MultiplyMatrixMatrix(unsigned int resultMatrixIndex, unsigned int leftMatrixIndex, 00146 unsigned int rightMatrixIndex); 00147 00148 virtual void MultiplyMatrixVector(unsigned int resultVectorIndex, unsigned int matrixIndex, unsigned int vectorIndex); 00149 00150 private: 00151 00153 // std::vector< vnl_sparse_matrix<Float>* > *m_Matrices; 00154 MatrixHolder *m_Matrices; 00155 00157 std::vector<vnl_vector<Float> *> *m_Vectors; 00158 00160 std::vector<vnl_vector<Float> *> *m_Solutions; 00161 }; 00162 } 00163 } // end namespace itk::fem 00164 00165 #endif 00166