All Classes Namespaces Functions Variables Typedefs Enumerations Enumerator Pages
finite_element_2s_13_mfr_mindlin.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_MFR_MINDLIN_H_
21 #define _FINITE_ELEMENT_2S_13_MFR_MINDLIN_H_
22 
23 #include "mbsimFlexibleBody/discretization_interface.h"
24 #include<cmath>
25 
26 namespace MBSimFlexibleBody {
27 
37  public:
46  FiniteElement2s13MFRMindlin(double E_,double nu_,double rho_,double d0_,double d1_,double d2_,const fmatvec::Vec &NodeCoordinates);
47 
52 
53  /* INTERFACE OF DISCRETIZATIONINTERFACE */
54  virtual const fmatvec::Vec& geth() const;
55  virtual const fmatvec::SqrMat& getdhdq() const;
56  virtual const fmatvec::SqrMat& getdhdu() const;
57  inline int getqSize() const { return RefDofs + 4*NodeDofs; }
58  inline int getuSize() const { return RefDofs + 4*NodeDofs; }
59  virtual void computeM(const fmatvec::Vec& q);
60  virtual void computeh(const fmatvec::Vec& q,const fmatvec::Vec& u);
61  virtual void computedhdz(const fmatvec::Vec& q,const fmatvec::Vec& u);
62  virtual double computeKineticEnergy(const fmatvec::Vec& q,const fmatvec::Vec& u);
63  virtual double computeGravitationalEnergy(const fmatvec::Vec& q);
64  virtual double computeElasticEnergy(const fmatvec::Vec& q);
65  virtual fmatvec::Vec3 getPosition(const fmatvec::Vec& qElement, const fmatvec::Vec2 &s);
66  virtual fmatvec::SqrMat3 getOrientation(const fmatvec::Vec& qElement, const fmatvec::Vec2 &s);
67  virtual fmatvec::Vec3 getVelocity (const fmatvec::Vec& qElement, const fmatvec::Vec& qpElement, const fmatvec::Vec2 &s);
68  virtual fmatvec::Vec3 getAngularVelocity(const fmatvec::Vec& qElement, const fmatvec::Vec& qpElement, const fmatvec::Vec2 &s);
69  virtual fmatvec::Mat getJacobianOfMotion(const fmatvec::Vec& qElement, const fmatvec::Vec2 &s);
70  /***************************************************/
71 
72  /* GETTER / SETTER */
73  const fmatvec::SymMat& getM() const;
74  const fmatvec::SymMat& getK() const { return *K; }
75  const fmatvec::SymMat& getM_RR() const { return *M_RR; }
76  const fmatvec::Mat& getN_compl() const { return *N_compl; }
77  const fmatvec::SqrMat& getN_ij(int i, int j) const { return *(N_ij[i][j]); }
78  const fmatvec::RowVec& getNR_ij(int i, int j) const { return *(NR_ij[i][j]); }
79  const fmatvec::Vec& getR_compl() const { return *R_compl; }
80  const fmatvec::SymMat& getR_ij() const { return *R_ij; }
81  void setEModul(double E_) { E = E_; }
82  void setPoissonRatio(double nu_) { nu = nu_; }
83  void setDensity(double rho_) { rho = rho_; }
84  void setShearCorrectionFactor(double alphaS_) { alphaS = alphaS_; }
85  /***************************************************/
86 
87  /* Freeer */
88  void freeK() { delete K; K=0; }
89  void freeM_RR() { delete M_RR; M_RR=0; }
90  void freeN_ij(int i, int j) { delete N_ij[i][j]; N_ij[i][j]=0; }
91  void freeN_compl() { delete N_compl; N_compl=0; }
92  void freeNR_ij(int i, int j) { delete NR_ij[i][j]; NR_ij[i][j]=0; }
93  void freeR_compl() { delete R_compl; R_compl=0; }
94  void freeR_ij() { delete R_ij; R_ij=0; }
95  /***************************************************/
96 
100  void computeStiffnessMatrix();
101 
105  void computeM_RR();
106 
110  void computeN_compl();
111 
116  void computeN_ij(int i, int j);
117  void computeN_11();
118  void computeN_12();
119  void computeN_13();
120  void computeN_21();
121  void computeN_22();
122  void computeN_23();
123  void computeN_31();
124  void computeN_32();
125  void computeN_33();
126 
131  void computeNR_ij(int i, int j);
132  void computeNR_11();
133  void computeNR_12();
134  void computeNR_13();
135  void computeNR_21();
136  void computeNR_22();
137  void computeNR_23();
138  void computeNR_31();
139  void computeNR_32();
140  void computeNR_33();
141 
145  void computeR_compl();
146 
150  void computeR_ij();
151 
161  fmatvec::Vector<fmatvec::Fixed<6>, double> getPositions(const fmatvec::Vec &NodeCoordinates, const fmatvec::Vec &qElement, const fmatvec::Vec2 &s, double d1, double d2);
162 
163  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);
164 
170  fmatvec::Mat JGeneralized(const fmatvec::Vec &NodeCoordinates,const fmatvec::Vec2 &s);
171 
172  private:
176  double E;
177 
181  double nu;
182 
186  double G;
187 
191  double d0, d1, d2;
192 
196  double rho;
197 
201  double alphaS;
202 
206  int RefDofs;
207 
211  int NodeDofs;
212 
216  int Nodes;
217 
222 
227 
232 
237 
242 
247 
252 
257  };
258 
259 }
260 
261 #endif /* _FINITE_ELEMENT_2S_13_MFR_MINDLIN_H_ */
262 
fmatvec::Vec * R_compl
part of the mass matrix
Definition: finite_element_2s_13_mfr_mindlin.h:251
const fmatvec::SymMat & getM() const
Definition: finite_element_2s_13_mfr_mindlin.cc:85
double G
shear modulus
Definition: finite_element_2s_13_mfr_mindlin.h:186
fmatvec::Mat JGeneralized(const fmatvec::Vec &NodeCoordinates, const fmatvec::Vec2 &s)
compute Jacobian of contact description at contour point
Definition: finite_element_2s_13_mfr_mindlin.cc:145
virtual const fmatvec::Vec & geth() const
Definition: finite_element_2s_13_mfr_mindlin.cc:89
virtual const fmatvec::SqrMat & getdhdu() const
Definition: finite_element_2s_13_mfr_mindlin.cc:97
fmatvec::SymMat * M_RR
part of the mass matrix
Definition: finite_element_2s_13_mfr_mindlin.h:231
fmatvec::SymMat * K
stiffness matrix
Definition: finite_element_2s_13_mfr_mindlin.h:226
virtual double computeElasticEnergy(const fmatvec::Vec &q)
compute elastic energy
Definition: finite_element_2s_13_mfr_mindlin.cc:121
int getqSize() const
Definition: finite_element_2s_13_mfr_mindlin.h:57
FiniteElement2s13MFRMindlin(double E_, double nu_, double rho_, double d0_, double d1_, double d2_, const fmatvec::Vec &NodeCoordinates)
constructor
Definition: finite_element_2s_13_mfr_mindlin.cc:33
int NodeDofs
elastic dof per node
Definition: finite_element_2s_13_mfr_mindlin.h:211
void computeStiffnessMatrix()
computes stiffnes matrix
Definition: K.cc:29
virtual double computeKineticEnergy(const fmatvec::Vec &q, const fmatvec::Vec &u)
compute kinetic energy
Definition: finite_element_2s_13_mfr_mindlin.cc:113
void computeR_ij()
computes a part of the mass matrix
Definition: R_ij.cc:29
fmatvec::Vector< fmatvec::Fixed< 6 >, double > getPositions(const fmatvec::Vec &NodeCoordinates, const fmatvec::Vec &qElement, const fmatvec::Vec2 &s, double d1, double d2)
virtual const fmatvec::SqrMat & getdhdq() const
Definition: finite_element_2s_13_mfr_mindlin.cc:93
void computeR_compl()
computes a part of the mass matrix
Definition: R_compl.cc:29
void computeN_compl()
computes a part of the mass matrix
Definition: N_compl.cc:29
double nu
Poisson ratio.
Definition: finite_element_2s_13_mfr_mindlin.h:181
virtual double computeGravitationalEnergy(const fmatvec::Vec &q)
compute gravitational energy
Definition: finite_element_2s_13_mfr_mindlin.cc:117
fmatvec::RowVec * NR_ij[3][3]
part of the mass matrix
Definition: finite_element_2s_13_mfr_mindlin.h:246
double rho
density
Definition: finite_element_2s_13_mfr_mindlin.h:196
int Nodes
number of nodes
Definition: finite_element_2s_13_mfr_mindlin.h:216
fmatvec::Mat * N_compl
part of the mass matrix
Definition: finite_element_2s_13_mfr_mindlin.h:236
virtual void computeh(const fmatvec::Vec &q, const fmatvec::Vec &u)
compute smooth right hand side
Definition: finite_element_2s_13_mfr_mindlin.cc:105
void computeN_ij(int i, int j)
computes a part of the mass matrix
Definition: finite_element_2s_13_mfr_mindlin.cc:68
void computeM_RR()
computes a part of the mass matrix
Definition: M_RR.cc:29
virtual void computedhdz(const fmatvec::Vec &q, const fmatvec::Vec &u)
compute Jacobian for implicit integration
Definition: finite_element_2s_13_mfr_mindlin.cc:109
void computeNR_ij(int i, int j)
computes a part of the mass matrix
Definition: finite_element_2s_13_mfr_mindlin.cc:149
fmatvec::Vec NodeCoordinates
radial and azimuthal coordinates of corner nodes
Definition: finite_element_2s_13_mfr_mindlin.h:221
fmatvec::SqrMat * N_ij[3][3]
part of the mass matrix
Definition: finite_element_2s_13_mfr_mindlin.h:241
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_mfr_mindlin.cc:101
FE for Reissner-Mindlin Plate using MFR.
Definition: finite_element_2s_13_mfr_mindlin.h:36
virtual ~FiniteElement2s13MFRMindlin()
destructor
Definition: finite_element_2s_13_mfr_mindlin.cc:55
fmatvec::SymMat * R_ij
part of the mass matrix
Definition: finite_element_2s_13_mfr_mindlin.h:256
int getuSize() const
Definition: finite_element_2s_13_mfr_mindlin.h:58
double E
Young&#39;s modulus.
Definition: finite_element_2s_13_mfr_mindlin.h:176
int RefDofs
reference dof
Definition: finite_element_2s_13_mfr_mindlin.h:206
double alphaS
shear correction factor
Definition: finite_element_2s_13_mfr_mindlin.h:201
double d0
geometric factors of the disk
Definition: finite_element_2s_13_mfr_mindlin.h:191

Impressum / Disclaimer / Datenschutz Generated by doxygen 1.8.5 Valid HTML