All Classes Namespaces Functions Variables Typedefs Enumerations Enumerator Pages
object.h
1 /* Copyright (C) 2004-2014 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: martin.o.foerg@googlemail.com
18  */
19 
20 #ifndef _OBJECT_H_
21 #define _OBJECT_H_
22 
23 #include "mbsim/element.h"
24 //#include "mbsim/object_interface.h"
25 
26 namespace MBSim {
27 
28  class DynamicSystem;
29  class Link;
30 
41  //class Object : public Element, public ObjectInterface {
42  class Object : public Element {
43  public:
47  Object(const std::string &name);
48 
52  virtual ~Object();
53 
54  /* INHERITED INTERFACE OF OBJECTINTERFACE */
55  virtual void updateT(double t) {};
56  virtual void updateh(double t, int j=0) {};
57  virtual void updateM(double t, int i=0) {};
58  virtual void updatedhdz(double t);
59  virtual void updatedq(double t, double dt) { qd = T * u * dt; }
60  virtual void updatedu(double t, double dt) { ud[0] = slvLLFac(LLM[0], h[0] * dt + r[0]); }
61  virtual void updateud(double t, int i=0) { ud[i] = slvLLFac(LLM[i], h[i] + r[i]); }
62  virtual void updateqd(double t) { qd = T * u; }
63  virtual void updatezd(double t) { updateqd(t); updateud(t); updatexd(t); }
64  virtual void sethSize(int hSize_, int i=0);
65  virtual int gethSize(int i=0) const { return hSize[i]; }
66  virtual int getqSize() const { return qSize; }
67  virtual int getuSize(int i=0) const { return uSize[i]; }
68  virtual void calcqSize() {};
69  virtual void calcuSize(int j) {};
70  //virtual int getqInd(DynamicSystem* sys);
71  virtual int getuInd(int i=0) { return uInd[i]; }
72  //virtual int getuInd(DynamicSystem* sys, int i=0);
73  virtual void setqInd(int qInd_) { qInd = qInd_; }
74  virtual void setuInd(int uInd_, int i=0) { uInd[i] = uInd_; }
75  //virtual int gethInd(DynamicSystem* sys,int i=0);
76  virtual const fmatvec::Vec& getq() const { return q; };
77  virtual const fmatvec::Vec& getu() const { return u; };
78  virtual H5::GroupBase *getPlotGroup() { return plotGroup; }
79  virtual PlotFeatureStatus getPlotFeature(PlotFeature fp) { return Element::getPlotFeature(fp); };
80  virtual PlotFeatureStatus getPlotFeatureForChildren(PlotFeature fp) { return Element::getPlotFeatureForChildren(fp); };
81  virtual void updateStateDependentVariables(double t) = 0;
82  virtual void updateStateDerivativeDependentVariables(double t) {};
83  virtual void updateJacobians(double t, int j=0) = 0;
84  virtual void updatehInverseKinetics(double t, int i=0) {};
85  /*******************************************************/
86 
87  /* INHERITED INTERFACE OF ELEMENT */
88  virtual void plot(double t, double dt = 1);
89  virtual void closePlot();
90  virtual std::string getType() const { return "Object"; }
91  //virtual void setDynamicSystemSolver(DynamicSystemSolver *sys);
92  /*******************************************************/
93 
94  virtual void updatedx(double t, double dt) {}
95  virtual void updatexd(double t) {}
96  virtual void calcxSize() { xSize = 0; }
97  virtual const fmatvec::Vec& getx() const { return x; }
98  virtual fmatvec::Vec& getx() { return x; }
99  virtual void setxInd(int xInd_) { xInd = xInd_; };
100  virtual int getxSize() const { return xSize; }
101  virtual void updatexRef(const fmatvec::Vec& ref);
102  virtual void updatexdRef(const fmatvec::Vec& ref);
103 
108  virtual void updateqRef(const fmatvec::Vec& qRef);
109 
114  virtual void updateqdRef(const fmatvec::Vec& qdRef);
115 
120  virtual void updateuRef(const fmatvec::Vec& uRef);
121 
126  virtual void updateuallRef(const fmatvec::Vec& uallRef);
127 
132  virtual void updateudRef(const fmatvec::Vec& udRef, int i=0);
133 
138  virtual void updateudallRef(const fmatvec::Vec& udallRef, int i=0);
139 
145  virtual void updatehRef(const fmatvec::Vec& hRef, int i=0);
146 
152  virtual void updatedhdqRef(const fmatvec::Mat& dhdqRef, int i=0);
153 
159  virtual void updatedhduRef(const fmatvec::SqrMat& dhduRef, int i=0);
160 
166  virtual void updatedhdtRef(const fmatvec::Vec& dhdtRef, int i=0);
167 
172  virtual void updaterRef(const fmatvec::Vec& ref, int i=0);
173 
178  virtual void updateTRef(const fmatvec::Mat &ref);
179 
185  virtual void updateMRef(const fmatvec::SymMat &ref, int i=0);
186 
192  virtual void updateLLMRef(const fmatvec::SymMat &ref, int i=0);
193 
197  virtual void init(InitStage stage);
198 
202  virtual void initz();
203 
207  virtual void writez(H5::GroupBase *group);
208 
212  virtual void readz0(H5::GroupBase *group);
213 
217  virtual void facLLM(int i=0) { LLM[i] = facLL(M[i]); }
218 
223  virtual void calcSize(int j) {}
224 
228  virtual double computeKineticEnergy() { return 0.5*u.T()*M[0]*u; }
229 
233  virtual double computePotentialEnergy() { return 0; }
234 
238  virtual void setUpInverseKinetics() {}
239  /*******************************************************/
240 
241  /* GETTER / SETTER */
242  void setqSize(int qSize_) { qSize = qSize_; }
243  void setuSize(int uSize_, int i=0) { uSize[i] = uSize_; }
244  int getzSize() const { return qSize + uSize[0]; }
245 
246  virtual void sethInd(int hInd_, int i=0);
247  int gethInd(int i=0) { return hInd[i]; }
248 
249  const fmatvec::Vec& geth(int i=0) const { return h[i]; };
250  fmatvec::Vec& geth(int i=0) { return h[i]; };
251  const fmatvec::Vec& getr(int i=0) const { return r[i]; };
252  fmatvec::Vec& getr(int i=0) { return r[i]; };
253  const fmatvec::SymMat& getM(int i=0) const { return M[i]; };
254  fmatvec::SymMat& getM(int i=0) { return M[i]; };
255  const fmatvec::Mat& getT() const { return T; };
256  fmatvec::Mat& getT() { return T; };
257  const fmatvec::SymMat& getLLM(int i=0) const { return LLM[i]; };
258  fmatvec::SymMat& getLLM(int i=0) { return LLM[i]; };
259 
260  fmatvec::Vec& getq() { return q; };
261  fmatvec::Vec& getu() { return u; };
262 
263  const fmatvec::Vec& getq0() const { return q0; };
264  const fmatvec::Vec& getu0() const { return u0; };
265  fmatvec::Vec& getq0() { return q0; };
266  fmatvec::Vec& getu0() { return u0; };
267 
268  const fmatvec::Vec& getqd() const { return qd; };
269  const fmatvec::Vec& getud(int i=0) const { return ud[i]; };
270  fmatvec::Vec& getqd() { return qd; };
271  fmatvec::Vec& getud(int i=0) { return ud[i]; };
272 
273  void setq(const fmatvec::Vec &q_) { q = q_; }
274  void setu(const fmatvec::Vec &u_) { u = u_; }
275 
276  /*******************************************************/
277 
278  virtual void initializeUsingXML(xercesc::DOMElement *element);
279  virtual xercesc::DOMElement* writeXMLFile(xercesc::DOMNode *element);
280 
281  protected:
285  int qSize;
286 
290  int uSize[2];
291 
295  int hSize[2];
296 
300  int qInd, uInd[2], hInd[2];
301 
305  int xSize, xInd;
306 
310  fmatvec::Vec q, u, uall, x;
311 
315  fmatvec::Vec q0, u0, x0;
316 
320  fmatvec::Vec qd, ud[2], udall[2], xd;
321 
325  fmatvec::Vec h[2], r[2];
326 
327  fmatvec::Mat W[2], V[2];
328 
333  fmatvec::SqrMat dhdu;
334  fmatvec::Vec dhdt;
335 
340 
345 
350  };
351 
352 }
353 
354 #endif /* _OBJECT_H_ */
355 
int qInd
indices of positions, velocities, right hand side
Definition: object.h:300
Object(const std::string &name)
constructor
Definition: object.cc:35
int xSize
size and local index of order one parameters
Definition: object.h:305
virtual void updateTRef(const fmatvec::Mat &ref)
references to linear transformation matrix between differentiated positions and velocities of dynamic...
Definition: object.cc:216
virtual double computeKineticEnergy()
Definition: object.h:228
virtual void updatedhdqRef(const fmatvec::Mat &dhdqRef, int i=0)
references to object Jacobian for implicit integration of dynamic system parent regarding positions ...
Definition: object.cc:200
fmatvec::Mat T
linear relation matrix of differentiated position and velocity parameters
Definition: object.h:339
int qSize
size of object positions
Definition: object.h:285
virtual void initz()
Definition: object.cc:266
fmatvec::Vec q
positions, velocities
Definition: object.h:310
basic class of MBSim mainly for plotting
Definition: element.h:58
virtual void updateudRef(const fmatvec::Vec &udRef, int i=0)
references to differentiated velocities of dynamic system parent
Definition: object.cc:180
virtual void writez(H5::GroupBase *group)
writes its z-Vector to a subgroup of the given group
Definition: object.cc:287
fmatvec::SymMat M[2]
mass matrix
Definition: object.h:344
int hSize[2]
size of object h-vector (columns of J)
Definition: object.h:295
virtual void updateLLMRef(const fmatvec::SymMat &ref, int i=0)
references to Cholesky decomposition of dynamic system parent
Definition: object.cc:224
virtual void updatedhduRef(const fmatvec::SqrMat &dhduRef, int i=0)
references to object Jacobian for implicit integration of dynamic system parent regarding velocities ...
Definition: object.cc:204
InitStage
The stages of the initialization.
Definition: element.h:97
virtual void updatedhdtRef(const fmatvec::Vec &dhdtRef, int i=0)
references to object Jacobian for implicit integration of dynamic system parent regarding time ...
Definition: object.cc:208
PlotFeatureStatus getPlotFeature(PlotFeature pf)
Definition: element.h:217
fmatvec::SymMat LLM[2]
LU-decomposition of mass matrix.
Definition: object.h:349
std::string name
name of element
Definition: element.h:290
virtual void updateuRef(const fmatvec::Vec &uRef)
references to velocities of dynamic system parent
Definition: object.cc:172
fmatvec::Mat dhdq
Jacobians of h.
Definition: object.h:332
virtual void init(InitStage stage)
initialize object at start of simulation with respect to contours and frames
Definition: object.cc:228
virtual void updateqRef(const fmatvec::Vec &qRef)
references to positions of dynamic system parent
Definition: object.cc:164
virtual void closePlot()
closes plot file
Definition: object.cc:154
int uSize[2]
size of object velocities
Definition: object.h:290
virtual void updatehRef(const fmatvec::Vec &hRef, int i=0)
references to smooth force vector of dynamic system parent
Definition: object.cc:196
virtual void readz0(H5::GroupBase *group)
reads the z-Vector of a subgroup of the given group
Definition: object.cc:293
fmatvec::Vec h[2]
complete and object smooth and nonsmooth right hand side
Definition: object.h:325
H5::GroupBase * plotGroup
associated plot group
Definition: element.h:324
virtual ~Object()
destructor
Definition: object.cc:49
virtual void calcSize(int j)
calculates size of right hand side
Definition: object.h:223
virtual double computePotentialEnergy()
Definition: object.h:233
virtual void plot(double t, double dt=1)
plots time dependent data
Definition: object.cc:122
PlotFeatureStatus
Plot feature status.
Definition: element.h:61
virtual std::string getType() const
Definition: object.h:90
virtual void updateuallRef(const fmatvec::Vec &uallRef)
references to velocities of dynamic system parent
Definition: object.cc:176
PlotFeature
Plot Features.
Definition: element.h:74
virtual void updaterRef(const fmatvec::Vec &ref, int i=0)
references to nonsmooth force vector of dynamic system parent
Definition: object.cc:212
virtual void updateudallRef(const fmatvec::Vec &udallRef, int i=0)
references to differentiated velocities of dynamic system parent
Definition: object.cc:184
virtual void updateMRef(const fmatvec::SymMat &ref, int i=0)
references to mass matrix of dynamic system parent
Definition: object.cc:220
virtual void updateqdRef(const fmatvec::Vec &qdRef)
references to differentiated positions of dynamic system parent
Definition: object.cc:168
fmatvec::Vec qd
differentiated positions, velocities
Definition: object.h:320
PlotFeatureStatus getPlotFeatureForChildren(PlotFeature pf)
Definition: element.h:222
virtual void facLLM(int i=0)
perform Cholesky decomposition of mass martix
Definition: object.h:217
class for all objects having own dynamics and mass
Definition: object.h:42
fmatvec::Vec q0
initial position, velocity
Definition: object.h:315
virtual void setUpInverseKinetics()
TODO.
Definition: object.h:238

Impressum / Disclaimer / Datenschutz Generated by doxygen 1.8.5 Valid HTML