Go to the documentation of this file.
16 #ifndef SURGSIM_PHYSICS_FEMELEMENT_H
17 #define SURGSIM_PHYSICS_FEMELEMENT_H
73 size_t getNodeId(
size_t elementNodeId)
const;
166 virtual void addMatVec(
double alphaM,
double alphaD,
double alphaK,
202 template <
typename DerivedSub,
typename T,
int Opt,
typename Index>
204 size_t blockSize, Eigen::SparseMatrix<T, Opt, Index>* matrix,
271 #endif // SURGSIM_PHYSICS_FEMELEMENT_H
double getMassDensity() const
Gets the mass density (in Kg.m-3)
Definition: FemElement.cpp:101
void setYoungModulus(double E)
Sets the Young modulus (in N.m-2)
Definition: FemElement.cpp:76
bool m_initializedFMDK
Flag to check in the f, M, D, K variables have been initialized.
Definition: FemElement.h:262
void setPoissonRatio(double nu)
Sets the Poisson ratio (unitless)
Definition: FemElement.cpp:86
virtual void doUpdateFMDK(const Math::OdeState &state, int options)=0
Update the FemElement based on the given state.
double m_E
Young modulus (in N.m-2)
Definition: FemElement.h:240
The state of an ode of 2nd order of the form with boundary conditions.
Definition: OdeState.h:38
void setNumDofPerNode(size_t numDofPerNode)
Sets the number of degrees of freedom per node.
Definition: FemElement.cpp:56
std::vector< size_t > m_nodeIds
Node ids connected by this element.
Definition: FemElement.h:234
FemElement()
Constructor.
Definition: FemElement.cpp:27
size_t getNumNodes() const
Gets the number of nodes connected by this element.
Definition: FemElement.cpp:61
SurgSim::Math::Vector m_f
The force vector.
Definition: FemElement.h:246
virtual ~FemElement()
Virtual destructor.
Definition: FemElement.cpp:31
bool isValidCoordinate(const SurgSim::Math::Vector &naturalCoordinate) const
Determines whether a given natural coordinate is valid.
Definition: FemElement.cpp:185
static FactoryType & getFactory()
Definition: FemElement.cpp:45
SurgSim::Math::Matrix m_M
The mass matrix.
Definition: FemElement.h:249
double m_rho
Mass density (in Kg.m-3)
Definition: FemElement.h:237
Definition: CompoundShapeToGraphics.cpp:29
void setMassDensity(double rho)
Sets the mass density (in Kg.m-3)
Definition: FemElement.cpp:96
double getPoissonRatio() const
Gets the Poisson ratio (unitless)
Definition: FemElement.cpp:91
An object factory, once a class is registered with the factory it can be used to create instances of ...
Definition: ObjectFactory.h:83
virtual void doInitializeFMDK()
Function to be overridden by the derived classes to initialize the f, M, D, K variables.
Definition: FemElement.cpp:208
void initializeFMDK()
Initialize f, M, D, K variables.
Definition: FemElement.cpp:199
size_t getNumDofPerNode() const
Gets the number of degree of freedom per node.
Definition: FemElement.cpp:51
virtual double getVolume(const SurgSim::Math::OdeState &state) const =0
Gets the element volume based on the input state (in m-3)
virtual void addForce(SurgSim::Math::Vector *F, double scale=1.0) const
Adds the element force (computed for a given state) to a complete system force vector F (assembly)
Definition: FemElement.cpp:111
Eigen::SparseMatrix< double > SparseMatrix
A sparse matrix.
Definition: SparseMatrix.h:32
virtual void initialize(const SurgSim::Math::OdeState &state)
Initialize the FemElement once everything has been set.
Definition: FemElement.cpp:34
virtual void addFMDK(SurgSim::Math::Vector *F, SurgSim::Math::SparseMatrix *M, SurgSim::Math::SparseMatrix *D, SurgSim::Math::SparseMatrix *K) const
Adds the element force vector, mass, stiffness and damping matrices (computed for a given state) into...
Definition: FemElement.cpp:134
virtual void addStiffness(SurgSim::Math::SparseMatrix *K, double scale=1.0) const
Adds the element stiffness matrix K (= -df/dx) (computed for a given state) to a complete system stif...
Definition: FemElement.cpp:129
virtual SurgSim::Math::Vector computeNaturalCoordinate(const SurgSim::Math::OdeState &state, const SurgSim::Math::Vector &cartesianCoordinate) const =0
Computes a natural coordinate given a global coordinate.
virtual SurgSim::Math::Vector computeCartesianCoordinate(const SurgSim::Math::OdeState &state, const SurgSim::Math::Vector &naturalCoordinate) const =0
Computes a given natural coordinate in cartesian coordinates.
SurgSim::Framework::ObjectFactory1< FemElement, std::shared_ptr< FemElementStructs::FemElementParameter > > FactoryType
Definition: FemElement.h:59
double m_nu
Poisson ratio (unitless)
Definition: FemElement.h:243
virtual void addDamping(SurgSim::Math::SparseMatrix *D, double scale=1.0) const
Adds the element damping matrix D (= -df/dv) (comuted for a given state) to a complete system damping...
Definition: FemElement.cpp:121
void updateFMDK(const Math::OdeState &state, int options)
Update the FemElement based on the given state.
Definition: FemElement.cpp:194
virtual void addMatVec(double alphaM, double alphaD, double alphaK, const SurgSim::Math::Vector &x, SurgSim::Math::Vector *F) const
Adds the element matrix-vector contribution F += (alphaM.M + alphaD.D + alphaK.K)....
Definition: FemElement.cpp:145
Eigen::Matrix< double, Eigen::Dynamic, 1 > Vector
A dynamic size column vector.
Definition: Vector.h:68
double getMass(const SurgSim::Math::OdeState &state) const
Gets the element mass based on the input state (in Kg)
Definition: FemElement.cpp:106
Base class for all Fem Element (1D, 2D, 3D) It handles the node ids to which it is connected and requ...
Definition: FemElement.h:45
const std::vector< size_t > & getNodeIds() const
Gets the node ids for this element.
Definition: FemElement.cpp:71
virtual void addMass(SurgSim::Math::SparseMatrix *M, double scale=1.0) const
Adds the element mass matrix M (computed for a given state) to a complete system mass matrix M (assem...
Definition: FemElement.cpp:116
size_t getNodeId(size_t elementNodeId) const
Gets the elementNodeId-th node id.
Definition: FemElement.cpp:66
bool m_useDamping
Flag to specify of the damping is used.
Definition: FemElement.h:255
SurgSim::Math::Matrix m_K
The stiffness matrix.
Definition: FemElement.h:258
SurgSim::Math::Matrix m_D
The damping matrix.
Definition: FemElement.h:252
void assembleMatrixBlocks(const DerivedSub &subMatrix, const std::vector< size_t > blockIds, size_t blockSize, Eigen::SparseMatrix< T, Opt, Index > *matrix, bool initialize=true) const
Helper method to add a sub-matrix made of squared-blocks into a matrix, for the sake of clarity.
Definition: FemElement-inl.h:31
size_t m_numDofPerNode
Number of degree of freedom per node for this element.
Definition: FemElement.h:231
Eigen::Matrix< double, Eigen::Dynamic, Eigen::Dynamic > Matrix
A dynamic size matrix.
Definition: Matrix.h:65
double getYoungModulus() const
Gets the Young modulus (in N.m-2)
Definition: FemElement.cpp:81