All Classes Namespaces Functions Variables Typedefs Enumerations Enumerator Pages
flexible_body_1s_33_rcm.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 _FLEXIBLE_BODY_1S_33_RCM_H_
21 #define _FLEXIBLE_BODY_1S_33_RCM_H_
22 
23 #include "mbsimFlexibleBody/flexible_body/flexible_body_1s.h"
24 #include "mbsimFlexibleBody/pointer.h"
25 
26 namespace MBSimFlexibleBody {
27 
41  public:
47  FlexibleBody1s33RCM(const std::string &name="",bool openStructure=false);
48 
49  /* INHERITED INTERFACE OF FLEXIBLE BODY */
50  virtual void BuildElements();
51 
52  virtual void GlobalVectorContribution(int n, const fmatvec::Vec& locVec, fmatvec::Vec& gloVec);
53  virtual void GlobalMatrixContribution(int n, const fmatvec::Mat& locMat, fmatvec::Mat& gloMat);
54  virtual void GlobalMatrixContribution(int n, const fmatvec::SymMat& locMat, fmatvec::SymMat& gloMat);
55 
56  virtual void exportPositionVelocity(const std::string & filenamePos, const std::string & filenameVel = std::string( ), const int & deg = 3, const bool & writePsFile = false);
57  virtual void importPositionVelocity(const std::string & filenamePos, const std::string & filenameVel = std::string( ));
58 
59  virtual void updatePositions(Frame1s* frame);
60  virtual void updateVelocities(Frame1s* frame);
61  virtual void updateAccelerations(Frame1s* frame);
62  virtual void updateJacobians(Frame1s* frame, int j=0);
63  virtual void updateGyroscopicAccelerations(Frame1s* frame);
64 
65  virtual void updatePositions(NodeFrame* frame);
66  virtual void updateVelocities(NodeFrame* frame);
67  virtual void updateAccelerations(NodeFrame* frame);
68  virtual void updateJacobians(NodeFrame* frame, int j=0);
69  virtual void updateGyroscopicAccelerations(NodeFrame* frame);
70 
71  virtual fmatvec::Vec3 getAngles(double s);
72  /***************************************************/
73 
74  /* INHERITED INTERFACE OF OBJECT */
75  virtual void init(InitStage stage);
76  /***************************************************/
77 
78  /* INHERITED INTERFACE OF ELEMENT */
79  virtual std::string getType() const { return "FlexibleBody1s33RCM"; }
80  void initializeUsingXML(xercesc::DOMElement * element);
81  /***************************************************/
82 
83  /* GETTER / SETTER */
84  void setNumberElements(int n);
85  int getNumberElements(){ return Elements; }
86  void setEGModuls(double E_, double G_) { E = E_; G = G_; }
87  void setDensity(double rho_) { rho = rho_; }
88  void setCrossSectionalArea(double A_) { A = A_; }
89  void setMomentsInertia(double I1_, double I2_, double I0_) { I1 = I1_; I2 = I2_; I0 = I0_; }
90  void setMaterialDamping(double epstD_,double k0D_);
91 
92  void setGauss(int nGauss_) { nGauss = nGauss_; }
93  void setCurlRadius(double R1_, double R2_);
94  void setLehrDamping(double epstL_, double k0L_);
95  /***************************************************/
96 
102 
108 
113  double computePhysicalStrain(const double sGlobal);
114 
118  void initInfo();
119 
120  private:
124  RevCardanPtr angle;
125 
129  int Elements;
130 
134  double l0;
135 
139  double E, G;
140 
144  double A;
145 
149  double I1, I2, I0;
150 
154  double rho;
155 
159  double R1, R2;
160 
164  double epstD, k0D, epstL, k0L;
165 
170 
174  int nGauss;
175 
182  void BuildElement(const double& sGlobal, double& sLocal, int& currentElement);
183  };
184 
185 }
186 
187 #endif /* _FLEXIBLE_BODY_1S_33_RCM_H_ */
virtual void BuildElements()
references finite element coordinates to assembled coordinates
Definition: flexible_body_1s_33_rcm.cc:51
bool initialised
initialised FLAG
Definition: flexible_body_1s_33_rcm.h:169
double A
area of cross-section
Definition: flexible_body_1s_33_rcm.h:144
void BuildElement(const double &sGlobal, double &sLocal, int &currentElement)
detect current finite element
Definition: flexible_body_1s_33_rcm.cc:323
virtual void exportPositionVelocity(const std::string &filenamePos, const std::string &filenameVel=std::string(), const int &deg=3, const bool &writePsFile=false)
interpolates the position and optional the velocity coordinates of the flexible body with Nurbs-packa...
Definition: flexible_body_1s_33_rcm.cc:390
double computePhysicalStrain(const double sGlobal)
compute the physical strain of the element
Definition: flexible_body_1s_33_rcm.cc:305
tbd
Definition: flexible_body_1s.h:33
Definition: frame_1s.h:27
virtual void GlobalVectorContribution(int n, const fmatvec::Vec &locVec, fmatvec::Vec &gloVec)
insert 'local' information in global vectors
Definition: flexible_body_1s_33_rcm.cc:97
int nGauss
number of Gauss points for rotational kinetic energy
Definition: flexible_body_1s_33_rcm.h:174
spatial beam using Redundant Coordinate Method (RCM)
Definition: flexible_body_1s_33_rcm.h:40
int Elements
number of elements
Definition: flexible_body_1s_33_rcm.h:129
fmatvec::Vector< fmatvec::Fixed< 6 >, double > getVelocities(double x)
compute velocities and differentiated angles at Lagrangian coordinate in local FE coordinates ...
Definition: flexible_body_1s_33_rcm.cc:298
bool openStructure
flag for open (cantilever beam) or closed (rings) structures
Definition: flexible_body_1s.h:73
RevCardanPtr angle
angle parametrisation
Definition: flexible_body_1s_33_rcm.h:124
double l0
length of entire beam and finite elements
Definition: flexible_body_1s_33_rcm.h:134
fmatvec::Vector< fmatvec::Fixed< 6 >, double > getPositions(double x)
compute positions and angle at Lagrangian coordinate in local FE coordinates
Definition: flexible_body_1s_33_rcm.cc:291
FlexibleBody1s33RCM(const std::string &name="", bool openStructure=false)
constructor
Definition: flexible_body_1s_33_rcm.cc:48
double R1
radius of undeformed shape
Definition: flexible_body_1s_33_rcm.h:159
std::vector< Frame * > frame
virtual void GlobalMatrixContribution(int n, const fmatvec::Mat &locMat, fmatvec::Mat &gloMat)
insert &#39;local&#39; information in global matrices
Definition: flexible_body_1s_33_rcm.cc:109
double rho
density
Definition: flexible_body_1s_33_rcm.h:154
virtual void importPositionVelocity(const std::string &filenamePos, const std::string &filenameVel=std::string())
imports the interpolated position and optional the velocity files (created with exportPositionVelocit...
Definition: flexible_body_1s_33_rcm.cc:454
double epstD
damping
Definition: flexible_body_1s_33_rcm.h:164
cartesian frame on nodes of flexible bodies
Definition: node_frame.h:31
double E
elastic modules
Definition: flexible_body_1s_33_rcm.h:139
double I1
area moment of inertia
Definition: flexible_body_1s_33_rcm.h:149
void initInfo()
Definition: flexible_body_1s_33_rcm.cc:312

Impressum / Disclaimer / Datenschutz Generated by doxygen 1.8.5 Valid HTML