19 #ifndef itkFEMSolver_h 20 #define itkFEMSolver_h 72 template <
unsigned int VDimension = 3>
118 itkSetMacro(Origin, InterpolationGridPointType);
119 itkGetMacro(Origin, InterpolationGridPointType);
123 itkSetMacro(Spacing, InterpolationGridSpacingType);
124 itkGetMacro(Spacing, InterpolationGridSpacingType);
128 itkSetMacro(
Region, InterpolationGridRegionType);
129 itkGetMacro(
Region, InterpolationGridRegionType);
133 itkSetMacro(Direction, InterpolationGridDirectionType);
134 itkGetMacro(Direction, InterpolationGridDirectionType);
148 Float
GetSolution(
unsigned int i,
unsigned int which = 0);
153 virtual void SetInput( FEMObjectType *fem);
155 virtual void SetInput(
unsigned int, FEMObjectType * fem);
159 FEMObjectType *
GetInput(
unsigned int idx);
215 const InterpolationGridPointType & bb2);
220 InterpolationGridPointType bb1;
225 InterpolationGridPointType bb2;
228 bb2[i] = size[i] - 1.0;
236 const InterpolationGridPointType & origin,
237 const InterpolationGridSpacingType & spacing,
238 const InterpolationGridDirectionType & direction);
260 virtual DataObjectPointer
MakeOutput(DataObjectPointerArraySizeType) ITK_OVERRIDE;
277 FEMObjectType * GetOutput(
unsigned int idx);
351 void ApplyBC(
int dim = 0,
unsigned int matrix = 0);
409 Solver(
const Self &) ITK_DELETE_FUNCTION;
410 void operator=(const Self &) ITK_DELETE_FUNCTION;
422 #ifndef ITK_MANUAL_INSTANTIATION 423 #include "itkFEMSolver.hxx" 426 #endif // #ifndef itkFEMSolver_h void ApplyBC(int dim=0, unsigned int matrix=0)
Float GetSolution(unsigned int i, unsigned int which=0)
Superclass::RegionType RegionType
void InitializeInterpolationGrid(const InterpolationGridSizeType &size)
InterpolationGridType::PointType InterpolationGridPointType
Implements N-dimensional Finite element (FE) models including elements, materials, and loads.
void FillInterpolationGrid()
void SetLinearSystemWrapper(LinearSystemWrapper::Pointer ls)
virtual void FinalizeMatrixAfterAssembly(void)
Represent the size (bounds) of a n-dimensional image.
InterpolationGridType::SpacingType InterpolationGridSpacingType
void InitializeInterpolationGrid(const InterpolationGridSizeType &size, const InterpolationGridPointType &bb1, const InterpolationGridPointType &bb2)
InterpolationGridType::DirectionType InterpolationGridDirectionType
InterpolationGridPointerType m_InterpolationGrid
InterpolationGridSpacingType m_Spacing
virtual void InitializeLinearSystemWrapper()
Element::VectorType VectorType
FEMObjectType::ConstPointer FEMObjectConstPointer
FEMObjectType::Pointer FEMObjectPointer
virtual DataObjectPointer MakeOutput(DataObjectPointerArraySizeType) override
InterpolationGridDirectionType m_Direction
The base class for all process objects (source, filters, mappers) in the Insight data processing pipe...
The "itk" namespace contains all Insight Segmentation and Registration Toolkit (ITK) classes...
itk::Image< Element::ConstPointer, VDimension > InterpolationGridType
ObjectType * GetPointer() const
void Fill(const ValueType &)
itk::fem::FEMObject< VDimension > FEMObjectType
Load::ArrayType LoadArray
SmartPointer< Self > Pointer
InterpolationGridType::RegionType InterpolationGridRegionType
FEMObjectPointer m_FEMObject
static ITK_CONSTEXPR double e
The base of the natural logarithm or Euler's number
static const unsigned int MaxDimensions
virtual void AssembleElementMatrix(Element::Pointer e)
FEMObjectType * GetOutput()
Defines all functions required by Solver class to allocate, assemble and solve a linear system of equ...
InterpolationGridType::IndexType InterpolationGridIndexType
LinearSystemWrapper::Pointer m_LinearSystem
virtual void AssembleLandmarkContribution(Element::ConstPointer e, float)
LinearSystemWrapper::Pointer GetLinearSystemWrapper()
DataObjectPointerArray::size_type DataObjectPointerArraySizeType
virtual void SetTimeStep(Float dt)
Superclass::IndexType IndexType
vnl_vector< Float > VectorType
virtual void SetInput(const DataObjectIdentifierType &key, DataObject *input)
Protected method for setting indexed and named inputs.
InterpolationGridRegionType m_Region
Element::ArrayType ElementArray
InterpolationGridType::Pointer InterpolationGridPointerType
InterpolationGridPointType m_Origin
FEM solver used to generate a solution for a FE formulation.
A region represents some portion or piece of data.
virtual void InitializeMatrixForAssembly(unsigned int N)
virtual DataObjectPointer MakeOutput(DataObjectPointerArraySizeType idx)
Make a DataObject of the correct type to used as the specified output.
Element::Node::ArrayType NodeArray
virtual void SetInput(FEMObjectType *fem)
virtual void PrintSelf(std::ostream &os, Indent indent) const override
void AssembleF(int dim=0)
Abstract base element class.
InterpolationGridType::SizeType InterpolationGridSizeType
static const unsigned int FEMDimension
void UpdateDisplacements()
ProcessObject::DataObjectPointerArraySizeType DataObjectPointerArraySizeType
DataObject::Pointer DataObjectPointer
virtual Float GetTimeStep() const
const InterpolationGridType * GetInterpolationGrid(void) const
Control indentation during Print() invocation.
Float GetDeformationEnergy(unsigned int SolutionIndex=0)
FEMObjectType * GetInput()
virtual void GenerateData() override
const Element * GetElementAtPoint(const VectorType &pt) const
SmartPointer< const Self > ConstPointer
LinearSystemWrapper class that uses VNL numeric library functions to define a sparse linear system of...
Templated n-dimensional image class.
LinearSystemWrapperVNL m_LinearSystemVNL
Material::ArrayType MaterialArray