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 
21 #include "SurgSim/Math/Matrix.h"
22 #include "SurgSim/Math/Vector.h"
23 
24 namespace SurgSim
25 {
26 
27 namespace Math
28 {
29 class OdeState;
30 };
31 
32 namespace Physics
33 {
34 
43 {
44 public:
46  FemElement();
47 
49  virtual ~FemElement();
50 
53  virtual void initialize(const SurgSim::Math::OdeState& state);
54 
57  size_t getNumDofPerNode() const;
58 
61  size_t getNumNodes() const;
62 
65  size_t getNodeId(size_t elementNodeId) const;
66 
69  const std::vector<size_t>& getNodeIds() const;
70 
73  void setYoungModulus(double E);
76  double getYoungModulus() const;
77 
80  void setPoissonRatio(double nu);
83  double getPoissonRatio() const;
84 
87  void setMassDensity(double rho);
90  double getMassDensity() const;
91 
95  double getMass(const SurgSim::Math::OdeState& state) const;
96 
100  virtual double getVolume(const SurgSim::Math::OdeState& state) const = 0;
101 
109  virtual void addForce(const SurgSim::Math::OdeState& state, SurgSim::Math::Vector* F, double scale = 1.0) = 0;
110 
118  virtual void addMass(const SurgSim::Math::OdeState& state, SurgSim::Math::Matrix* M, double scale = 1.0) = 0;
119 
128  virtual void addDamping(const SurgSim::Math::OdeState& state, SurgSim::Math::Matrix* D,
129  double scale = 1.0) = 0;
130 
139  virtual void addStiffness(const SurgSim::Math::OdeState& state, SurgSim::Math::Matrix* K,
140  double scale = 1.0) = 0;
141 
151  virtual void addFMDK(const SurgSim::Math::OdeState& state,
155  SurgSim::Math::Matrix* K) = 0;
156 
167  virtual void addMatVec(const SurgSim::Math::OdeState& state, double alphaM, double alphaD, double alphaK,
169 
174  virtual bool update(const SurgSim::Math::OdeState& state);
175 
179  bool isValidCoordinate(const SurgSim::Math::Vector& naturalCoordinate) const;
180 
185  virtual SurgSim::Math::Vector computeCartesianCoordinate(
186  const SurgSim::Math::OdeState& state,
187  const SurgSim::Math::Vector& naturalCoordinate) const = 0;
188 
193  virtual SurgSim::Math::Vector computeNaturalCoordinate(
194  const SurgSim::Math::OdeState& state,
195  const SurgSim::Math::Vector& cartesianCoordinate) const = 0;
196 
197 protected:
202  void setNumDofPerNode(size_t numDofPerNode);
203 
206 
208  std::vector<size_t> m_nodeIds;
209 
211  double m_rho;
212 
214  double m_E;
215 
217  double m_nu;
218 };
219 
220 } // namespace Physics
221 
222 } // namespace SurgSim
223 
224 #endif // SURGSIM_PHYSICS_FEMELEMENT_H
Definition: DriveElementFromInputBehavior.cpp:27
size_t m_numDofPerNode
Number of degree of freedom per node for this element.
Definition: FemElement.h:205
Eigen::Matrix< double, Eigen::Dynamic, Eigen::Dynamic > Matrix
A dynamic size matrix.
Definition: Matrix.h:65
double m_nu
Poisson ratio (unitless)
Definition: FemElement.h:217
OdeState defines the state y of an ode of 2nd order of the form M(x,v).a = F(x, v) with boundary cond...
Definition: OdeState.h:34
Base class for all Fem Element (1D, 2D, 3D) It handles the node ids to which it is connected and requ...
Definition: FemElement.h:42
double m_rho
Mass density (in Kg.m-3)
Definition: FemElement.h:211
Eigen::Matrix< double, Eigen::Dynamic, 1 > Vector
A dynamic size column vector.
Definition: Vector.h:67
double m_E
Young modulus (in N.m-2)
Definition: FemElement.h:214
Definitions of small fixed-size square matrix types.
Definitions of small fixed-size vector types.
std::vector< size_t > m_nodeIds
Node ids connected by this element.
Definition: FemElement.h:208