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

itkFEMElementBase.h

Go to the documentation of this file.
00001 /*=========================================================================
00002 
00003   Program:   Insight Segmentation & Registration Toolkit
00004   Module:    $RCSfile: itkFEMElementBase.h,v $
00005   Language:  C++
00006   Date:      $Date: 2002/10/01 21:49:54 $
00007   Version:   $Revision: 1.37 $
00008 
00009   Copyright (c) 2002 Insight 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 __itkFEMElementBase_h
00019 #define __itkFEMElementBase_h
00020 
00021 #include "itkFEMLightObject.h"
00022 #include "itkFEMPArray.h"
00023 #include "itkFEMMaterialBase.h"
00024 #include "itkFEMSolution.h"
00025 #include "itkVisitorDispatcher.h"
00026 #include "vnl/vnl_matrix.h"
00027 #include "vnl/vnl_vector.h"
00028 #include <set>
00029 #include <vector>
00030 
00031 namespace itk {
00032 namespace fem {
00033 
00034 
00035 
00036 // FIXME: Write better documentation
00069 #define HANDLE_ELEMENT_LOADS() \
00070  \
00071   typedef void (*LoadImplementationFunctionPointer)(ConstPointer,Element::LoadPointer, Element::VectorType& ); \
00072   virtual void GetLoadVector( Element::LoadPointer l, Element::VectorType& Fe ) const \
00073   { VisitorDispatcher<Self,Element::LoadType, LoadImplementationFunctionPointer>::Visit(l)(this,l,Fe); }
00074 
00075 class Element : public FEMLightObject
00076 {
00077 FEM_ABSTRACT_CLASS(Element,FEMLightObject)
00078 public:
00079 
00083   typedef double Float;
00084 
00088   typedef FEMPArray<Element> ArrayType;
00089 
00093   typedef vnl_matrix<Float> MatrixType;
00094 
00098   typedef vnl_vector<Float> VectorType;
00099 
00111   typedef FEMLightObject LoadType;
00112   typedef LoadType::Pointer LoadPointer;
00113 
00117   typedef unsigned int DegreeOfFreedomIDType;
00118 
00124   enum{ InvalidDegreeOfFreedomID = 0xffffffff };
00125 
00134   class Node : public FEMLightObject
00135   {
00136   FEM_CLASS(Node,FEMLightObject)
00137   public:
00138 
00142     typedef double Float;
00143 
00147     typedef FEMPArray<Self> ArrayType;
00148 
00149 
00150     /* Windows visualization */
00151   #ifdef FEM_BUILD_VISUALIZATION
00152 
00153     void Draw(CDC* pDC, Solution::ConstPointer sol) const;
00155     static double& DC_Scale;
00156   #endif
00157 
00161     Node() {}
00162 
00166     Node(Float x, Float y) : m_coordinates(VectorType(2))
00167     { m_coordinates[0]=x; m_coordinates[1]=y; }
00168 
00172     Node(Float x, Float y, Float z) : m_coordinates(VectorType(x,y,z)) {}
00173 
00178     const VectorType& GetCoordinates( void ) const
00179     { return m_coordinates; }
00180 
00184     void SetCoordinates( const VectorType& coords )
00185     { m_coordinates=coords; }
00186 
00190     DegreeOfFreedomIDType GetDegreeOfFreedom(unsigned int i) const
00191     {
00192       if( i>=m_dof.size() ) { return InvalidDegreeOfFreedomID; }
00193       return m_dof[i];
00194     }
00195 
00199     void SetDegreeOfFreedom(unsigned int i, DegreeOfFreedomIDType dof) const
00200     {
00201       if( i>=m_dof.size() ) { m_dof.resize(i+1, InvalidDegreeOfFreedomID); }
00202       m_dof[i]=dof;
00203     }
00204 
00205     virtual void ClearDegreesOfFreedom( void ) const
00206     {
00207       m_dof.clear();
00208     }
00209 
00210     virtual void Read(  std::istream& f, void* info );
00211     virtual void Write( std::ostream& f ) const;
00212 
00213   public:
00218     typedef std::set<Element*> SetOfElements;
00219     mutable SetOfElements m_elements;
00220 
00221   private:
00225     VectorType m_coordinates;
00226 
00231     mutable std::vector<DegreeOfFreedomIDType> m_dof;
00232 
00233   }; // end class Node
00234 
00235 
00236 
00237 
00239   /*
00240    * Methods related to the physics of the problem.
00241    */
00242 
00263   virtual void GetStiffnessMatrix( MatrixType& Ke ) const;
00264 
00274   virtual Float GetElementDeformationEnergy( MatrixType& LocalSolution ) const;
00275 
00287   virtual void GetMassMatrix( MatrixType& Me ) const;
00288 
00318   virtual void GetLoadVector( LoadPointer l, VectorType& Fe ) const = 0;
00319 
00327   virtual void GetStrainDisplacementMatrix( MatrixType& B, const MatrixType& shapeDgl ) const = 0;
00328 
00334   virtual void GetMaterialMatrix( MatrixType& D ) const = 0;
00335 
00347   virtual VectorType InterpolateSolution( const VectorType& pt, const Solution& sol , unsigned int solutionIndex=0 ) const;
00348 
00362   virtual Float InterpolateSolutionN( const VectorType& pt, const Solution& sol, unsigned int f , unsigned int solutionIndex=0 ) const;
00363 
00370   DegreeOfFreedomIDType GetDegreeOfFreedom( unsigned int local_dof ) const
00371   {
00372     if(local_dof>this->GetNumberOfDegreesOfFreedom()) { return InvalidDegreeOfFreedomID; }
00373     return this->GetNode(local_dof/this->GetNumberOfDegreesOfFreedomPerNode())->GetDegreeOfFreedom(local_dof%this->GetNumberOfDegreesOfFreedomPerNode());
00374   }
00375 
00392   virtual Material::ConstPointer GetMaterial(void) const { return 0; }
00393 
00402   virtual void SetMaterial(Material::ConstPointer) {} // FIXME: maybe we should throw an exception instead
00403 
00404 
00405 
00406 
00408   /*
00409    * Methods related to numeric integration
00410    */
00411 
00434   virtual void GetIntegrationPointAndWeight( unsigned int i, VectorType& pt, Float& w, unsigned int order=0 ) const = 0;
00435 
00446   virtual unsigned int GetNumberOfIntegrationPoints( unsigned int order=0 ) const = 0;
00447 
00456   enum { gaussMaxOrder=10 };
00457 
00469   static const Float gaussPoint[gaussMaxOrder+1][gaussMaxOrder];
00470 
00476   static const Float gaussWeight[gaussMaxOrder+1][gaussMaxOrder];
00477 
00478 
00480   /*
00481    * Methods related to the geometry of an element
00482    */
00483 
00488   typedef Node::ConstPointer NodeIDType;
00489 
00493   virtual unsigned int GetNumberOfNodes( void ) const = 0;
00494 
00498   virtual NodeIDType GetNode(unsigned int n) const = 0;
00499 
00503   virtual void SetNode(unsigned int n, NodeIDType node) = 0;
00504 
00510   virtual const VectorType& GetNodeCoordinates( unsigned int n ) const = 0;
00511 
00517   virtual VectorType GetGlobalFromLocalCoordinates( const VectorType& pt ) const;
00518 
00525   virtual bool GetLocalFromGlobalCoordinates( const VectorType& globalPt , VectorType& localPt ) const = 0;
00526 
00532   virtual unsigned int GetNumberOfSpatialDimensions() const = 0;
00533 
00541   virtual VectorType ShapeFunctions( const VectorType& pt ) const = 0;
00542 
00558   virtual void ShapeFunctionDerivatives( const VectorType& pt, MatrixType& shapeD ) const = 0;
00559 
00579   virtual void ShapeFunctionGlobalDerivatives( const VectorType& pt, MatrixType& shapeDgl, const MatrixType* pJ=0, const MatrixType* pshapeD=0 ) const;
00580 
00602   virtual void Jacobian( const VectorType& pt, MatrixType& J, const MatrixType* pshapeD = 0 ) const;
00603 
00613   virtual Float JacobianDeterminant( const VectorType& pt, const MatrixType* pJ = 0 ) const;
00614 
00626   virtual void JacobianInverse( const VectorType& pt, MatrixType& invJ, const MatrixType* pJ = 0 ) const;
00627 
00633   virtual unsigned int GetNumberOfDegreesOfFreedom( void ) const
00634   {
00635     return this->GetNumberOfNodes() * this->GetNumberOfDegreesOfFreedomPerNode();
00636   }
00637 
00645   virtual unsigned int GetNumberOfDegreesOfFreedomPerNode( void ) const = 0;
00646 
00647 
00648 
00649 
00651   /*
00652    * Methods and classes related to IO and drawing
00653    */
00654 
00655 #ifdef FEM_BUILD_VISUALIZATION
00656 
00659   virtual void Draw(CDC* pDC, Solution::ConstPointer sol) const {}
00661   static double DC_Scale;
00662 #endif
00663 
00664 };
00665 
00666 // Make sure that Element::Node class is registered with the object factory.
00667 static INITClass Initializer_ElementNode(Element::Node::CLID());
00668 
00669 // Alias for Element::Node class
00670 typedef Element::Node Node;
00671 
00672 
00673 
00674 
00686 class ReadInfoType
00687 {
00688 public:
00689 
00690   typedef Node::ArrayType::ConstPointer NodeArrayPointer;
00691   typedef Element::ArrayType::ConstPointer ElementArrayPointer;
00692   typedef Material::ArrayType::ConstPointer MaterialArrayPointer;
00693 
00695   NodeArrayPointer m_node;
00696 
00698   ElementArrayPointer m_el;
00699 
00701   MaterialArrayPointer m_mat;
00702 
00704   ReadInfoType( NodeArrayPointer node_, ElementArrayPointer el_, MaterialArrayPointer mat_) :
00705     m_node(node_), m_el(el_), m_mat(mat_) {}
00706 };
00707 
00708 
00709 
00710 
00711 }} // end namespace itk::fem
00712 
00713 #endif // #ifndef __itkFEMElementBase_h

Generated at Wed Mar 12 01:12:54 2003 for ITK by doxygen 1.2.15 written by Dimitri van Heesch, © 1997-2000