Go to the documentation of this file.00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018 #ifndef __itkFEMSolver_h
00019 #define __itkFEMSolver_h
00020
00021 #include "itkFEMElementBase.h"
00022 #include "itkFEMMaterialBase.h"
00023 #include "itkFEMLoadBase.h"
00024
00025 #include "itkFEMLinearSystemWrapper.h"
00026 #include "itkFEMLinearSystemWrapperVNL.h"
00027
00028 #include "itkImage.h"
00029
00030 namespace itk {
00031 namespace fem {
00032
00041 class Solver
00042 {
00043 public:
00044
00048 typedef Element::Float Float;
00049
00054 typedef Element::ArrayType ElementArray;
00055 ElementArray el;
00056
00060 typedef Node::ArrayType NodeArray;
00061 NodeArray node;
00062
00066 typedef Load::ArrayType LoadArray;
00067 LoadArray load;
00068
00072 typedef Material::ArrayType MaterialArray;
00073 MaterialArray mat;
00074
00078 typedef Element::VectorType VectorType;
00079
00090 itkStaticConstMacro(MaxGridDimensions, unsigned int, 3);
00091
00095 typedef itk::Image<Element::ConstPointer,MaxGridDimensions> InterpolationGridType;
00096
00111 void InitializeInterpolationGrid(const VectorType& size, const VectorType& bb1, const VectorType& bb2);
00112
00116 void InitializeInterpolationGrid(const VectorType& size)
00117 {
00118 InitializeInterpolationGrid(size, VectorType(size.size(),0.0), size-1.0);
00119 }
00121
00132 const InterpolationGridType * GetInterpolationGrid(void) const
00133 { return m_InterpolationGrid.GetPointer(); }
00134
00143 const Element * GetElementAtPoint(const VectorType& pt) const;
00144
00148 void Read( std::istream& f );
00149
00153 void Write( std::ostream& f );
00154
00158 virtual void Clear( void );
00159
00160
00169 void GenerateGFN( void );
00170
00174 void AssembleK( void );
00175
00182 virtual void InitializeMatrixForAssembly(unsigned int N);
00183
00189 virtual void FinalizeMatrixAfterAssembly( void )
00190 {
00191
00192 this->ApplyBC();
00193 }
00194
00201 virtual void AssembleElementMatrix(Element::Pointer e);
00202
00210 virtual void AssembleLandmarkContribution(Element::Pointer e, float);
00211
00225 void ApplyBC(int dim=0, unsigned int matrix=0);
00226
00234 void AssembleF(int dim=0);
00235
00239 void DecomposeK( void );
00240
00244 virtual void Solve( void );
00245
00250 void UpdateDisplacements( void );
00251
00252 Float GetSolution(unsigned int i,unsigned int which=0)
00253 {
00254 return m_ls->GetSolutionValue(i,which);
00255 }
00256
00257 unsigned int GetNumberOfDegreesOfFreedom( void )
00258 {
00259 return NGFN;
00260 }
00261
00263 Float GetDeformationEnergy(unsigned int SolutionIndex=0);
00264
00265 public:
00270 Solver();
00271
00275 virtual ~Solver() {}
00276
00291 void SetLinearSystemWrapper(LinearSystemWrapper::Pointer ls);
00292
00298 LinearSystemWrapper::Pointer GetLinearSystemWrapper() { return m_ls; }
00299
00304 virtual void InitializeLinearSystemWrapper(void);
00305
00309 virtual Float GetTimeStep( void ) const { return 0.0; }
00310
00316 virtual void SetTimeStep(Float) {}
00317
00318 protected:
00319
00323 unsigned int NGFN;
00324
00329 unsigned int NMFC;
00330
00332 LinearSystemWrapper::Pointer m_ls;
00333
00334 private:
00335
00339 LinearSystemWrapperVNL m_lsVNL;
00340
00346 InterpolationGridType::Pointer m_InterpolationGrid;
00347
00348 };
00349
00350 }}
00351
00352 #endif // #ifndef __itkFEMSolver_h
00353