Fem1DElementBeam.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_FEM1DELEMENTBEAM_H
17 #define SURGSIM_PHYSICS_FEM1DELEMENTBEAM_H
18 
19 #include <array>
20 
22 
23 namespace SurgSim
24 {
25 
26 namespace Physics
27 {
28 
36 {
37 public:
41  explicit Fem1DElementBeam(std::array<size_t, 2> nodeIds);
42 
45  void setRadius(double radius);
46 
49  double getRadius() const;
50 
54  virtual void initialize(const SurgSim::Math::OdeState& state) override;
55 
59  virtual double getVolume(const SurgSim::Math::OdeState& state) const override;
60 
63  bool getShearingEnabled() const;
64 
69  void setShearingEnabled(bool enabled);
70 
78  virtual void addForce(const SurgSim::Math::OdeState& state, SurgSim::Math::Vector* F,
79  double scale = 1.0) override;
80 
88  virtual void addMass(const SurgSim::Math::OdeState& state, SurgSim::Math::Matrix* M,
89  double scale = 1.0) override;
90 
100  virtual void addDamping(const SurgSim::Math::OdeState& state, SurgSim::Math::Matrix* D,
101  double scale = 1.0) override;
102 
111  virtual void addStiffness(const SurgSim::Math::OdeState& state, SurgSim::Math::Matrix* K,
112  double scale = 1.0) override;
113 
125 
136  virtual void addMatVec(const SurgSim::Math::OdeState& state, double alphaM, double alphaD, double alphaK,
138 
140  const SurgSim::Math::OdeState& state,
141  const SurgSim::Math::Vector& naturalCoordinate) const;
142 
144  const SurgSim::Math::OdeState& state,
145  const SurgSim::Math::Vector& cartesianCoordinate) const override;
146 
147 protected:
152 
156  void computeStiffness(const SurgSim::Math::OdeState& state,
157  Eigen::Matrix<double, 12, 12>* k);
158 
162  void computeMass(const SurgSim::Math::OdeState& state, Eigen::Matrix<double, 12, 12>* m);
163 
165  Eigen::Matrix<double, 12, 1> m_x0;
166 
168  Eigen::Matrix<double, 12, 12> m_R0;
169 
171  Eigen::Matrix<double, 12, 12> m_M;
173  Eigen::Matrix<double, 12, 12> m_MLocal;
175  Eigen::Matrix<double, 12, 12> m_K;
177  Eigen::Matrix<double, 12, 12> m_KLocal;
178 
180  double m_G;
181 
183  double m_restLength;
185  double m_radius;
187  double m_A;
193  double m_Asy, m_Asz;
197  double m_Phi_y, m_Phi_z;
199  double m_Iy, m_Iz;
201  double m_J;
202 };
203 
204 } // namespace Physics
205 
206 } // namespace SurgSim
207 
208 #endif // SURGSIM_PHYSICS_FEM1DELEMENTBEAM_H
Eigen::Matrix< double, 12, 12 > m_M
Mass matrix (in global coordinate frame)
Definition: Fem1DElementBeam.h:171
Definition: DriveElementFromInputBehavior.cpp:27
double m_J
Polar moment of inertia.
Definition: Fem1DElementBeam.h:201
Eigen::Matrix< double, Eigen::Dynamic, Eigen::Dynamic > Matrix
A dynamic size matrix.
Definition: Matrix.h:65
void setRadius(double radius)
Sets the beam's circular cross-section radius.
Definition: Fem1DElementBeam.cpp:56
void computeStiffness(const SurgSim::Math::OdeState &state, Eigen::Matrix< double, 12, 12 > *k)
Computes the beam's stiffness matrix.
Definition: Fem1DElementBeam.cpp:254
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's force vector, mass, stiffness and damping matrices (computed for a given state) in...
Definition: Fem1DElementBeam.cpp:139
Eigen::Matrix< double, 12, 12 > m_R0
Initial rotation matrix for the element.
Definition: Fem1DElementBeam.h:168
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
virtual void addDamping(const SurgSim::Math::OdeState &state, SurgSim::Math::Matrix *D, double scale=1.0) override
Adds the element's damping matrix D (= -df/dv) (computed for a given state) to a complete system damp...
Definition: Fem1DElementBeam.cpp:130
virtual void addForce(const SurgSim::Math::OdeState &state, SurgSim::Math::Vector *F, double scale=1.0) override
Adds the element's force (computed for a given state) to a complete system force vector F (assembly) ...
Definition: Fem1DElementBeam.cpp:113
Eigen::Matrix< double, Eigen::Dynamic, 1 > Vector
A dynamic size column vector.
Definition: Vector.h:67
bool m_haveShear
Does this beam element have shear.
Definition: Fem1DElementBeam.h:189
double m_restLength
Rest length.
Definition: Fem1DElementBeam.h:183
virtual void addMass(const SurgSim::Math::OdeState &state, SurgSim::Math::Matrix *M, double scale=1.0) override
Adds the element's mass matrix M (computed for a given state) to a complete system mass matrix M (ass...
Definition: Fem1DElementBeam.cpp:125
Eigen::Matrix< double, 12, 12 > m_KLocal
Stiffness matrix (in local coordinate frame)
Definition: Fem1DElementBeam.h:177
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: Fem1DElementBeam.cpp:380
virtual void addStiffness(const SurgSim::Math::OdeState &state, SurgSim::Math::Matrix *K, double scale=1.0) override
Adds the element's stiffness matrix K (= -df/dx) (computed for a given state) to a complete system st...
Definition: Fem1DElementBeam.cpp:134
1D FemElement based on a beam volume discretization with a fixed cross section
Definition: Fem1DElementBeam.h:35
Fem1DElementBeam(std::array< size_t, 2 > nodeIds)
Constructor.
Definition: Fem1DElementBeam.cpp:35
double m_Iy
Cross sectional moment of inertia.
Definition: Fem1DElementBeam.h:199
double m_shearFactor
Shear factor (usually 5/8)
Definition: Fem1DElementBeam.h:191
virtual void initialize(const SurgSim::Math::OdeState &state) override
Initializes the FemElement once everything has been set.
Definition: Fem1DElementBeam.cpp:87
double m_radius
radius for a circular Beam
Definition: Fem1DElementBeam.h:185
double m_Asz
Definition: Fem1DElementBeam.h:193
double m_Asy
The shear area in the y and z directions (=0 => no shear) http://en.wikipedia.org/wiki/Timoshenko_bea...
Definition: Fem1DElementBeam.h:193
virtual void addMatVec(const SurgSim::Math::OdeState &state, double alphaM, double alphaD, double alphaK, const SurgSim::Math::Vector &x, SurgSim::Math::Vector *F)
Adds the element's 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: Fem1DElementBeam.cpp:154
void setShearingEnabled(bool enabled)
Enables or disables shearing for the element.
Definition: Fem1DElementBeam.cpp:74
double getRadius() const
Gets the beam's circular cross-section radius.
Definition: Fem1DElementBeam.cpp:64
void computeMass(const SurgSim::Math::OdeState &state, Eigen::Matrix< double, 12, 12 > *m)
Computes the beam's mass matrix.
Definition: Fem1DElementBeam.cpp:185
double m_Phi_z
Definition: Fem1DElementBeam.h:197
void computeInitialRotation(const SurgSim::Math::OdeState &state)
Computes the beam element's initial rotation.
Definition: Fem1DElementBeam.cpp:337
Eigen::Matrix< double, 12, 12 > m_MLocal
Stiffness matrix (in local coordinate frame)
Definition: Fem1DElementBeam.h:173
virtual double getVolume(const SurgSim::Math::OdeState &state) const override
Gets the element's volume based on the input state.
Definition: Fem1DElementBeam.cpp:79
Eigen::Matrix< double, 12, 1 > m_x0
The element's rest state.
Definition: Fem1DElementBeam.h:165
Eigen::Matrix< double, 12, 12 > m_K
Stiffness matrix (in global coordinate frame)
Definition: Fem1DElementBeam.h:175
double m_G
Physical shear modulus G = E/( 2(1+mu) )
Definition: Fem1DElementBeam.h:180
double m_A
Cross sectional area = PI.radius.radius if circular.
Definition: Fem1DElementBeam.h:187
virtual SurgSim::Math::Vector computeCartesianCoordinate(const SurgSim::Math::OdeState &state, const SurgSim::Math::Vector &naturalCoordinate) const
Computes a given natural coordinate in cartesian coordinates.
Definition: Fem1DElementBeam.cpp:362
double m_Phi_y
Shear deformation parameters Phi_y=12.E.Iz/(G.Asy.L^2) or 0 if As?=0 Phi_z=12.E.Iy/(G.Asz.L^2) or 0 if As?=0.
Definition: Fem1DElementBeam.h:197
double m_Iz
Definition: Fem1DElementBeam.h:199
bool getShearingEnabled() const
Gets whether shearing is enabled for the element.
Definition: Fem1DElementBeam.cpp:69