Main Page   Groups   Namespace List   Class Hierarchy   Alphabetical List   Compound List   File List   Namespace Members   Compound Members   File Members   Concepts

itkFEMSolverCrankNicolson.h

Go to the documentation of this file.
00001 /*========================================================================= 00002 00003 Program: Insight Segmentation & Registration Toolkit 00004 Module: $RCSfile: itkFEMSolverCrankNicolson.h,v $ 00005 Language: C++ 00006 Date: $Date: 2003/09/10 14:29:44 $ 00007 Version: $Revision: 1.17 $ 00008 00009 Copyright (c) Insight Software Consortium. All rights reserved. 00010 See ITKCopyright.txt or http://www.itk.org/HTML/Copyright.htm for details. 00011 00012 This software is distributed WITHOUT ANY WARRANTY; without even 00013 the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR 00014 PURPOSE. See the above copyright notices for more information. 00015 00016 =========================================================================*/ 00017 00018 #ifndef __itkFEMSolverCrankNicolson_h 00019 #define __itkFEMSolverCrankNicolson_h 00020 00021 #include "itkFEMSolver.h" 00022 #include "itkFEMElementBase.h" 00023 #include "itkFEMMaterialBase.h" 00024 #include "itkFEMLoadBase.h" 00025 #include "itkFEMLinearSystemWrapperVNL.h" 00026 00027 #include "vnl/vnl_sparse_matrix.h" 00028 #include "vnl/vnl_matrix.h" 00029 #include "vnl/vnl_vector.h" 00030 #include "vnl/algo/vnl_svd.h" 00031 #include "vnl/algo/vnl_cholesky.h" 00032 #include <vnl/vnl_sparse_matrix_linear_system.h> 00033 #include <vnl/algo/vnl_lsqr.h> 00034 #include <math.h> 00035 00036 00037 namespace itk { 00038 namespace fem { 00039 00040 00063 class SolverCrankNicolson : public Solver 00064 { 00065 public: 00066 00067 /* 00068 * helper initialization function before assembly but after generate GFN. 00069 */ 00070 void InitializeForSolution(); 00075 void AssembleKandM(); 00076 00084 void AssembleFforTimeStep(int dim=0); 00085 00089 void Solve(); 00090 00095 void AddToDisplacements(Float optimum=1.0); 00096 void AverageLastTwoDisplacements(Float t=0.5); 00097 void ZeroVector(int which=0); 00098 void PrintDisplacements(); 00099 void PrintForce(); 00100 00102 inline void SetAlpha(Float a = 0.5) { m_alpha=a; } 00103 00105 inline void SetDeltatT(Float T) { m_deltaT=T; } 00106 00108 inline void SetRho(Float rho) { m_rho=rho; } 00109 00113 void RecomputeForceVector(unsigned int index); 00114 00115 /* Finds a triplet that brackets the energy minimum. From Numerical Recipes.*/ 00116 void FindBracketingTriplet(Float* a,Float* b,Float* c); 00120 Float GoldenSection(Float tol=0.01,unsigned int MaxIters=25); 00121 /* Brents method from Numerical Recipes. */ 00122 Float BrentsMethod(Float tol=0.01,unsigned int MaxIters=25); 00123 Float EvaluateResidual(Float t=1.0); 00124 Float GetDeformationEnergy(Float t=1.0); 00125 inline Float GSSign(Float a,Float b) { return (b > 0.0 ? fabs(a) : -1.*fabs(a)); } 00126 inline Float GSMax(Float a,Float b) { return (a > b ? a : b); } 00127 00128 void SetEnergyToMin(Float xmin); 00129 inline LinearSystemWrapper* GetLS(){ return m_ls;} 00130 00131 Float GetCurrentMaxSolution() { return m_CurrentMaxSolution; } 00132 00134 void PrintMinMaxOfSolution(); 00139 SolverCrankNicolson() 00140 { 00141 m_deltaT=0.5; 00142 m_rho=1.; 00143 m_alpha=0.5; 00144 // BUG FIXME NOT SURE IF SOLVER IS USING VECTOR INDEX 1 FOR BCs 00145 ForceTIndex=0; // vector 00146 ForceTMinus1Index=2; // vector 00147 SolutionVectorTMinus1Index=3; // vector 00148 DiffMatrixBySolutionTMinus1Index=4; // vector 00149 ForceTotalIndex=5; // vector 00150 SolutionTIndex=0; // solution 00151 TotalSolutionIndex=1; // solution 00152 SolutionTMinus1Index=2; // solution 00153 SumMatrixIndex=0; // matrix 00154 DifferenceMatrixIndex=1; // matrix 00155 m_CurrentMaxSolution=1.0; 00156 } 00157 00158 00159 ~SolverCrankNicolson() { } 00160 00161 Float m_deltaT; 00162 Float m_rho; 00163 Float m_alpha; 00164 Float m_CurrentMaxSolution; 00165 00166 unsigned int ForceTIndex; 00167 unsigned int ForceTotalIndex; 00168 unsigned int ForceTMinus1Index; 00169 unsigned int SolutionTIndex; 00170 unsigned int SolutionTMinus1Index; 00171 unsigned int SolutionVectorTMinus1Index; 00172 unsigned int TotalSolutionIndex; 00173 unsigned int DifferenceMatrixIndex; 00174 unsigned int SumMatrixIndex; 00175 unsigned int DiffMatrixBySolutionTMinus1Index; 00176 00177 }; 00178 00179 00180 00181 00182 }} // end namespace itk::fem 00183 00184 #endif // #ifndef __itkFEMSolverCrankNicolson_h

Generated at Sat Mar 31 02:18:31 2007 for ITK by doxygen 1.3.8 written by Dimitri van Heesch, © 1997-2000