All Classes Namespaces Functions Variables Typedefs Enumerations Enumerator Pages
dynamic_system_solver.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  * rzander@users.berlios.de
19  */
20 
21 #ifndef _DYNAMIC_SYSTEM_SOLVER_H_
22 #define _DYNAMIC_SYSTEM_SOLVER_H_
23 
24 #include "mbsim/group.h"
25 #include "fmatvec/sparse_matrix.h"
26 
27 namespace MBSim {
28 
29  class Graph;
30 
49  class DynamicSystemSolver : public Group {
50  public:
51 
55  enum Solver { FixedPointTotal, FixedPointSingle, GaussSeidel, LinearEquations, RootFinding };
56 
60  enum Strategy { global, local };
61 
65  enum LinAlg { LUDecomposition, LevenbergMarquardt, PseudoInverse };
66 
71  DynamicSystemSolver(const std::string &name="");
72 
76  virtual ~DynamicSystemSolver();
77 
116  void initialize();
117 
118  /*
119  * If true (the default) then
120  * the simulation output files (h5 files) are deleted/truncated/regenerated.
121  * If false, then these files left are untouched from a previous run.
122  * This is useful to regenerate the
123  * e.g. OpenMBV XML files without doing a new integration.
124  */
125  void setTruncateSimulationFiles(bool trunc) { truncateSimulationFiles=trunc; }
126 
127  /* INHERITED INTERFACE OF GROUP */
128  void init(InitStage stage);
129  using Group::plot;
130  /***************************************************/
131 
132  /* INHERITED INTERFACE OF DYNAMICSYSTEM */
133  virtual int solveConstraintsFixpointSingle();
134  virtual int solveImpactsFixpointSingle(double dt = 0);
135  virtual int solveConstraintsGaussSeidel();
136  virtual int solveImpactsGaussSeidel(double dt = 0);
137  virtual int solveConstraintsRootFinding();
138  virtual int solveImpactsRootFinding(double dt = 0);
139  virtual void checkConstraintsForTermination();
140  virtual void checkImpactsForTermination(double dt = 0);
141  /***************************************************/
142 
143  /* INHERITED INTERFACE OF OBJECTINTERFACE */
144  virtual void updateh(double t, int i=0);
145  virtual void updateM(double t, int i=0);
146  virtual void updateStateDependentVariables(double t); // this function is called once every time step by every integrator
147  /***************************************************/
148 
149  /* INHERITED INTERFACE OF LINKINTERFACE */
150 
161  virtual void updater(double t, int j=0);
162  virtual void updatewb(double t, int j=0);
163  virtual void updateW(double t, int j=0);
164  virtual void updateV(double t, int j=0);
165  /***************************************************/
166 
167  /* INHERITED INTERFACE OF ELEMENT */
169  virtual std::string getType() const { return "DynamicSystemSolver"; }
170  virtual void plot(const fmatvec::Vec& z, double t, double dt=1); // TODO completely rearrange
171  virtual void plot2(const fmatvec::Vec& z, double t, double dt=1); // TODO completely rearrange
172 
173  virtual void closePlot();
174  /***************************************************/
175 
176  /* INHERITED INTERFACE FOR DERIVED CLASS */
181  virtual int solveConstraints();
182 
188  virtual int solveImpacts(double dt = 0);
189  /***************************************************/
190 
191  /* GETTER / SETTER */
192 
193  void setConstraintSolver(Solver solver_) { contactSolver = solver_; }
194  void setImpactSolver(Solver solver_) { impactSolver = solver_; }
195  const Solver& getConstraintSolver() { return contactSolver; }
196  const Solver& getImpactSolver() { return impactSolver; }
197  void setTermination(bool term_) { term = term_; }
198  void setStrategy(Strategy strategy_) { strategy = strategy_; }
199  void setMaxIter(int iter) { maxIter = iter; }
200  void setHighIter(int iter) { highIter = iter; }
201  void setNumJacProj(bool numJac_) { numJac = numJac_; }
202  void setMaxDampingSteps(int maxDSteps) { maxDampingSteps = maxDSteps; }
203  void setLevenbergMarquardtParam(double lmParm_) { lmParm = lmParm_; }
204  void setLinAlg(LinAlg linAlg_) { linAlg = linAlg_; }
205 
206  void setUseOldla(bool flag) { useOldla = flag; }
207  void setDecreaseLevels(const fmatvec::VecInt &decreaseLevels_) { decreaseLevels = decreaseLevels_; }
208  void setCheckTermLevels(const fmatvec::VecInt &checkTermLevels_) { checkTermLevels = checkTermLevels_; }
209  void setCheckGSize(bool checkGSize_) { checkGSize = checkGSize_; }
210  void setLimitGSize(int limitGSize_) { limitGSize = limitGSize_; checkGSize = false; }
211 
212  const fmatvec::SparseMat& getGs() const { return Gs; }
213  fmatvec::SparseMat& getGs() { return Gs; }
214  const fmatvec::SqrMat& getG() const { return G; }
215  fmatvec::SqrMat& getG() { return G; }
216  const fmatvec::Vec& getb() const { return b; }
217  fmatvec::Vec& getb() { return b; }
218  const fmatvec::SqrMat& getJprox() const { return Jprox; }
219  fmatvec::SqrMat& getJprox() { return Jprox; }
220 
221  const fmatvec::Mat& getWParent(int i=0) const { return WParent[i]; }
222  const fmatvec::Mat& getVParent(int i=0) const { return VParent[i]; }
223  const fmatvec::Vec& getlaParent() const { return laParent; }
224  const fmatvec::Vec& getgdParent() const { return gdParent; }
225  const fmatvec::Vec& getresParent() const { return resParent; }
226  const fmatvec::Vec& getrFactorParent() const { return rFactorParent; }
227 
228  DynamicSystemSolver* getDynamicSystemSolver() { return this; }
229  bool getIntegratorExitRequest() { return integratorExitRequest; }
230  int getMaxIter() {return maxIter;}
231  /***************************************************/
232 
237 
244  fmatvec::Vec deltau(const fmatvec::Vec &zParent, double t, double dt);
245 
252  fmatvec::Vec deltaq(const fmatvec::Vec &zParent, double t, double dt);
253 
260  fmatvec::Vec deltax(const fmatvec::Vec &zParent, double t, double dt);
261 
265  virtual void initz(fmatvec::Vec& z0);
266 
272 
278  int (DynamicSystemSolver::*solveImpacts_)(double dt);
279 
286 
293  int solveImpactsLinearEquations(double dt = 0);
294 
299  virtual void updateG(double t, int i=0);
300 
304  void decreaserFactors();
305 
313  void update(const fmatvec::Vec &z, double t, int options=0);
314 
321  virtual void shift(fmatvec::Vec& z, const fmatvec::VecInt& jsv, double t);
322 
329  virtual void zdot(const fmatvec::Vec& z, fmatvec::Vec& zd, double t);
330 
336  virtual fmatvec::Vec zdot(const fmatvec::Vec &zParent, double t);
337 
344  virtual void getsv(const fmatvec::Vec& z, fmatvec::Vec& svExt, double t);
345 
350  void getLinkStatus(fmatvec::VecInt &LinkStatusExt, double t);
351 
357  void getLinkStatusReg(fmatvec::VecInt &LinkStatusRegExt, double t);
362  void projectGeneralizedPositions(double t, int mode, bool fullUpdate=false);
363 
368  void projectGeneralizedVelocities(double t, int mode);
369 
374  void savela(double dt=1.0);
375 
380  void initla(double dt=1.0);
381 
385  double computeKineticEnergy() { return 0.5*u.T()*M[0]*u; }
386 
391  double computePotentialEnergy();
392 
397  void addElement(Element *element_);
398 
404  Element* getElement(const std::string &name);
405 
409  std::string getSolverInfo();
410 
415  void setStopIfNoConvergence(bool flag, bool dropInfo = false) { stopIfNoConvergence = flag; dropContactInfo=dropInfo; }
416 
420  void dropContactMatrices();
421 
425  static void sigInterruptHandler(int);
426 
430  static void sigAbortHandler(int);
431 
432  static void sigSegfaultHandler(int);
433 
434  // TODO just for testing
435  void setPartialEventDrivenSolver(bool peds_) { peds = peds_; }
436 
442  void writez(std::string fileName, bool formatH5=true);
443 
448  void readz0(std::string fileName);
449 
450  virtual void initializeUsingXML(xercesc::DOMElement *element);
451  virtual xercesc::DOMElement* writeXMLFile(xercesc::DOMNode *element);
452 
453  static DynamicSystemSolver* readXMLFile(const std::string &filename);
454  void writeXMLFile(const std::string &name);
455  void writeXMLFile() { writeXMLFile(getName()); }
456 
457 
462  void setProjectionTolerance(double tol) { tolProj = tol; }
463 
467  void setgTol(double tol) {gTol = tol; Group::setgTol(tol);}
468 
472  void setgdTol(double tol) {gdTol = tol; Group::setgdTol(tol);}
473 
477  void setgddTol(double tol) {gddTol = tol; Group::setgddTol(tol);}
478 
482  void setlaTol(double tol) {laTol = tol; Group::setlaTol(tol);}
483 
487  void setLaTol(double tol) {LaTol = tol; Group::setLaTol(tol);}
488 
493  void updatezRef(const fmatvec::Vec &ext);
494 
499  void setFlushEvery(unsigned int every) {flushEvery = every;}
500 
501  void setAlwaysConsiderContact(bool alwaysConsiderContact_) {alwaysConsiderContact = alwaysConsiderContact_;}
502 
503  void setInverseKinetics(bool inverseKinetics_) {inverseKinetics = inverseKinetics_;}
504 
505  void setInitialProjection(bool initialProjection_) {initialProjection = initialProjection_;}
506 
507  void setUseConstraintSolverForPlot(bool useConstraintSolverForPlot_) {useConstraintSolverForPlot = useConstraintSolverForPlot_;}
508 
509  fmatvec::Mat dhdq(double t, int lb=0, int ub=0);
510  fmatvec::Mat dhdu(double t, int lb=0, int ub=0);
511  fmatvec::Mat dhdx(double t);
512  fmatvec::Vec dhdt(double t);
513 
514  void setRootID(int ID) {rootID = ID;}
515  int getRootID() const {return rootID;}
516  void setq(const fmatvec::Vec& q_){ q = q_;}
517  void setLa(const fmatvec::Vec& la_){ la = la_;}
518 
519  protected:
524 
529 
534 
539 
544 
549 
554 
559 
564 
569 
574 
579 
584 
589 
590  fmatvec::Vec udParent1;
591 
596 
601 
606 
611 
615  fmatvec::VecInt jsvParent;
616 
620  fmatvec::VecInt LinkStatusParent;
621 
626  fmatvec::VecInt LinkStatusRegParent;
627 
631  fmatvec::SparseMat Gs;
632 
637 
642 
647 
651  bool term;
652 
656  int maxIter, highIter, maxDampingSteps;
657 
661  double lmParm;
662 
666  Solver contactSolver, impactSolver;
667 
672 
677 
682 
687 
691  bool useOldla;
692 
696  bool numJac;
697 
701  fmatvec::VecInt decreaseLevels;
702 
706  fmatvec::VecInt checkTermLevels;
707 
712 
717 
722 
726  bool peds;
727 
731  unsigned int driftCount;
732 
737  unsigned int flushEvery;
738 
743  unsigned int flushCount;
744 
748  double tolProj;
749 
754  void updatezdRef(const fmatvec::Vec &ext);
755 
760  virtual void updaterFactors();
761 
766  void computeConstraintForces(double t);
767 
771  fmatvec::Vec laInverseKineticsParent;
772  fmatvec::Mat bInverseKineticsParent;
773 
777  fmatvec::Mat WInverseKineticsParent[2];
778 
779  bool alwaysConsiderContact;
780  bool inverseKinetics;
781  bool initialProjection;
782  bool useConstraintSolverForPlot;
783 
784  fmatvec::Vec corrParent;
785 
786  int rootID;
787 
788  double gTol, gdTol, gddTol, laTol, LaTol;
789 
790  private:
794  void constructor();
795 
800 
804  static bool exitRequest;
805 
809  bool READZ0;
810 
819  void addToGraph(Graph* graph, fmatvec::SqrMat &A, int i, std::vector<Object*> &objList);
820 
821  bool truncateSimulationFiles;
822 
823  // Holds the dynamic systems before the "reorganize hierarchie" takes place.
824  // This is required since all elements of all other containers from DynamicSystem are readded to DynamicSystemSolver,
825  // except this container (which is a "pure" container = no calculation is done in DynamicSystem)
826  // However, we must hold this container until the dtor of DynamicSystemSolver is called to avoid the deallocation of other
827  // elements hold by DynamicSystem elements (especially (currently only) hdf5File)
828  std::vector<DynamicSystem*> dynamicsystemPreReorganize;
829 
830  double facSizeGs;
831  };
832 
833 }
834 
835 #endif /* _DYNAMIC_SYSTEM_SOLVER_H_ */
836 
bool useOldla
flag if contac force parameter of last time step should be used
Definition: dynamic_system_solver.h:691
fmatvec::Vec fParent
right hand side of order one parameters
Definition: dynamic_system_solver.h:605
int solveConstraintsLinearEquations()
solution of contact equations with Cholesky decomposition
Definition: dynamic_system_solver.cc:978
solver-interface for dynamic systems
Definition: solver.h:38
fmatvec::Vec deltau(const fmatvec::Vec &zParent, double t, double dt)
Definition: dynamic_system_solver.cc:935
virtual int solveConstraintsRootFinding()
solve contact equations with Newton scheme
Definition: dynamic_system_solver.cc:533
virtual void setlaTol(double tol)
Definition: dynamic_system.cc:1292
fmatvec::Vec q
positions, differentiated positions, initial positions
Definition: dynamic_system.h:784
fmatvec::Vec gParent
relative distances
Definition: dynamic_system_solver.h:573
void updatezRef(const fmatvec::Vec &ext)
references to external state
Definition: dynamic_system_solver.cc:1332
bool checkGSize
boolean if force action matrix should be resized in each step
Definition: dynamic_system_solver.h:711
virtual void updateG(double t, int i=0)
updates mass action matrix
Definition: dynamic_system_solver.cc:989
fmatvec::Vec rParent[2]
nonsmooth right hand side
Definition: dynamic_system_solver.h:600
fmatvec::SymMat M[2]
mass matrix
Definition: dynamic_system.h:774
void getLinkStatus(fmatvec::VecInt &LinkStatusExt, double t)
Definition: dynamic_system_solver.cc:1049
bool term
boolean to check for termination of contact equations solution
Definition: dynamic_system_solver.h:651
void readz0(std::string fileName)
reads state from a file
Definition: dynamic_system_solver.cc:1324
bool dropContactInfo
flag if contact matrices for debugging should be dropped in no-convergence case
Definition: dynamic_system_solver.h:686
double tolProj
Tolerance for projection of generalized position.
Definition: dynamic_system_solver.h:748
solver interface for modelling and simulation of dynamic systems
Definition: dynamic_system_solver.h:49
void savela(double dt=1.0)
save contact force parameter for use as starting value in next time step
Definition: dynamic_system_solver.cc:1165
fmatvec::SparseMat Gs
sparse mass action matrix
Definition: dynamic_system_solver.h:631
virtual void zdot(const fmatvec::Vec &z, fmatvec::Vec &zd, double t)
update for event driven integrator during smooth phase
Definition: dynamic_system_solver.cc:1042
Strategy strategy
relaxarion strategy for solution of fixed-point scheme
Definition: dynamic_system_solver.h:671
void setgddTol(double tol)
Definition: dynamic_system_solver.h:477
fmatvec::VecInt checkTermLevels
TODO.
Definition: dynamic_system_solver.h:706
DynamicSystemSolver(const std::string &name="")
constructor
Definition: dynamic_system_solver.cc:69
Solver contactSolver
solver for contact equations and impact equations
Definition: dynamic_system_solver.h:666
void setgdTol(double tol)
Definition: dynamic_system_solver.h:472
void projectGeneralizedPositions(double t, int mode, bool fullUpdate=false)
drift projection for positions
Definition: dynamic_system_solver.cc:1065
bool stopIfNoConvergence
flag if the contact equations should be stopped if there is no convergence
Definition: dynamic_system_solver.h:681
fmatvec::Vec deltax(const fmatvec::Vec &zParent, double t, double dt)
return x-state difference for current time
Definition: dynamic_system_solver.cc:952
void setgTol(double tol)
Definition: dynamic_system_solver.h:467
fmatvec::Mat TParent
matrix of linear relation between differentiated positions and velocities
Definition: dynamic_system_solver.h:533
fmatvec::Vec sParent
TODO.
Definition: dynamic_system_solver.h:563
fmatvec::Vec u
velocities, differentiated velocities, initial velocities
Definition: dynamic_system.h:789
double computePotentialEnergy()
compute potential energy of entire dynamic system change? TODO
Definition: dynamic_system_solver.cc:1175
fmatvec::VecInt LinkStatusRegParent
status vector of single valued links
Definition: dynamic_system_solver.h:626
fmatvec::Vec gdParent
relative velocities
Definition: dynamic_system_solver.h:578
virtual int solveConstraints()
solves prox-functions for contact forces using sparsity structure
Definition: dynamic_system_solver.cc:846
fmatvec::VecInt LinkStatusParent
status vector of set valued links with piecewise link equation (which piece is valid) ...
Definition: dynamic_system_solver.h:620
group ingredients do not depend on each other
Definition: group.h:35
LinAlg
linear algebra for Newton scheme in solution of contact equations
Definition: dynamic_system_solver.h:65
void projectGeneralizedVelocities(double t, int mode)
drift projection for positions
Definition: dynamic_system_solver.cc:1124
virtual void closePlot()
closes plot file
Definition: dynamic_system_solver.cc:840
unsigned int flushCount
counts plot-calls until files to be flushed TODO
Definition: dynamic_system_solver.h:743
std::string getSolverInfo()
Definition: dynamic_system_solver.cc:1227
virtual void setgTol(double tol)
Definition: dynamic_system.cc:1271
basic class of MBSim mainly for plotting
Definition: element.h:58
static void sigInterruptHandler(int)
handler for user interrupt signal
Definition: dynamic_system_solver.cc:1285
double lmParm
Levenberg-Marquard parameter.
Definition: dynamic_system_solver.h:661
void addElement(Element *element_)
Definition: dynamic_system_solver.cc:1188
void addToGraph(Graph *graph, fmatvec::SqrMat &A, int i, std::vector< Object * > &objList)
adds list of objects to tree
Definition: dynamic_system_solver.cc:1542
int solveImpactsLinearEquations(double dt=0)
solution of contact equations with Cholesky decomposition on velocity level
Definition: dynamic_system_solver.cc:984
virtual void setLaTol(double tol)
Definition: dynamic_system.cc:1299
fmatvec::Vec rFactorParent
relaxation parameters for contact equations
Definition: dynamic_system_solver.h:558
virtual int solveImpactsGaussSeidel(double dt=0)
solve impact equations with Gauss-Seidel scheme on velocity level
Definition: dynamic_system_solver.cc:513
virtual void plot(double t, double dt)
plots time dependent data
Definition: dynamic_system.cc:368
fmatvec::VecInt decreaseLevels
decreasing relaxation factors is done in levels containing the number of contact iterations as condit...
Definition: dynamic_system_solver.h:701
Strategy
relaxation strategies in solution of contact equations
Definition: dynamic_system_solver.h:60
void setStopIfNoConvergence(bool flag, bool dropInfo=false)
Definition: dynamic_system_solver.h:415
void getLinkStatusReg(fmatvec::VecInt &LinkStatusRegExt, double t)
Definition: dynamic_system_solver.cc:1057
fmatvec::SymMat MParent[2]
mass matrix
Definition: dynamic_system_solver.h:523
fmatvec::Vec b
TODO.
Definition: dynamic_system_solver.h:646
virtual void checkImpactsForTermination(double dt=0)
validate force laws concerning given tolerances on velocity level
Definition: dynamic_system_solver.cc:695
void dropContactMatrices()
Definition: dynamic_system_solver.cc:1266
void initialize()
Initialize the system.
Definition: dynamic_system_solver.cc:84
double computeKineticEnergy()
compute kinetic energy of entire dynamic system
Definition: dynamic_system_solver.h:385
bool numJac
flag if Jacobian for Newton scheme should be calculated numerically
Definition: dynamic_system_solver.h:696
virtual void updater(double t, int j=0)
update smooth link force law
Definition: dynamic_system_solver.cc:821
unsigned int flushEvery
flushes all hdf5-files every x-times the plot-routine is called TODO
Definition: dynamic_system_solver.h:737
unsigned int driftCount
TODO, additional stop in event driven solver for drift correction.
Definition: dynamic_system_solver.h:731
fmatvec::Mat WParent[2]
contact force directions
Definition: dynamic_system_solver.h:538
int warnLevel
level for warning output (0-2)
Definition: dynamic_system_solver.h:721
virtual void initz(fmatvec::Vec &z0)
initialises state variables
Definition: dynamic_system_solver.cc:960
class for tree-structured mechanical systems with recursive and flat memory mechanism ...
Definition: graph.h:37
InitStage
The stages of the initialization.
Definition: element.h:97
virtual int solveImpacts(double dt=0)
solves prox-functions for impacts on velocity level using sparsity structure
Definition: dynamic_system_solver.cc:880
int(DynamicSystemSolver::* solveImpacts_)(double dt)
function pointer for election of prox-solver for impact equations on velocity level ...
Definition: dynamic_system_solver.h:278
void setProjectionTolerance(double tol)
set tolerance for projection of generalized position
Definition: dynamic_system_solver.h:462
int(DynamicSystemSolver::* solveConstraints_)()
Definition: dynamic_system_solver.h:271
Element * getElement(const std::string &name)
Definition: dynamic_system_solver.cc:1200
fmatvec::VecInt jsvParent
boolean evaluation of stopvector
Definition: dynamic_system_solver.h:615
virtual void getsv(const fmatvec::Vec &z, fmatvec::Vec &svExt, double t)
evaluation of stop vector
Definition: dynamic_system_solver.cc:1694
fmatvec::SymMat LLMParent[2]
Cholesky decomposition of mass matrix.
Definition: dynamic_system_solver.h:528
fmatvec::Vec deltaq(const fmatvec::Vec &zParent, double t, double dt)
Definition: dynamic_system_solver.cc:944
LinAlg linAlg
linear system solver used for Newton scheme in contact equations
Definition: dynamic_system_solver.h:676
virtual void setgdTol(double tol)
Definition: dynamic_system.cc:1278
std::string name
name of element
Definition: element.h:290
int limitGSize
TODO.
Definition: dynamic_system_solver.h:716
virtual int solveConstraintsFixpointSingle()
solve contact equations with single step fixed point scheme
Definition: dynamic_system_solver.cc:429
virtual void updaterFactors()
update relaxation factors for contact equations
Definition: dynamic_system_solver.cc:1356
static bool exitRequest
boolean signal evaluation for end integration set by user
Definition: dynamic_system_solver.h:804
bool peds
TODO, flag for occuring impact and sticking in event driven solver.
Definition: dynamic_system_solver.h:726
const std::string & getName() const
Definition: element.h:161
fmatvec::Vec hParent[2]
smooth, smooth with respect to objects, smooth with respect to links right hand side ...
Definition: dynamic_system_solver.h:595
void constructor()
set plot feature default values
Definition: dynamic_system_solver.cc:1382
virtual int solveImpactsRootFinding(double dt=0)
solve impact equations with Newton scheme on velocity level
Definition: dynamic_system_solver.cc:609
fmatvec::SqrMat G
mass action matrix
Definition: dynamic_system_solver.h:641
void setlaTol(double tol)
Definition: dynamic_system_solver.h:482
bool integratorExitRequest
boolean signal evaluation for end integration set by program
Definition: dynamic_system_solver.h:799
fmatvec::Vec wbParent
TODO.
Definition: dynamic_system_solver.h:548
virtual void checkConstraintsForTermination()
validate force laws concerning given tolerances
Definition: dynamic_system_solver.cc:684
void init(InitStage stage)
plots time series header
Definition: dynamic_system_solver.cc:112
void update(const fmatvec::Vec &z, double t, int options=0)
update of dynamic system for time-stepping integrator
Definition: dynamic_system_solver.cc:1008
virtual void shift(fmatvec::Vec &z, const fmatvec::VecInt &jsv, double t)
update for event driven integrator for event
Definition: dynamic_system_solver.cc:1552
void initla(double dt=1.0)
load contact force parameter for use as starting value
Definition: dynamic_system_solver.cc:1170
fmatvec::Vec la
contact force parameters
Definition: dynamic_system.h:809
int maxIter
maximum number of contact iterations, high number of contact iterations for warnings, maximum number of damping steps for Newton scheme
Definition: dynamic_system_solver.h:656
Solver
solver for contact equations
Definition: dynamic_system_solver.h:55
fmatvec::Vec zParent
state
Definition: dynamic_system_solver.h:583
fmatvec::Vec laParent
contact force parameters
Definition: dynamic_system_solver.h:553
virtual std::string getType() const
Definition: dynamic_system_solver.h:169
fmatvec::Vec resParent
residuum of contact equations
Definition: dynamic_system_solver.h:568
void setLaTol(double tol)
Definition: dynamic_system_solver.h:487
void setFlushEvery(unsigned int every)
set the number of plot-routine-calls after which all hdf5-files will be flushed
Definition: dynamic_system_solver.h:499
virtual ~DynamicSystemSolver()
destructor
Definition: dynamic_system_solver.cc:74
fmatvec::Vec svParent
stopvector (rootfunctions for event driven integration
Definition: dynamic_system_solver.h:610
fmatvec::SqrMat Jprox
JACOBIAN of contact equations for Newton scheme.
Definition: dynamic_system_solver.h:636
fmatvec::VecInt jsv
boolean evaluation of stop vector concerning roots
Definition: dynamic_system.h:839
fmatvec::Mat VParent[2]
condensed contact force directions
Definition: dynamic_system_solver.h:543
void writez(std::string fileName, bool formatH5=true)
writes state to a file
Definition: dynamic_system_solver.cc:1304
virtual int solveConstraintsGaussSeidel()
solve contact equations with Gauss-Seidel scheme
Definition: dynamic_system_solver.cc:493
static void sigAbortHandler(int)
handler for abort signals
Definition: dynamic_system_solver.cc:1290
void decreaserFactors()
decrease relaxation factors if mass action matrix is not diagonal dominant
Definition: dynamic_system_solver.cc:1002
void computeInitialCondition()
compute initial condition for links for event driven integrator
Definition: dynamic_system_solver.cc:917
void updatezdRef(const fmatvec::Vec &ext)
references to differentiated external state
Definition: dynamic_system_solver.cc:1343
bool READZ0
is a state read from a file
Definition: dynamic_system_solver.h:809
void computeConstraintForces(double t)
compute inverse kinetics constraint forces
Definition: dynamic_system_solver.cc:1378
virtual void setgddTol(double tol)
Definition: dynamic_system.cc:1285
fmatvec::Vec zdParent
differentiated state
Definition: dynamic_system_solver.h:588
virtual int solveImpactsFixpointSingle(double dt=0)
solve impact equations with single step fixed point scheme on velocity level
Definition: dynamic_system_solver.cc:461

Impressum / Disclaimer / Datenschutz Generated by doxygen 1.8.5 Valid HTML