FemElement.h
Go to the documentation of this file.
1 // This file is a part of the OpenSurgSim project.
2 // Copyright 2013, SimQuest Solutions Inc.
3 //
4 // Licensed under the Apache License, Version 2.0 (the "License");
5 // you may not use this file except in compliance with the License.
6 // You may obtain a copy of the License at
7 //
8 // http://www.apache.org/licenses/LICENSE-2.0
9 //
10 // Unless required by applicable law or agreed to in writing, software
11 // distributed under the License is distributed on an "AS IS" BASIS,
12 // WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
13 // See the License for the specific language governing permissions and
14 // limitations under the License.
15 
16 #ifndef SURGSIM_PHYSICS_FEMELEMENT_H
17 #define SURGSIM_PHYSICS_FEMELEMENT_H
18 
19 #include <vector>
20 
22 #include "SurgSim/Math/Matrix.h"
24 #include "SurgSim/Math/Vector.h"
25 #include "SurgSim/Physics/Fem.h"
26 
27 namespace SurgSim
28 {
29 
30 namespace Math
31 {
32 class OdeState;
33 };
34 
35 namespace Physics
36 {
37 
46 {
47 public:
49  FemElement();
50 
52  virtual ~FemElement();
53 
56  virtual void initialize(const SurgSim::Math::OdeState& state);
57 
60 
61  static FactoryType& getFactory();
62 
65  size_t getNumDofPerNode() const;
66 
69  size_t getNumNodes() const;
70 
73  size_t getNodeId(size_t elementNodeId) const;
74 
77  const std::vector<size_t>& getNodeIds() const;
78 
81  void setYoungModulus(double E);
84  double getYoungModulus() const;
85 
88  void setPoissonRatio(double nu);
91  double getPoissonRatio() const;
92 
95  void setMassDensity(double rho);
98  double getMassDensity() const;
99 
103  double getMass(const SurgSim::Math::OdeState& state) const;
104 
108  virtual double getVolume(const SurgSim::Math::OdeState& state) const = 0;
109 
116  virtual void addForce(SurgSim::Math::Vector* F, double scale = 1.0) const;
117 
124  virtual void addMass(SurgSim::Math::SparseMatrix* M, double scale = 1.0) const;
125 
133  virtual void addDamping(SurgSim::Math::SparseMatrix* D, double scale = 1.0) const;
134 
142  virtual void addStiffness(SurgSim::Math::SparseMatrix* K, double scale = 1.0) const;
143 
152  virtual void addFMDK(SurgSim::Math::Vector* F,
155  SurgSim::Math::SparseMatrix* K) const;
156 
166  virtual void addMatVec(double alphaM, double alphaD, double alphaK,
167  const SurgSim::Math::Vector& x, SurgSim::Math::Vector* F) const;
168 
172  bool isValidCoordinate(const SurgSim::Math::Vector& naturalCoordinate) const;
173 
179  const SurgSim::Math::OdeState& state,
180  const SurgSim::Math::Vector& naturalCoordinate) const = 0;
181 
187  const SurgSim::Math::OdeState& state,
188  const SurgSim::Math::Vector& cartesianCoordinate) const = 0;
189 
202  template <typename DerivedSub, typename T, int Opt, typename Index>
203  void assembleMatrixBlocks(const DerivedSub& subMatrix, const std::vector<size_t> blockIds,
204  size_t blockSize, Eigen::SparseMatrix<T, Opt, Index>* matrix,
205  bool initialize = true) const;
206 
210  void updateFMDK(const Math::OdeState& state, int options);
211 
212 protected:
217  void setNumDofPerNode(size_t numDofPerNode);
218 
222  virtual void doUpdateFMDK(const Math::OdeState& state, int options) = 0;
223 
225  void initializeFMDK();
226 
228  virtual void doInitializeFMDK();
229 
232 
234  std::vector<size_t> m_nodeIds;
235 
237  double m_rho;
238 
240  double m_E;
241 
243  double m_nu;
244 
247 
250 
253 
256 
259 
260 private:
263 };
264 
265 } // namespace Physics
266 
267 } // namespace SurgSim
268 
270 
271 #endif // SURGSIM_PHYSICS_FEMELEMENT_H
SurgSim::Physics::FemElement::getMassDensity
double getMassDensity() const
Gets the mass density (in Kg.m-3)
Definition: FemElement.cpp:101
SurgSim::Physics::FemElement::setYoungModulus
void setYoungModulus(double E)
Sets the Young modulus (in N.m-2)
Definition: FemElement.cpp:76
SurgSim::Physics::FemElement::m_initializedFMDK
bool m_initializedFMDK
Flag to check in the f, M, D, K variables have been initialized.
Definition: FemElement.h:262
SurgSim::Physics::FemElement::setPoissonRatio
void setPoissonRatio(double nu)
Sets the Poisson ratio (unitless)
Definition: FemElement.cpp:86
SurgSim::Physics::FemElement::doUpdateFMDK
virtual void doUpdateFMDK(const Math::OdeState &state, int options)=0
Update the FemElement based on the given state.
SurgSim::Physics::FemElement::m_E
double m_E
Young modulus (in N.m-2)
Definition: FemElement.h:240
Vector.h
SurgSim::Math::OdeState
The state of an ode of 2nd order of the form with boundary conditions.
Definition: OdeState.h:38
SurgSim::Physics::FemElement::setNumDofPerNode
void setNumDofPerNode(size_t numDofPerNode)
Sets the number of degrees of freedom per node.
Definition: FemElement.cpp:56
SurgSim::Physics::FemElement::m_nodeIds
std::vector< size_t > m_nodeIds
Node ids connected by this element.
Definition: FemElement.h:234
SurgSim::Physics::FemElement::FemElement
FemElement()
Constructor.
Definition: FemElement.cpp:27
SurgSim::Physics::FemElement::getNumNodes
size_t getNumNodes() const
Gets the number of nodes connected by this element.
Definition: FemElement.cpp:61
Matrix.h
SurgSim::Physics::FemElement::m_f
SurgSim::Math::Vector m_f
The force vector.
Definition: FemElement.h:246
SurgSim::Physics::FemElement::~FemElement
virtual ~FemElement()
Virtual destructor.
Definition: FemElement.cpp:31
SurgSim::Physics::FemElement::isValidCoordinate
bool isValidCoordinate(const SurgSim::Math::Vector &naturalCoordinate) const
Determines whether a given natural coordinate is valid.
Definition: FemElement.cpp:185
SurgSim::Physics::FemElement::getFactory
static FactoryType & getFactory()
Definition: FemElement.cpp:45
SurgSim::Physics::FemElement::m_M
SurgSim::Math::Matrix m_M
The mass matrix.
Definition: FemElement.h:249
SurgSim::Physics::FemElement::m_rho
double m_rho
Mass density (in Kg.m-3)
Definition: FemElement.h:237
SparseMatrix.h
ObjectFactory.h
SurgSim
Definition: CompoundShapeToGraphics.cpp:29
SurgSim::Physics::FemElement::setMassDensity
void setMassDensity(double rho)
Sets the mass density (in Kg.m-3)
Definition: FemElement.cpp:96
SurgSim::Physics::FemElement::getPoissonRatio
double getPoissonRatio() const
Gets the Poisson ratio (unitless)
Definition: FemElement.cpp:91
SurgSim::Framework::ObjectFactory1
An object factory, once a class is registered with the factory it can be used to create instances of ...
Definition: ObjectFactory.h:83
SurgSim::Physics::FemElement::doInitializeFMDK
virtual void doInitializeFMDK()
Function to be overridden by the derived classes to initialize the f, M, D, K variables.
Definition: FemElement.cpp:208
SurgSim::Physics::FemElement::initializeFMDK
void initializeFMDK()
Initialize f, M, D, K variables.
Definition: FemElement.cpp:199
SurgSim::Physics::FemElement::getNumDofPerNode
size_t getNumDofPerNode() const
Gets the number of degree of freedom per node.
Definition: FemElement.cpp:51
SurgSim::Physics::FemElement::getVolume
virtual double getVolume(const SurgSim::Math::OdeState &state) const =0
Gets the element volume based on the input state (in m-3)
SurgSim::Physics::FemElement::addForce
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
SurgSim::Math::SparseMatrix
Eigen::SparseMatrix< double > SparseMatrix
A sparse matrix.
Definition: SparseMatrix.h:32
SurgSim::Physics::FemElement::initialize
virtual void initialize(const SurgSim::Math::OdeState &state)
Initialize the FemElement once everything has been set.
Definition: FemElement.cpp:34
SurgSim::Physics::FemElement::addFMDK
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
SurgSim::Physics::FemElement::addStiffness
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
SurgSim::Physics::FemElement::computeNaturalCoordinate
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.
SurgSim::Physics::FemElement::computeCartesianCoordinate
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::Physics::FemElement::FactoryType
SurgSim::Framework::ObjectFactory1< FemElement, std::shared_ptr< FemElementStructs::FemElementParameter > > FactoryType
Definition: FemElement.h:59
SurgSim::Physics::FemElement::m_nu
double m_nu
Poisson ratio (unitless)
Definition: FemElement.h:243
SurgSim::Physics::FemElement::addDamping
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
SurgSim::Physics::FemElement::updateFMDK
void updateFMDK(const Math::OdeState &state, int options)
Update the FemElement based on the given state.
Definition: FemElement.cpp:194
SurgSim::Physics::FemElement::addMatVec
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
SurgSim::Math::Vector
Eigen::Matrix< double, Eigen::Dynamic, 1 > Vector
A dynamic size column vector.
Definition: Vector.h:68
SurgSim::Physics::FemElement::getMass
double getMass(const SurgSim::Math::OdeState &state) const
Gets the element mass based on the input state (in Kg)
Definition: FemElement.cpp:106
SurgSim::Physics::FemElement
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
FemElement-inl.h
SurgSim::Physics::FemElement::getNodeIds
const std::vector< size_t > & getNodeIds() const
Gets the node ids for this element.
Definition: FemElement.cpp:71
SurgSim::Physics::FemElement::addMass
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
SurgSim::Physics::FemElement::getNodeId
size_t getNodeId(size_t elementNodeId) const
Gets the elementNodeId-th node id.
Definition: FemElement.cpp:66
Fem.h
SurgSim::Physics::FemElement::m_useDamping
bool m_useDamping
Flag to specify of the damping is used.
Definition: FemElement.h:255
SurgSim::Physics::FemElement::m_K
SurgSim::Math::Matrix m_K
The stiffness matrix.
Definition: FemElement.h:258
SurgSim::Physics::FemElement::m_D
SurgSim::Math::Matrix m_D
The damping matrix.
Definition: FemElement.h:252
SurgSim::Physics::FemElement::assembleMatrixBlocks
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
SurgSim::Physics::FemElement::m_numDofPerNode
size_t m_numDofPerNode
Number of degree of freedom per node for this element.
Definition: FemElement.h:231
SurgSim::Math::Matrix
Eigen::Matrix< double, Eigen::Dynamic, Eigen::Dynamic > Matrix
A dynamic size matrix.
Definition: Matrix.h:65
SurgSim::Physics::FemElement::getYoungModulus
double getYoungModulus() const
Gets the Young modulus (in N.m-2)
Definition: FemElement.cpp:81