All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Fem3DElementCube.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_FEM3DELEMENTCUBE_H
17 #define SURGSIM_PHYSICS_FEM3DELEMENTCUBE_H
18 
19 #include <array>
20 
23 
24 namespace SurgSim
25 {
26 
27 namespace Physics
28 {
29 
41 {
42 public:
51  explicit Fem3DElementCube(std::array<size_t, 8> nodeIds);
52 
63  virtual void initialize(const SurgSim::Math::OdeState& state) override;
64 
67  virtual double getVolume(const SurgSim::Math::OdeState& state) const override;
68 
76  virtual void addForce(const SurgSim::Math::OdeState& state, SurgSim::Math::Vector* F,
77  double scale = 1.0) override;
78 
86  virtual void addMass(const SurgSim::Math::OdeState& state, SurgSim::Math::Matrix* M,
87  double scale = 1.0) override;
88 
98  virtual void addDamping(const SurgSim::Math::OdeState& state, SurgSim::Math::Matrix* D,
99  double scale = 1.0) override;
100 
109  virtual void addStiffness(const SurgSim::Math::OdeState& state, SurgSim::Math::Matrix* K,
110  double scale = 1.0) override;
111 
121  virtual void addFMDK(const SurgSim::Math::OdeState& state,
125  SurgSim::Math::Matrix* K) override;
126 
137  virtual void addMatVec(const SurgSim::Math::OdeState& state,
138  double alphaM, double alphaD, double alphaK,
139  const SurgSim::Math::Vector& x, SurgSim::Math::Vector* F) override;
140 
142  const SurgSim::Math::OdeState& state,
143  const SurgSim::Math::Vector& naturalCoordinate) const override;
144 
146  const SurgSim::Math::OdeState& state,
147  const SurgSim::Math::Vector& cartesianCoordinate) const override;
148 
149 protected:
152  void buildConstitutiveMaterialMatrix(Eigen::Matrix<double, 6, 6>* constitutiveMatrix);
153 
157  void computeStiffness(const SurgSim::Math::OdeState& state,
158  Eigen::Matrix<double, 6, 24>* strain,
159  Eigen::Matrix<double, 6, 24>* stress,
160  Eigen::Matrix<double, 24, 24>* k);
161 
165  void computeMass(const SurgSim::Math::OdeState& state, Eigen::Matrix<double, 24, 24>* m);
166 
176  void addForce(const SurgSim::Math::OdeState& state, const Eigen::Matrix<double, 24, 24>& k,
177  SurgSim::Math::Vector* F, double scale = 1.0);
178 
188  Eigen::Matrix<double, 6, 24>* strain,
189  Eigen::Matrix<double, 6, 24>* stress,
190  Eigen::Matrix<double, 24, 24>* k);
191 
200  Eigen::Matrix<double, 24, 24>* m);
201 
207  void evaluateJ(const SurgSim::Math::OdeState& state, double epsilon, double eta, double mu,
210  double *detJ) const;
211 
217  void evaluateStrainDisplacement(double epsilon, double eta, double mu, const SurgSim::Math::Matrix33d& Jinv,
218  Eigen::Matrix<double, 6, 24> *B) const;
219 
221  double m_restVolume;
222 
231  std::array<double, 8> m_shapeFunctionsEpsilonSign;
232  std::array<double, 8> m_shapeFunctionsEtaSign;
233  std::array<double, 8> m_shapeFunctionsMuSign;
235 
238 
253  double shapeFunction(size_t i, double epsilon, double eta, double mu) const;
260 
268  double dShapeFunctiondepsilon(size_t i, double epsilon, double eta, double mu) const;
269 
277  double dShapeFunctiondeta(size_t i, double epsilon, double eta, double mu) const;
278 
286  double dShapeFunctiondmu(size_t i, double epsilon, double eta, double mu) const;
287 
289  Eigen::Matrix<double, 24, 1> m_elementRestPosition;
290 
292  Eigen::Matrix<double, 6, 24> m_strain;
294  Eigen::Matrix<double, 6, 24> m_stress;
296  Eigen::Matrix<double, 6, 6> m_constitutiveMaterial;
297 
299  Eigen::Matrix<double, 24, 24> m_mass;
301  Eigen::Matrix<double, 24, 24> m_stiffness;
302 };
303 
304 } // namespace Physics
305 
306 } // namespace SurgSim
307 
308 #endif // SURGSIM_PHYSICS_FEM3DELEMENTCUBE_H
void computeStiffness(const SurgSim::Math::OdeState &state, Eigen::Matrix< double, 6, 24 > *strain, Eigen::Matrix< double, 6, 24 > *stress, Eigen::Matrix< double, 24, 24 > *k)
Computes the cube stiffness matrix along with the strain and stress matrices.
Definition: Fem3DElementCube.cpp:131
void evaluateStrainDisplacement(double epsilon, double eta, double mu, const SurgSim::Math::Matrix33d &Jinv, Eigen::Matrix< double, 6, 24 > *B) const
Helper method to evaluate the strain-displacement matrix at a given 3D parametric location c...
Definition: Fem3DElementCube.cpp:200
void addMassMatrixAtPoint(const SurgSim::Math::OdeState &state, const SurgSim::Math::gaussQuadraturePoint &epsilon, const SurgSim::Math::gaussQuadraturePoint &eta, const SurgSim::Math::gaussQuadraturePoint &mu, Eigen::Matrix< double, 24, 24 > *m)
Helper method to evaluate mass integral terms with a discrete sum using a Gauss quadrature rule...
Definition: Fem3DElementCube.cpp:277
double dShapeFunctiondmu(size_t i, double epsilon, double eta, double mu) const
Shape functions derivative .
Definition: Fem3DElementCube.cpp:432
Eigen::Matrix< double, 24, 24 > m_mass
Mass matrix (usually noted )
Definition: Fem3DElementCube.h:299
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: Fem3DElementCube.cpp:86
double shapeFunction(size_t i, double epsilon, double eta, double mu) const
Shape functions .
Definition: Fem3DElementCube.cpp:408
double dShapeFunctiondeta(size_t i, double epsilon, double eta, double mu) const
Shape functions derivative .
Definition: Fem3DElementCube.cpp:424
Eigen::Matrix< double, Eigen::Dynamic, Eigen::Dynamic > Matrix
A dynamic size matrix.
Definition: Matrix.h:65
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: Fem3DElementCube.cpp:459
void addStrainStressStiffnessAtPoint(const SurgSim::Math::OdeState &state, const SurgSim::Math::gaussQuadraturePoint &epsilon, const SurgSim::Math::gaussQuadraturePoint &eta, const SurgSim::Math::gaussQuadraturePoint &mu, Eigen::Matrix< double, 6, 24 > *strain, Eigen::Matrix< double, 6, 24 > *stress, Eigen::Matrix< double, 24, 24 > *k)
Helper method to evaluate strain-stress and stiffness integral terms with a discrete sum using a Gaus...
Definition: Fem3DElementCube.cpp:254
Eigen::Matrix< double, 24, 1 > m_elementRestPosition
The cube rest state (nodes ordered by m_nodeIds)
Definition: Fem3DElementCube.h:289
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
void evaluateJ(const SurgSim::Math::OdeState &state, double epsilon, double eta, double mu, SurgSim::Math::Matrix33d *J, SurgSim::Math::Matrix33d *Jinv, double *detJ) const
Helper method to evaluate matrix J = d(x,y,z)/d(epsilon,eta,mu) at a given 3D parametric location J e...
Definition: Fem3DElementCube.cpp:156
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
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: Fem3DElementCube.cpp:122
std::array< double, 8 > m_shapeFunctionsMuSign
Definition: Fem3DElementCube.h:233
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: Fem3DElementCube.cpp:313
Definitions of a n-point Gaussian quadrature rule (a.k.a.
Eigen::Matrix< double, 24, 24 > m_stiffness
Stiffness matrix (usually noted )
Definition: Fem3DElementCube.h:301
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: Fem3DElementCube.cpp:308
Eigen::Matrix< double, Eigen::Dynamic, 1 > Vector
A dynamic size column vector.
Definition: Vector.h:67
Class for Fem Element 3D based on a cube volume discretization.
Definition: Fem3DElementCube.h:40
std::array< double, 8 > m_shapeFunctionsEpsilonSign
Definition: Fem3DElementCube.h:231
Eigen::Matrix< double, 6, 6 > m_constitutiveMaterial
Constitutive material matrix (Hooke's law in this case) defines the relationship between stress and s...
Definition: Fem3DElementCube.h:296
std::array< double, 8 > m_shapeFunctionsEtaSign
Definition: Fem3DElementCube.h:232
void buildConstitutiveMaterialMatrix(Eigen::Matrix< double, 6, 6 > *constitutiveMatrix)
Build the constitutive material 6x6 matrix.
Definition: Fem3DElementCube.cpp:240
virtual void initialize(const SurgSim::Math::OdeState &state) override
Initializes the element once everything has been set.
Definition: Fem3DElementCube.cpp:50
void computeMass(const SurgSim::Math::OdeState &state, Eigen::Matrix< double, 24, 24 > *m)
Computes the cube mass matrix.
Definition: Fem3DElementCube.cpp:91
double dShapeFunctiondepsilon(size_t i, double epsilon, double eta, double mu) const
Shape functions derivative .
Definition: Fem3DElementCube.cpp:416
virtual double getVolume(const SurgSim::Math::OdeState &state) const override
Gets the element volume based on the input state.
Definition: Fem3DElementCube.cpp:363
Eigen::Matrix< double, 3, 3, Eigen::RowMajor > Matrix33d
A 3x3 matrix of doubles.
Definition: Matrix.h:51
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: Fem3DElementCube.cpp:440
double m_restVolume
Cube rest volume.
Definition: Fem3DElementCube.h:221
Eigen::Matrix< double, 6, 24 > m_stress
Stress matrix (usually noted )
Definition: Fem3DElementCube.h:294
Definition: GaussLegendreQuadrature.h:31
Eigen::Matrix< double, 6, 24 > m_strain
Strain matrix (usually noted )
Definition: Fem3DElementCube.h:292
Fem3DElementCube(std::array< size_t, 8 > nodeIds)
Constructor.
Definition: Fem3DElementCube.cpp:33
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: Fem3DElementCube.cpp:331
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) (computed for a given state) to a complete system dampin...
Definition: Fem3DElementCube.cpp:127