ITK  4.8.0
Insight Segmentation and Registration Toolkit
itkFEMElementBase.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 itkFEMElementBase_h
20 #define itkFEMElementBase_h
21 
22 #include "itkFEMLightObject.h"
23 #include "itkFEMPArray.h"
24 #include "itkFEMMaterialBase.h"
25 #include "itkFEMSolution.h"
26 #include "itkVectorContainer.h"
27 #include "vnl/vnl_matrix.h"
28 #include "vnl/vnl_vector.h"
29 #include "ITKFEMExport.h"
30 
31 #include <set>
32 #include <vector>
33 
34 namespace itk
35 {
36 namespace fem
37 {
73 class ITKFEM_EXPORT Element : public FEMLightObject
74 {
75 public:
77  typedef Element Self;
81 
83  itkTypeMacro(Element, FEMLightObject);
84 
88  typedef double Float;
89  typedef unsigned long ElementIdentifier;
90 
94  // FIXME - Remove FEMPArray Type and replace with VectorContainer version
97 
101  typedef vnl_matrix<Float> MatrixType;
102 
106  typedef vnl_vector<Float> VectorType;
107 
122 
127  typedef unsigned int DegreeOfFreedomIDType;
128 
134  enum { InvalidDegreeOfFreedomID = 0xffffffff };
135 
136 
148  class ITKFEM_EXPORT Node : public FEMLightObject
149  {
150  public:
151 
153  typedef Node Self;
157 
159  // itkNewMacro(Self);
160  static Pointer New(void)
161  {
164 
165  if( smartPtr.IsNull() )
166  {
167  smartPtr = static_cast<Pointer>(new Self);
168  }
169  smartPtr->UnRegister();
170  return smartPtr;
171  }
172 
174  itkTypeMacro(Node, FEMLightObject);
175 
178  virtual::itk::LightObject::Pointer CreateAnother(void) const ITK_OVERRIDE;
179 
183  typedef double Float;
184 
189 
194  {
195  }
196 
201  {
202  this->ClearDegreesOfFreedom();
203  this->m_elements.clear();
204  }
206 
211  const VectorType & GetCoordinates(void) const
212  {
213  return m_coordinates;
214  }
215 
219  void SetCoordinates(const VectorType & coords)
220  {
221  m_coordinates = coords;
222  }
223 
228  {
229  if( i >= m_dof.size() )
230  {
231  return InvalidDegreeOfFreedomID;
232  }
233  return m_dof[i];
234  }
236 
240  void SetDegreeOfFreedom(unsigned int i, DegreeOfFreedomIDType dof) const
241  {
242  if( i >= m_dof.size() )
243  {
244  m_dof.resize(i + 1, InvalidDegreeOfFreedomID);
245  }
246  m_dof[i] = dof;
247  }
249 
250  virtual void ClearDegreesOfFreedom(void) const;
251 
252  public:
257  typedef std::set<Element *> SetOfElements;
259  protected:
260  virtual void PrintSelf(std::ostream& os, Indent indent) const ITK_OVERRIDE;
261  private:
262 
267 
272  mutable std::vector<DegreeOfFreedomIDType> m_dof;
273  }; // end class Node
274 
275 // ////////////////////////////////////////////////////////////////////////
276 /*
277  * Methods related to the physics of the problem.
278  */
279 
280  virtual VectorType GetStrainsAtPoint(const VectorType & pt, const Solution & sol, unsigned int index) const;
281 
282  virtual VectorType GetStressesAtPoint(const VectorType & pt, const VectorType & e, const Solution & sol,
283  unsigned int index) const;
284 
305  virtual void GetStiffnessMatrix(MatrixType & Ke) const;
306 
316  virtual Float GetElementDeformationEnergy(MatrixType & LocalSolution) const;
317 
329  virtual void GetMassMatrix(MatrixType & Me) const;
330 
342  virtual void GetLandmarkContributionMatrix(float eta, MatrixType & Le) const;
343 
351  virtual void GetStrainDisplacementMatrix(MatrixType & B, const MatrixType & shapeDgl) const = 0;
352 
358  virtual void GetMaterialMatrix(MatrixType & D) const = 0;
359 
371  virtual VectorType InterpolateSolution(const VectorType & pt,
372  const Solution & sol,
373  unsigned int solutionIndex = 0) const;
374 
388  virtual Float InterpolateSolutionN(const VectorType & pt, const Solution & sol, unsigned int f,
389  unsigned int solutionIndex = 0) const;
390 
397  DegreeOfFreedomIDType GetDegreeOfFreedom(unsigned int local_dof) const
398  {
399  if( local_dof > this->GetNumberOfDegreesOfFreedom() )
400  {
401  return InvalidDegreeOfFreedomID;
402  }
403  return this->GetNode(local_dof /
404  this->GetNumberOfDegreesOfFreedomPerNode() )
405  ->GetDegreeOfFreedom(local_dof % this->GetNumberOfDegreesOfFreedomPerNode() );
406  }
408 
425  virtual Material::ConstPointer GetMaterial(void) const;
426 
435  virtual void SetMaterial(Material::ConstPointer);
436 
437  // ////////////////////////////////////////////////////////////////////////
464  virtual void GetIntegrationPointAndWeight(unsigned int i,
465  VectorType & pt,
466  Float & w,
467  unsigned int order = 0) const = 0;
468 
479  virtual unsigned int GetNumberOfIntegrationPoints(unsigned int order = 0) const = 0;
480 
489  enum { gaussMaxOrder = 10 };
490 
502  static const Float gaussPoint[gaussMaxOrder + 1][gaussMaxOrder];
503 
509  static const Float gaussWeight[gaussMaxOrder + 1][gaussMaxOrder];
510 
511 // ////////////////////////////////////////////////////////////////////////
512 /*
513  * Methods related to the geometry of an element
514  */
515 
521 
525  virtual unsigned int GetNumberOfNodes(void) const = 0;
526 
530  virtual NodeIDType GetNode(unsigned int n) const = 0;
531 
535  virtual void SetNode(unsigned int n, NodeIDType node) = 0;
536  virtual void SetNode(unsigned int n, Node::Pointer node);
537 
543  virtual const VectorType & GetNodeCoordinates(unsigned int n) const = 0;
544 
550  virtual VectorType GetGlobalFromLocalCoordinates(const VectorType & pt) const;
551 
558  virtual bool GetLocalFromGlobalCoordinates(const VectorType & globalPt, VectorType & localPt) const = 0;
559 
565  virtual unsigned int GetNumberOfSpatialDimensions() const = 0;
566 
574  virtual VectorType ShapeFunctions(const VectorType & pt) const = 0;
575 
591  virtual void ShapeFunctionDerivatives(const VectorType & pt, MatrixType & shapeD) const = 0;
592 
612  virtual void ShapeFunctionGlobalDerivatives(const VectorType & pt, MatrixType & shapeDgl, const MatrixType *pJ = ITK_NULLPTR,
613  const MatrixType *pshapeD = ITK_NULLPTR) const;
614 
636  virtual void Jacobian(const VectorType & pt, MatrixType & J, const MatrixType *pshapeD = ITK_NULLPTR) const;
637 
647  virtual Float JacobianDeterminant(const VectorType & pt, const MatrixType *pJ = ITK_NULLPTR) const;
648 
660  virtual void JacobianInverse(const VectorType & pt, MatrixType & invJ, const MatrixType *pJ = ITK_NULLPTR) const;
661 
667  virtual unsigned int GetNumberOfDegreesOfFreedom(void) const;
668 
672  virtual std::vector<std::vector<int> > GetEdgeIds(void) const;
673 
681  virtual unsigned int GetNumberOfDegreesOfFreedomPerNode(void) const = 0;
682 
684  virtual void PopulateEdgeIds(void) = 0;
685 
686 protected:
687 
688  // to store edge connectivity data
689  std::vector<std::vector<int> > m_EdgeIds;
690 
691  virtual void PrintSelf(std::ostream& os, Indent indent) const ITK_OVERRIDE;
692 
693 };
694 
695 }
696 } // end namespace itk::fem
697 
698 #endif // #ifndef itkFEMElementBase_h
DegreeOfFreedomIDType GetDegreeOfFreedom(unsigned int i) const
DegreeOfFreedomIDType GetDegreeOfFreedom(unsigned int local_dof) const
static Pointer New(void)
Light weight base class for most itk classes.
LoadType::Pointer LoadPointer
FEMLightObject LoadType
std::vector< DegreeOfFreedomIDType > m_dof
unsigned long ElementIdentifier
std::vector< std::vector< int > > m_EdgeIds
std::set< Element * > SetOfElements
SmartPointer< const Self > ConstPointer
static T::Pointer Create()
void SetCoordinates(const VectorType &coords)
static const double e
The base of the natural logarithm or Euler&#39;s number
Definition: itkMath.h:45
Class that stores information required to define a node.
void UnRegister() ITK_SP_NOEXCEPT
vnl_matrix< Float > MatrixType
unsigned int DegreeOfFreedomIDType
Provides functions to access the values of the solution vector.
VectorContainer< ElementIdentifier, Element::Pointer > ArrayType1
Base class for all classes that define the FEM system.
vnl_vector< Float > VectorType
SmartPointer< Self > Pointer
SmartPointer< Self > Pointer
bool IsNull() const
const VectorType & GetCoordinates(void) const
Array for FEMP objects.
Definition: itkFEMPArray.h:40
Abstract base element class.
Node::ConstPointer NodeIDType
Define a front-end to the STL &quot;vector&quot; container that conforms to the IndexedContainerInterface.
SmartPointer< const Self > ConstPointer
Control indentation during Print() invocation.
Definition: itkIndent.h:49
FEMPArray< Element > ArrayType
void SetDegreeOfFreedom(unsigned int i, DegreeOfFreedomIDType dof) const
FEMLightObject Superclass