All Classes Namespaces Functions Variables Typedefs Enumerations Enumerator Pages
finite_element_2s_13_disk.h
1 /* Copyright (C) 2004-2015 MBSim Development Team
2  *
3  * This library is free software; you can redistribute it and/or
4  * modify it under the terms of the GNU Lesser General Public
5  * License as published by the Free Software Foundation; either
6  * version 2.1 of the License, or (at your option) any later version.
7  *
8  * This library is distributed in the hope that it will be useful,
9  * but WITHOUT ANY WARRANTY; without even the implied warranty of
10  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
11  * Lesser General Public License for more details.
12  *
13  * You should have received a copy of the GNU Lesser General Public
14  * License along with this library; if not, write to the Free Software
15  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
16  *
17  * Contact: thorsten.schindler@mytum.de
18  */
19 
20 #ifndef _FINITE_ELEMENT_2S_13_DISK_H_
21 #define _FINITE_ELEMENT_2S_13_DISK_H_
22 
23 #include "mbsimFlexibleBody/discretization_interface.h"
24 #include <cmath>
25 
26 namespace MBSimFlexibleBody {
27 
39  public:
46  FiniteElement2s13Disk(double E_, double nu_, double rho_);
47 
51  virtual ~FiniteElement2s13Disk() { }
52 
53  /* INTERFACE OF DISCRETIZATIONINTERFACE */
54  virtual const fmatvec::SymMat& getM() const { return M; }
55  virtual const fmatvec::Vec& geth() const;
56  virtual const fmatvec::SqrMat& getdhdq() const;
57  virtual const fmatvec::SqrMat& getdhdu() const;
58  virtual int getqSize() const { return RefDofs + 4*NodeDofs; }
59  virtual int getuSize() const { return RefDofs + 4*NodeDofs; }
60  virtual void computeM(const fmatvec::Vec& q);
61  virtual void computeh(const fmatvec::Vec& q,const fmatvec::Vec& u);
62  virtual void computedhdz(const fmatvec::Vec& q,const fmatvec::Vec& u);
63  virtual double computeKineticEnergy(const fmatvec::Vec& q,const fmatvec::Vec& u);
64  virtual double computeGravitationalEnergy(const fmatvec::Vec& q);
65  virtual double computeElasticEnergy(const fmatvec::Vec& q);
66  virtual fmatvec::Vec3 getPosition(const fmatvec::Vec& qElement, const fmatvec::Vec2 &s);
67  virtual fmatvec::SqrMat3 getOrientation(const fmatvec::Vec& qElement, const fmatvec::Vec2 &s);
68  virtual fmatvec::Vec3 getVelocity (const fmatvec::Vec& qElement, const fmatvec::Vec& qpElement, const fmatvec::Vec2 &s);
69  virtual fmatvec::Vec3 getAngularVelocity(const fmatvec::Vec& qElement, const fmatvec::Vec& qpElement, const fmatvec::Vec2 &s);
70  virtual fmatvec::Mat getJacobianOfMotion(const fmatvec::Vec& qElement, const fmatvec::Vec2 &s);
71  /***************************************************/
72 
73  /* GETTER / SETTER */
74  const fmatvec::SymMat& getK() const { return K; }
75  void setEModul(double E_) { E = E_; }
76  void setPoissonRatio(double nu_) { nu = nu_; }
77  void setDensity(double rho_) { rho = rho_; }
78  void setShearCorrectionFactor(double alphaS_) { alphaS = alphaS_; }
79  /***************************************************/
80 
87  void computeConstantSystemMatrices(const fmatvec::Vec &NodeCoordinates, double d1, double d2);
88 
98  fmatvec::Vector<fmatvec::Fixed<6>, double> getPositions(const fmatvec::Vec &NodeCoordinates, const fmatvec::Vec &qElement, const fmatvec::Vec2 &s, double d1, double d2);
99 
100  fmatvec::Vector<fmatvec::Fixed<6>, double> getVelocities(const fmatvec::Vec &NodeCoordinates, const fmatvec::Vec &qElement, const fmatvec::Vec &qpElement, const fmatvec::Vec2 &s, double d1, double d2);
101 
107  fmatvec::Mat JGeneralized(const fmatvec::Vec &NodeCoordinates, const fmatvec::Vec2 &s);
108 
109  private:
113  double E;
114 
118  double nu;
119 
123  double G;
124 
128  double rho;
129 
133  double alphaS;
134 
138  int RefDofs;
139 
143  int NodeDofs;
144 
148  int Nodes;
149 
154  };
155 
156 }
157 
158 #endif /* _FINITE_ELEMENT_2S_13_RCM_H_ */
virtual const fmatvec::SqrMat & getdhdu() const
Definition: finite_element_2s_13_disk.cc:440
void computeConstantSystemMatrices(const fmatvec::Vec &NodeCoordinates, double d1, double d2)
computes mass and stiffness matrix
Definition: finite_element_2s_13_disk.cc:37
double G
shear modulus
Definition: finite_element_2s_13_disk.h:123
virtual int getqSize() const
Definition: finite_element_2s_13_disk.h:58
virtual double computeKineticEnergy(const fmatvec::Vec &q, const fmatvec::Vec &u)
compute kinetic energy
Definition: finite_element_2s_13_disk.cc:456
virtual const fmatvec::Vec & geth() const
Definition: finite_element_2s_13_disk.cc:432
double alphaS
shear correction factor
Definition: finite_element_2s_13_disk.h:133
double rho
density
Definition: finite_element_2s_13_disk.h:128
virtual void computedhdz(const fmatvec::Vec &q, const fmatvec::Vec &u)
compute Jacobian for implicit integration
Definition: finite_element_2s_13_disk.cc:452
int NodeDofs
elastic dof per node
Definition: finite_element_2s_13_disk.h:143
int RefDofs
reference dof
Definition: finite_element_2s_13_disk.h:138
virtual const fmatvec::SymMat & getM() const
Definition: finite_element_2s_13_disk.h:54
double nu
Poisson ratio.
Definition: finite_element_2s_13_disk.h:118
virtual int getuSize() const
Definition: finite_element_2s_13_disk.h:59
virtual ~FiniteElement2s13Disk()
destructor
Definition: finite_element_2s_13_disk.h:51
fmatvec::Mat JGeneralized(const fmatvec::Vec &NodeCoordinates, const fmatvec::Vec2 &s)
compute Jacobian of contact description at contour point
Definition: finite_element_2s_13_disk.cc:354
fmatvec::SymMat M
mass and stiffness matrix
Definition: finite_element_2s_13_disk.h:153
int Nodes
number of nodes
Definition: finite_element_2s_13_disk.h:148
double E
Young&#39;s modulus.
Definition: finite_element_2s_13_disk.h:113
virtual double computeGravitationalEnergy(const fmatvec::Vec &q)
compute gravitational energy
Definition: finite_element_2s_13_disk.cc:460
fmatvec::Vector< fmatvec::Fixed< 6 >, double > getPositions(const fmatvec::Vec &NodeCoordinates, const fmatvec::Vec &qElement, const fmatvec::Vec2 &s, double d1, double d2)
Definition: finite_element_2s_13_disk.cc:274
virtual double computeElasticEnergy(const fmatvec::Vec &q)
compute elastic energy
Definition: finite_element_2s_13_disk.cc:464
FE for Reissner-Mindlin Plate using MFR.
Definition: finite_element_2s_13_disk.h:38
discretization interface for flexible systems
Definition: discretization_interface.h:36
virtual void computeM(const fmatvec::Vec &q)
compute mass matrix
Definition: finite_element_2s_13_disk.cc:444
virtual void computeh(const fmatvec::Vec &q, const fmatvec::Vec &u)
compute smooth right hand side
Definition: finite_element_2s_13_disk.cc:448
virtual const fmatvec::SqrMat & getdhdq() const
Definition: finite_element_2s_13_disk.cc:436
FiniteElement2s13Disk(double E_, double nu_, double rho_)
constructor
Definition: finite_element_2s_13_disk.cc:35

Impressum / Disclaimer / Datenschutz Generated by doxygen 1.8.5 Valid HTML