Fem3DElementTetrahedron.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_FEM3DELEMENTTETRAHEDRON_H
17 #define SURGSIM_PHYSICS_FEM3DELEMENTTETRAHEDRON_H
18 
19 #include <array>
20 
21 #include "SurgSim/Math/Matrix.h"
22 #include "SurgSim/Math/Vector.h"
24 
25 namespace SurgSim
26 {
27 
28 namespace Physics
29 {
30 
39 {
40 public:
47  explicit Fem3DElementTetrahedron(std::array<size_t, 4> nodeIds);
48 
56  virtual void initialize(const SurgSim::Math::OdeState& state) override;
57 
60  virtual double getVolume(const SurgSim::Math::OdeState& state) const override;
61 
69  virtual void addForce(const SurgSim::Math::OdeState& state, SurgSim::Math::Vector* F,
70  double scale = 1.0) override;
71 
79  virtual void addMass(const SurgSim::Math::OdeState& state, SurgSim::Math::Matrix* M,
80  double scale = 1.0) override;
81 
91  virtual void addDamping(const SurgSim::Math::OdeState& state, SurgSim::Math::Matrix* D,
92  double scale = 1.0) override;
93 
102  virtual void addStiffness(const SurgSim::Math::OdeState& state, SurgSim::Math::Matrix* K,
103  double scale = 1.0) override;
104 
114  virtual void addFMDK(const SurgSim::Math::OdeState& state,
118  SurgSim::Math::Matrix* K) override;
119 
130  virtual void addMatVec(const SurgSim::Math::OdeState& state,
131  double alphaM, double alphaD, double alphaK,
132  const SurgSim::Math::Vector& x, SurgSim::Math::Vector* F) override;
133 
135  const SurgSim::Math::OdeState& state,
136  const SurgSim::Math::Vector& naturalCoordinate) const override;
137 
139  const SurgSim::Math::OdeState& state,
140  const SurgSim::Math::Vector& cartesianCoordinate) const override;
141 
142 protected:
151  double* volume,
152  std::array<double, 4>* ai,
153  std::array<double, 4>* bi,
154  std::array<double, 4>* ci,
155  std::array<double, 4>* di) const;
156 
160  void computeStiffness(const SurgSim::Math::OdeState& state,
161  Eigen::Matrix<double, 12, 12>* k);
162 
166  void computeMass(const SurgSim::Math::OdeState& state,
167  Eigen::Matrix<double, 12, 12>* m);
168 
178  void addForce(const SurgSim::Math::OdeState& state, const Eigen::Matrix<double, 12, 12>& k,
179  SurgSim::Math::Vector* F, double scale = 1.0);
180 
182  double m_restVolume;
184  std::array<double, 4> m_ai, m_bi, m_ci, m_di;
185 
187  Eigen::Matrix<double, 12, 1> m_x0;
188 
190  Eigen::Matrix<double, 6, 6> m_Em;
192  Eigen::Matrix<double, 6, 12> m_strain;
194  Eigen::Matrix<double, 6, 12> m_stress;
195 
197  Eigen::Matrix<double, 12, 12> m_M;
199  Eigen::Matrix<double, 12, 12> m_K;
200 };
201 
202 } // namespace Physics
203 
204 } // namespace SurgSim
205 
206 #endif // SURGSIM_PHYSICS_FEM3DELEMENTTETRAHEDRON_H
Definition: DriveElementFromInputBehavior.cpp:27
Eigen::Matrix< double, 6, 12 > m_strain
Strain matrix.
Definition: Fem3DElementTetrahedron.h:192
std::array< double, 4 > m_bi
Definition: Fem3DElementTetrahedron.h:184
Class for Fem Element 3D based on a tetrahedron volume discretization.
Definition: Fem3DElementTetrahedron.h:38
std::array< double, 4 > m_ai
Shape functions coefficients Ni(x,y,z) = 1/6V ( ai + x.bi + y.ci + z.di )
Definition: Fem3DElementTetrahedron.h:184
Eigen::Matrix< double, Eigen::Dynamic, Eigen::Dynamic > Matrix
A dynamic size matrix.
Definition: Matrix.h:65
void computeShapeFunctions(const SurgSim::Math::OdeState &state, double *volume, std::array< double, 4 > *ai, std::array< double, 4 > *bi, std::array< double, 4 > *ci, std::array< double, 4 > *di) const
Computes the tetrahedron shape functions.
Definition: Fem3DElementTetrahedron.cpp:277
virtual void addFMDK(const SurgSim::Math::OdeState &state, SurgSim::Math::Vector *F, SurgSim::Math::Matrix *M, SurgSim::Math::Matrix *D, SurgSim::Math::Matrix *K) override
Adds the element force vector, mass, stiffness and damping matrices (computed for a given state) into...
Definition: Fem3DElementTetrahedron.cpp:208
Eigen::Matrix< double, 12, 1 > m_x0
The tetrahedon rest state.
Definition: Fem3DElementTetrahedron.h:187
virtual void addMass(const SurgSim::Math::OdeState &state, SurgSim::Math::Matrix *M, double scale=1.0) override
Adds the element mass matrix M (computed for a given state) to a complete system mass matrix M (assem...
Definition: Fem3DElementTetrahedron.cpp:154
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
Eigen::Matrix< double, 6, 12 > m_stress
Stress matrix.
Definition: Fem3DElementTetrahedron.h:194
virtual SurgSim::Math::Vector computeCartesianCoordinate(const SurgSim::Math::OdeState &state, const SurgSim::Math::Vector &naturalCoordinate) const override
Computes a given natural coordinate in cartesian coordinates.
Definition: Fem3DElementTetrahedron.cpp:401
virtual void addDamping(const SurgSim::Math::OdeState &state, SurgSim::Math::Matrix *D, double scale=1.0) override
Adds the element damping matrix D (= -df/dv) (comuted for a given state) to a complete system damping...
Definition: Fem3DElementTetrahedron.cpp:159
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
Eigen::Matrix< double, 6, 6 > m_Em
Elasticity material matrix (contains the elastic properties of the material)
Definition: Fem3DElementTetrahedron.h:190
Eigen::Matrix< double, Eigen::Dynamic, 1 > Vector
A dynamic size column vector.
Definition: Vector.h:67
Eigen::Matrix< double, 12, 12 > m_M
Mass matrix.
Definition: Fem3DElementTetrahedron.h:197
std::array< double, 4 > m_ci
Definition: Fem3DElementTetrahedron.h:184
Eigen::Matrix< double, 12, 12 > m_K
Stiffness matrix.
Definition: Fem3DElementTetrahedron.h:199
virtual void addForce(const SurgSim::Math::OdeState &state, SurgSim::Math::Vector *F, double scale=1.0) override
Adds the element force (computed for a given state) to a complete system force vector F (assembly) ...
Definition: Fem3DElementTetrahedron.cpp:102
Definitions of small fixed-size square matrix types.
std::array< double, 4 > m_di
Definition: Fem3DElementTetrahedron.h:184
Definitions of small fixed-size vector types.
double m_restVolume
Shape functions: Tetrahedron rest volume.
Definition: Fem3DElementTetrahedron.h:182
virtual void addStiffness(const SurgSim::Math::OdeState &state, SurgSim::Math::Matrix *K, double scale=1.0) override
Adds the element stiffness matrix K (= -df/dx) (computed for a given state) to a complete system stif...
Definition: Fem3DElementTetrahedron.cpp:202
void computeStiffness(const SurgSim::Math::OdeState &state, Eigen::Matrix< double, 12, 12 > *k)
Computes the tetrahedron stiffness matrix.
Definition: Fem3DElementTetrahedron.cpp:163
virtual SurgSim::Math::Vector computeNaturalCoordinate(const SurgSim::Math::OdeState &state, const SurgSim::Math::Vector &cartesianCoordinate) const override
Computes a natural coordinate given a global coordinate.
Definition: Fem3DElementTetrahedron.cpp:420
virtual void initialize(const SurgSim::Math::OdeState &state) override
Initialize the FemElement once everything has been set.
Definition: Fem3DElementTetrahedron.cpp:55
virtual double getVolume(const SurgSim::Math::OdeState &state) const override
Get the element volume based on the input state.
Definition: Fem3DElementTetrahedron.cpp:258
Fem3DElementTetrahedron(std::array< size_t, 4 > nodeIds)
Constructor.
Definition: Fem3DElementTetrahedron.cpp:48
virtual void addMatVec(const SurgSim::Math::OdeState &state, double alphaM, double alphaD, double alphaK, const SurgSim::Math::Vector &x, SurgSim::Math::Vector *F) override
Adds the element matrix-vector contribution F += (alphaM.M + alphaD.D + alphaK.K).x (computed for a given state) into a complete system data structure F (assembly)
Definition: Fem3DElementTetrahedron.cpp:226
void computeMass(const SurgSim::Math::OdeState &state, Eigen::Matrix< double, 12, 12 > *m)
Computes the tetrahedron mass matrix.
Definition: Fem3DElementTetrahedron.cpp:108