All Classes Namespaces Functions Variables Typedefs Enumerations Enumerator Pages
flexible_body_1S_reference_curve.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  * Created on: Jul 30, 2013
18  * Contact: kilian.grundl@gmail.com
19  */
20 
21 #ifndef FLEXIBLE_BODY_1S_REFERENCE_CURVE_H_
22 #define FLEXIBLE_BODY_1S_REFERENCE_CURVE_H_
23 
24 #include <mbsimFlexibleBody/contours/neutral_contour/contour_1s_neutral_reference_curve.h>
25 
26 #include <mbsimFlexibleBody/flexible_body.h>
27 
28 #include <mbsim/functions/function.h>
29 
30 namespace MBSimFlexibleBody {
31 
36  public:
38  virtual ~ReferenceCurve();
39 
43  virtual void computeReference() = 0;
44 
48  virtual fmatvec::Vec3 computeVecAt(double xi, double Theta, int derXi, int derTheta) = 0;
49 
53  virtual fmatvec::VecV computeC1Positions(const fmatvec::Vec & qRef) = 0;
54 
55  protected:
59  double length;
60  };
61 
67  friend class FlexibleBody1SReferenceCurveFE;
69  friend class funcForWgamma;
70  friend class gethFunc;
71 
72  public:
73  FlexibleBody1SReferenceCurve(const std::string & name, ReferenceCurve * refCurve);
74 
76 
77  virtual void init(MBSim::Element::InitStage stage);
78  virtual void initInfo(fmatvec::Vec q0 = fmatvec::Vec(0,fmatvec::NONINIT), fmatvec::Vec u0 = fmatvec::Vec(0,fmatvec::NONINIT));
79 
80  /*INHERITED INTERFACE*/
81  virtual void BuildElements();
82  virtual void updateh(int k = 0);
83  virtual void updateM();
84  virtual void updatedq();
85  virtual void updatedu();
86  virtual void updateud();
87  virtual void updateqd();
88  virtual void initz();
89  virtual void plot(double t, double dt = 1.);
90 
91  virtual void GlobalVectorContribution(int, const fmatvec::Vec&, fmatvec::Vec&);
92 
93  virtual void GlobalMatrixContribution(int, const fmatvec::Mat3xV&, fmatvec::Mat3xV&);
94 
95  virtual void GlobalMatrixContribution(int, const fmatvec::Mat&, fmatvec::Mat&) {
96  throw MBSim::MBSimError("NOT IMPLEMENTED: " + std::string(__func__));
97  }
98  virtual void GlobalMatrixContribution(int, const fmatvec::SymMat&, fmatvec::SymMat&);
99 
100 // virtual void updateKinematicsForFrame(MBSim::ContourPointData&, MBSim::Frame::Feature, MBSim::Frame*);
101 // virtual void updateJacobiansForFrame(MBSim::ContourPointData&, MBSim::Frame*) {
102 // throw MBSim::MBSimError("NOT IMPLEMENTED: " + std::string(__func__));
103 // }
104 //
105 // virtual void updateKinematicsAtNode(MBSimFlexibleBody::NodeFrame *frame, MBSim::Frame::Feature ff);
106 
107  virtual void setq0(fmatvec::Vec q0_) {
109  qF << q0;
110  }
111  virtual void setu0(fmatvec::Vec u0_) {
112  FlexibleBodyContinuum<double>::setu0(u0_);
113  uF << u0;
114  }
115 
116  /*GETTER and SETTER*/
117  void setLength(double length_) {
118  length = length_;
119  }
120  double getlength() const {
121  return length;
122  }
123  void setDensity(double rho_) {
124  rho = rho_;
125  }
126  double getrho() const {
127  return rho;
128  }
129  void setCrossSectionalArea(double A_) {
130  A = A_;
131  }
132  double getA() const {
133  return A;
134  }
135  void setEModul(double E_) {
136  E = E_;
137  }
138  double getE() const {
139  return E;
140  }
141  void setNu(double nu_) {
142  nu = nu_;
143  }
144  double getNu() const {
145  return nu;
146  }
147  void setIn(double In_) {
148  In = In_;
149  }
150  double getIn() const {
151  return In;
152  }
153  void setIb(double Ib_) {
154  Ib = Ib_;
155  }
156  double getIb() const {
157  return Ib;
158  }
159  void setIt(double It_) {
160  It = It_;
161  }
162  double getIt() const {
163  return It;
164  }
165  void setdPulley(double dPulley_) {
166  dPulley = dPulley_;
167  }
168  double getdPulley() const {
169  return dPulley;
170  }
171 
172  void setNumberElements(int elements_) {
173  elements = elements_;
174  }
175 
176  void setElementOrder(int elementOrder_) {
177  elementOrder = elementOrder_;
178  }
179 
180  void setNodeDofs(int nodeDofs_) {
181  nodeDoFs = nodeDofs_;
182  }
183 
184  void setUseSpatialReferenceKinematics(bool switch_) {
186  }
187 
188  void setThetaDamping(double dTheta_) {
189  dTheta = dTheta_;
190  }
191 
192  void setLockedDofs(std::set<int> lockedDofs_) {
193  lockedDofsFull = lockedDofs_;
194  }
195 
196  void setLambdaqSwitches(const fmatvec::VecV & lambdaqSwitches_) {
197  lambdaqFSwitches = lambdaqSwitches_;
198  }
199 
200  void setLambdauSwitches(const fmatvec::VecV & lambdauSwitches_) {
201  lambdauFSwitches = lambdauSwitches_;
202  }
203 
204  virtual void setelongationActive(bool val) {
205  elongationActive = val;
206  }
207 
208  virtual void setnormalBendingActive(bool val) {
209  normalBendingActive = val;
210  }
211 
212  virtual double gets() {
213  if (not usePT1)
214  return q(0);
215  else
216  return qF(0);
217  }
218 
219 // virtual void updateStateDependentVariables(double t);
220 
224  Contour1sNeutralFlexibleBody1SReferenceCurve* createNeutralPhase(const std::string & contourName);
225 
229  double computePhysicalStrain(double xi);
230 
235 
236  private:
241 
245  int elements;
246 
250  int nodeDoFs;
251 
256 
261 
265  double length;
266 
270  double rho;
271 
275  double A;
276 
280  double E;
281 
285  double nu;
286 
290  double In;
291 
295  double Ib;
296 
300  double It;
301 
305  double dPulley;
306 
310  fmatvec::Vec3 b;
311 
315  double dTheta;
316 
320  std::set<int> lockedDofsFull;
321 
325  fmatvec::VecV qF;
326 
330  fmatvec::VecV uF;
331 
342  fmatvec::VecV lambdaqFSwitches;
343 
348  fmatvec::VecV lambdauFSwitches;
349 
353  fmatvec::VecV lambdaqF;
354 
358  fmatvec::VecV lambdauF;
359 
363  bool usePT1 = false;
364 
368  std::vector<fmatvec::Vec> qElementAll;
369 
373  std::vector<fmatvec::Vec> uElementAll;
374 
379 
384 
389 
393  bool elongationActive = true;
394 
398  bool normalBendingActive = true;
399 
403  double computeG() const {
404  return E / (2 * (1 + nu));
405  }
406 
410  fmatvec::Vec3 computer(double xi, int derXi = 0, int derTheta = 0);
411 
415  fmatvec::Vec3 computedrdqk(double xi, int derXi, int qInd);
416 
420  fmatvec::Vec3 computev(double xi);
421 
425  fmatvec::Vec3 computerRef(double xi, int derXi = 0, int derTheta = 0);
426 
430  fmatvec::SqrMat3 computeARef(double xi, int derXi = 0, int derTheta = 0);
431 
435  fmatvec::Vec3 computenRef(double xi, int derXi = 0, int derTheta = 0);
436 
440  fmatvec::Vec3 computeSElement(int globalDOF, double xi, int derXi);
441 
445  fmatvec::Vec3 computeSqfElement(double xi, int derXi);
446 
450  fmatvec::Vec3 computeSufElement(double xi, int derXi);
451 
455  fmatvec::Mat3xV computeP(double xi, int derXi = 0);
456 
460  fmatvec::Mat3xV computedPdqk(double xi, int qInd);
461 
465  fmatvec::SymMatV integratePTP();
466 
470  fmatvec::SqrMatV integratePTdPdxi();
471 
475  fmatvec::SqrMatV integratePTdPdt();
476 
480  fmatvec::SqrMatV integratePTdPdqk(int qInd);
481 
485  fmatvec::Vec3 integrateForWgammaWnWtau(int qInd);
486 
490  double integrateForWb(int qInd);
491 
496  int findElement(double & xi);
497 
502  std::vector<std::pair<int, int> > getElementAndLocalDoFNo(int globalColumn);
503 
507  int getEleDofs() const;
508  };
509 
510  class gethFunc : public MBSim::Function<fmatvec::Vec(fmatvec::Vec)> {
511  public:
513  refBody(refBody) {
514  }
515 
516  ~gethFunc() {
517  }
518 
519  fmatvec::Vec operator()(const fmatvec::Vec & q) {
520  refBody->setq(q);
521  refBody->updateh(0, 0);
522  return refBody->geth(0);
523  }
524  private:
526  };
527 
528 }
529 
530 #endif /* FLEXIBLE_BODY_1S_REFERENCE_CURVE_H_ */
fmatvec::Vec3 computeSqfElement(double xi, int derXi)
compute the deformation vector, i.e. S*qf, for the ring
Definition: flexible_body_1S_reference_curve.cc:878
std::vector< fmatvec::Vec > uElementAll
vector storing the values of the locked and not locked Dofs generelaized velocities ...
Definition: flexible_body_1S_reference_curve.h:373
std::vector< std::pair< int, int > > getElementAndLocalDoFNo(int globalColumn)
find the elemental DOF Number out of the global DOF Number
Definition: flexible_body_1S_reference_curve.cc:929
fmatvec::Mat3xV computeP(double xi, int derXi=0)
compute the P matrix
Definition: flexible_body_1S_reference_curve.cc:892
bool elongationActive
switch to enable / disable the elongation energy in the h-vector
Definition: flexible_body_1S_reference_curve.h:393
virtual fmatvec::VecV computeC1Positions(const fmatvec::Vec &qRef)=0
compute the only C1-continuous position in case the integral over the FE has to be split ...
bool useSpatialReferenceKinematics
switch to only use kinematical reference deformation
Definition: flexible_body_1S_reference_curve.h:260
std::vector< fmatvec::Vec > qElementAll
vector storing the values of the locked and not locked Dofs
Definition: flexible_body_1S_reference_curve.h:368
double rho
density of the body
Definition: flexible_body_1S_reference_curve.h:270
Contour1sNeutralFlexibleBody1SReferenceCurve * createNeutralPhase(const std::string &contourName)
create a neutral phase as a basis for overlaid contours
Definition: flexible_body_1S_reference_curve.cc:723
double tLastRefUpdate
marker since when the last update has been done on the reference
Definition: flexible_body_1S_reference_curve.h:388
int findElement(double &xi)
find the element number of the given global coordinate
Definition: flexible_body_1S_reference_curve.cc:914
fmatvec::Vec3 b
Binormal of the curve.
Definition: flexible_body_1S_reference_curve.h:310
this function returns basically three values. As these share some values it makes sense to not comput...
Definition: finite_element_1S_reference_curve_functions.h:157
fmatvec::VecV lambdaqF
lambda factors for qFs
Definition: flexible_body_1S_reference_curve.h:353
class that sets up the neutral contour for the body FlexibleBody1SReferenceCurve Not that as interfac...
Definition: contour_1s_neutral_reference_curve.h:37
bool usePT1
switch to use either the lambda PT1 version or not
Definition: flexible_body_1S_reference_curve.h:363
fmatvec::VecV qF
vector for PT1 in positions
Definition: flexible_body_1S_reference_curve.h:325
fmatvec::Vec q
virtual void GlobalVectorContribution(int, const fmatvec::Vec &, fmatvec::Vec &)
insert &#39;local&#39; information in global vectors
Definition: flexible_body_1S_reference_curve.cc:588
virtual void computeReference()=0
get the reference curve at a certain geometric ratio
fmatvec::VecV lambdaqFSwitches
these are the switches to the the lambdaqSwitches The first entrie is always the entro for &quot;s&quot; the ...
Definition: flexible_body_1S_reference_curve.h:342
fmatvec::VecV lambdauFSwitches
these are the switches to the the lambdauSwitches See the explanations for the qSwitches, however for the generalized velocities
Definition: flexible_body_1S_reference_curve.h:348
double E
Youngs modulus of the body.
Definition: flexible_body_1S_reference_curve.h:280
double length
length of belt
Definition: flexible_body_1S_reference_curve.h:59
fmatvec::Vec3 computedrdqk(double xi, int derXi, int qInd)
computes the derivative wrt the generalized positions
Definition: flexible_body_1S_reference_curve.cc:784
int referenceNotUpdated
counter of how many updates have been done since the last time
Definition: flexible_body_1S_reference_curve.h:383
flexible body entirely described within MBSim holding all informations about continuum approximations...
Definition: flexible_body.h:240
this is a interface definition for possible implementations of the reference curve used in the Flexib...
Definition: flexible_body_1S_reference_curve.h:35
Definition: flexible_body_1S_reference_curve.h:510
fmatvec::Vec3 computerRef(double xi, int derXi=0, int derTheta=0)
compute the position and its derivatives w.r.t. xi or theta at position xi and ratio Theta ...
Definition: flexible_body_1S_reference_curve.cc:808
virtual void BuildElements()
references finite element coordinates to assembled coordinates
Definition: flexible_body_1S_reference_curve.cc:429
double Ib
area moment of inertia around binormal axis
Definition: flexible_body_1S_reference_curve.h:295
fmatvec::SymMatV integratePTP()
computes the integral needed for the mass matrix ( P^T P)
double dPulley
distance between pulleys
Definition: flexible_body_1S_reference_curve.h:305
double length
length of ring
Definition: flexible_body_1S_reference_curve.h:265
fmatvec::Vec3 integrateForWgammaWnWtau(int qInd)
computes the integral for the first part of the h vector coming from the potential energy (intForWgam...
double It
area moment of inertia around torsional axis
Definition: flexible_body_1S_reference_curve.h:300
int updateReferenceEvery
number of how often the reference movement should be updated compared to the local movements ...
Definition: flexible_body_1S_reference_curve.h:378
double In
area moment of inertia around normal axis
Definition: flexible_body_1S_reference_curve.h:290
class that implements an one dimensional nonlinear beam that uses a reference curve for the nonlinear...
Definition: flexible_body_1S_reference_curve.h:66
double dTheta
damping coefficient for theta
Definition: flexible_body_1S_reference_curve.h:315
virtual void GlobalMatrixContribution(int, const fmatvec::Mat &, fmatvec::Mat &)
insert &#39;local&#39; information in global matrices
Definition: flexible_body_1S_reference_curve.h:95
Finite-Element class for the body FlexibleBody1SReferenceCurve.
Definition: finite_element_1S_reference_curve.h:37
fmatvec::Vec3 computer(double xi, int derXi=0, int derTheta=0)
compute the position and/or its derivative w.r.t. xi and/or theta
Definition: flexible_body_1S_reference_curve.cc:778
fmatvec::VecV lambdauF
lambda factors for uFs
Definition: flexible_body_1S_reference_curve.h:358
virtual fmatvec::Vec3 computeVecAt(double xi, double Theta, int derXi, int derTheta)=0
computes the vector/point on the surface for specific derivatives
fmatvec::Vec computeNeutralState(const fmatvec::Vec &q0)
compute the neutral state, i.e. h = 0
Definition: flexible_body_1S_reference_curve.cc:757
double computePhysicalStrain(double xi)
compute the physical strain at the Eulerian coordinate xi
Definition: flexible_body_1S_reference_curve.cc:740
bool normalBendingActive
switch to enable / disable the normal bending energy in the h-vector
Definition: flexible_body_1S_reference_curve.h:398
fmatvec::Vec3 computenRef(double xi, int derXi=0, int derTheta=0)
compute the normal vector at the given position
Definition: flexible_body_1S_reference_curve.cc:851
ReferenceCurve * refCurve
the reference curve
Definition: flexible_body_1S_reference_curve.h:240
fmatvec::Vec3 computev(double xi)
compute the velocity and/or its derivative w.r.t. xi and/or theta
Definition: flexible_body_1S_reference_curve.cc:802
int nodeDoFs
DoFs per element.
Definition: flexible_body_1S_reference_curve.h:250
fmatvec::Vec3 computeSufElement(double xi, int derXi)
compute the deformation vector, i.e. S*uf, for the ring
Definition: flexible_body_1S_reference_curve.cc:885
double A
cross section Area of the body
Definition: flexible_body_1S_reference_curve.h:275
int elementOrder
DoFs per element.
Definition: flexible_body_1S_reference_curve.h:255
double nu
Poisson-ratio of the body.
Definition: flexible_body_1S_reference_curve.h:285
double integrateForWb(int qInd)
computes the integral for the second part of the h vector coming from the potential energy (intForWga...
fmatvec::SqrMatV integratePTdPdt()
computes the integral for the second part of the h vector coming from the kinetic energy ( P^T * dPdt...
int getEleDofs() const
get the number of DOFs for one element
Definition: flexible_body_1S_reference_curve.cc:945
fmatvec::SqrMatV integratePTdPdxi()
computes the integral for the first part of the h vector coming from the kinetic energy ( P^T * P&#39;) ...
fmatvec::SqrMatV integratePTdPdqk(int qInd)
computes the integral for the third part of the h vector coming from the kinetic energy ( P^T * dPdqk...
fmatvec::VecV uF
vector for PT1 in velocities
Definition: flexible_body_1S_reference_curve.h:330
double computeG() const
computes the shear modulus G
Definition: flexible_body_1S_reference_curve.h:403
fmatvec::Vec3 computeSElement(int globalDOF, double xi, int derXi)
compute the column of the elements S matrix
Definition: flexible_body_1S_reference_curve.cc:866
int elements
number of elements used
Definition: flexible_body_1S_reference_curve.h:245
fmatvec::Mat3xV computedPdqk(double xi, int qInd)
compute the P matrix derived wrt the generalized position
Definition: flexible_body_1S_reference_curve.cc:902
fmatvec::SqrMat3 computeARef(double xi, int derXi=0, int derTheta=0)
compute the local transformation Matrix A and its derivatives w.r.t. xi or theta at position xi and r...
Definition: flexible_body_1S_reference_curve.cc:823
std::set< int > lockedDofsFull
this vector saves the locked DOFs of the system
Definition: flexible_body_1S_reference_curve.h:320
fmatvec::Vec q0

Impressum / Disclaimer / Datenschutz Generated by doxygen 1.8.5 Valid HTML