ITK  4.13.0
Insight Segmentation and Registration Toolkit
itkFEMSolver.h
Go to the documentation of this file.
1 /*=========================================================================
2  *
3  * Copyright Insight Software Consortium
4  *
5  * Licensed under the Apache License, Version 2.0 (the "License");
6  * you may not use this file except in compliance with the License.
7  * You may obtain a copy of the License at
8  *
9  * http://www.apache.org/licenses/LICENSE-2.0.txt
10  *
11  * Unless required by applicable law or agreed to in writing, software
12  * distributed under the License is distributed on an "AS IS" BASIS,
13  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14  * See the License for the specific language governing permissions and
15  * limitations under the License.
16  *
17  *=========================================================================*/
18 
19 #ifndef itkFEMSolver_h
20 #define itkFEMSolver_h
21 
22 #include "itkProcessObject.h"
23 #include "itkFEMObject.h"
24 
27 
28 #include "itkImage.h"
29 
30 namespace itk
31 {
32 namespace fem
33 {
71 
72 template <unsigned int VDimension = 3>
73 class ITK_TEMPLATE_EXPORT Solver : public ProcessObject
74 {
75 public:
77  typedef Solver Self;
81 
83  itkNewMacro(Self);
84 
86  itkTypeMacro(Solver, ProcessObject);
87 
88  itkStaticConstMacro(FEMDimension, unsigned int, VDimension);
89  itkStaticConstMacro(MaxDimensions, unsigned int, 3);
90 
96 
104 
116 
118  itkSetMacro(Origin, InterpolationGridPointType);
119  itkGetMacro(Origin, InterpolationGridPointType);
121 
123  itkSetMacro(Spacing, InterpolationGridSpacingType);
124  itkGetMacro(Spacing, InterpolationGridSpacingType);
126 
128  itkSetMacro(Region, InterpolationGridRegionType);
129  itkGetMacro(Region, InterpolationGridRegionType);
131 
133  itkSetMacro(Direction, InterpolationGridDirectionType);
134  itkGetMacro(Direction, InterpolationGridDirectionType);
136 
138  virtual Float GetTimeStep() const;
139 
145  virtual void SetTimeStep(Float dt);
146 
148  Float GetSolution(unsigned int i, unsigned int which = 0);
149 
152  using Superclass::SetInput;
153  virtual void SetInput( FEMObjectType *fem);
154 
155  virtual void SetInput( unsigned int, FEMObjectType * fem);
156 
157  FEMObjectType * GetInput();
158 
159  FEMObjectType * GetInput(unsigned int idx);
160 
169  const Element * GetElementAtPoint(const VectorType & pt) const;
170 
172  Float GetDeformationEnergy(unsigned int SolutionIndex = 0);
173 
188  void SetLinearSystemWrapper(LinearSystemWrapper::Pointer ls);
189 
196  {
197  return m_LinearSystem;
198  }
199 
214  void InitializeInterpolationGrid(const InterpolationGridSizeType & size, const InterpolationGridPointType & bb1,
215  const InterpolationGridPointType & bb2);
216 
219  {
222 
223  bb1.Fill(0.0);
224 
226  for( unsigned int i = 0; i < FEMDimension; i++ )
227  {
228  bb2[i] = size[i] - 1.0;
229  }
230  InitializeInterpolationGrid(size, bb1, bb2);
231  }
232 
235  void InitializeInterpolationGrid(const InterpolationGridRegionType & region,
236  const InterpolationGridPointType & origin,
237  const InterpolationGridSpacingType & spacing,
238  const InterpolationGridDirectionType & direction);
239 
251  {
252  return m_InterpolationGrid.GetPointer();
253  }
254 
255 
259  using Superclass::MakeOutput;
260  virtual DataObjectPointer MakeOutput(DataObjectPointerArraySizeType) ITK_OVERRIDE;
261 
275  FEMObjectType * GetOutput();
276 
277  FEMObjectType * GetOutput(unsigned int idx);
278 
279 protected:
280  Solver();
281  virtual ~Solver() ITK_OVERRIDE;
282  virtual void PrintSelf(std::ostream& os, Indent indent) const ITK_OVERRIDE;
283 
286  virtual void GenerateData() ITK_OVERRIDE;
287 
288 
297  // void GenerateGFN();
298 
300  void AssembleK();
301 
308  virtual void InitializeMatrixForAssembly(unsigned int N);
309 
315  virtual void FinalizeMatrixAfterAssembly(void)
316  {
317  // Apply the boundary conditions to the K matrix
318  this->ApplyBC();
319  }
320 
327  virtual void AssembleElementMatrix(Element::Pointer e);
328 
336  virtual void AssembleLandmarkContribution(Element::ConstPointer e, float);
337 
351  void ApplyBC(int dim = 0, unsigned int matrix = 0);
352 
360  void AssembleF(int dim = 0);
361 
363  void DecomposeK();
364 
367  virtual void RunSolver();
368 
373  void UpdateDisplacements();
374 
376  void FillInterpolationGrid();
377 
382  virtual void InitializeLinearSystemWrapper();
383 
385  unsigned int m_NGFN;
386 
391  unsigned int m_NMFC;
392 
395 
398 
405 
407 
408 private:
409  ITK_DISALLOW_COPY_AND_ASSIGN(Solver);
410 
416 
417 };
418 } // end namespace fem
419 } // end namespace itk
420 
421 #ifndef ITK_MANUAL_INSTANTIATION
422 #include "itkFEMSolver.hxx"
423 #endif
424 
425 #endif // #ifndef itkFEMSolver_h
Superclass::RegionType RegionType
Definition: itkImage.h:137
unsigned int m_NMFC
Definition: itkFEMSolver.h:391
void InitializeInterpolationGrid(const InterpolationGridSizeType &size)
Definition: itkFEMSolver.h:218
Light weight base class for most itk classes.
InterpolationGridType::PointType InterpolationGridPointType
Definition: itkFEMSolver.h:112
Implements N-dimensional Finite element (FE) models including elements, materials, and loads.
Definition: itkFEMObject.h:75
Represent the size (bounds) of a n-dimensional image.
Definition: itkSize.h:52
InterpolationGridType::SpacingType InterpolationGridSpacingType
Definition: itkFEMSolver.h:113
InterpolationGridType::DirectionType InterpolationGridDirectionType
Definition: itkFEMSolver.h:115
unsigned int m_NGFN
Definition: itkFEMSolver.h:385
InterpolationGridPointerType m_InterpolationGrid
Definition: itkFEMSolver.h:404
InterpolationGridSpacingType m_Spacing
Definition: itkFEMSolver.h:414
Element::VectorType VectorType
Definition: itkFEMSolver.h:99
FEMObjectType::ConstPointer FEMObjectConstPointer
Definition: itkFEMSolver.h:94
FEMObjectType::Pointer FEMObjectPointer
Definition: itkFEMSolver.h:93
InterpolationGridDirectionType m_Direction
Definition: itkFEMSolver.h:415
The base class for all process objects (source, filters, mappers) in the Insight data processing pipe...
itk::Image< Element::ConstPointer, VDimension > InterpolationGridType
Definition: itkFEMSolver.h:108
void Fill(const ValueType &)
itk::fem::FEMObject< VDimension > FEMObjectType
Definition: itkFEMSolver.h:92
Load::ArrayType LoadArray
Definition: itkFEMSolver.h:102
SmartPointer< Self > Pointer
Definition: itkFEMSolver.h:79
InterpolationGridType::RegionType InterpolationGridRegionType
Definition: itkFEMSolver.h:111
FEMObjectPointer m_FEMObject
Definition: itkFEMSolver.h:406
Defines all functions required by Solver class to allocate, assemble and solve a linear system of equ...
InterpolationGridType::IndexType InterpolationGridIndexType
Definition: itkFEMSolver.h:114
LinearSystemWrapper::Pointer m_LinearSystem
Definition: itkFEMSolver.h:394
LinearSystemWrapper::Pointer GetLinearSystemWrapper()
Definition: itkFEMSolver.h:195
DataObjectPointerArray::size_type DataObjectPointerArraySizeType
vnl_vector< Float > VectorType
InterpolationGridRegionType m_Region
Definition: itkFEMSolver.h:409
Element::ArrayType ElementArray
Definition: itkFEMSolver.h:101
InterpolationGridType::Pointer InterpolationGridPointerType
Definition: itkFEMSolver.h:109
InterpolationGridPointType m_Origin
Definition: itkFEMSolver.h:413
FEM solver used to generate a solution for a FE formulation.
Definition: itkFEMSolver.h:73
A region represents some portion or piece of data.
Definition: itkRegion.h:64
Element::Node::ArrayType NodeArray
Definition: itkFEMSolver.h:100
Array for FEMP objects.
Definition: itkFEMPArray.h:40
ProcessObject Superclass
Definition: itkFEMSolver.h:78
Abstract base element class.
InterpolationGridType::SizeType InterpolationGridSizeType
Definition: itkFEMSolver.h:110
ProcessObject::DataObjectPointerArraySizeType DataObjectPointerArraySizeType
Definition: itkFEMSolver.h:258
DataObject::Pointer DataObjectPointer
Definition: itkFEMSolver.h:95
const InterpolationGridType * GetInterpolationGrid(void) const
Definition: itkFEMSolver.h:250
Control indentation during Print() invocation.
Definition: itkIndent.h:49
Element::Float Float
Definition: itkFEMSolver.h:98
static ITK_CONSTEXPR_VAR double e
The base of the natural logarithm or Euler&#39;s number
Definition: itkMath.h:56
SmartPointer< const Self > ConstPointer
Definition: itkFEMSolver.h:80
LinearSystemWrapper class that uses VNL numeric library functions to define a sparse linear system of...
Templated n-dimensional image class.
Definition: itkImage.h:75
LinearSystemWrapperVNL m_LinearSystemVNL
Definition: itkFEMSolver.h:397
Material::ArrayType MaterialArray
Definition: itkFEMSolver.h:103