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  */
19 
20 #ifndef _DYNAMIC_SYSTEM_SOLVER_H_
21 #define _DYNAMIC_SYSTEM_SOLVER_H_
22 
23 #include "mbsim/group.h"
24 #include "fmatvec/sparse_matrix.h"
25 
26 namespace MBSim {
27 
28  class Graph;
29 
48  class DynamicSystemSolver : public Group {
49  public:
50 
54  enum Solver { FixedPointTotal, FixedPointSingle, GaussSeidel, LinearEquations, RootFinding };
55 
59  enum Strategy { global, local };
60 
64  enum LinAlg { LUDecomposition, LevenbergMarquardt, PseudoInverse };
65 
70  DynamicSystemSolver(const std::string &name="");
71 
75  virtual ~DynamicSystemSolver();
76 
115  void initialize();
116 
117  /*
118  * If true (the default) then
119  * the simulation output files (h5 files) are deleted/truncated/regenerated.
120  * If false, then these files left are untouched from a previous run.
121  * This is useful to regenerate the
122  * e.g. OpenMBV XML files without doing a new integration.
123  */
124  void setTruncateSimulationFiles(bool trunc) { truncateSimulationFiles=trunc; }
125 
126  /* INHERITED INTERFACE OF GROUP */
127  void init(InitStage stage);
128  using Group::plot;
129  /***************************************************/
130 
131  /* INHERITED INTERFACE OF DYNAMICSYSTEM */
132  virtual int solveConstraintsFixpointSingle();
133  virtual int solveImpactsFixpointSingle();
134  virtual int solveConstraintsGaussSeidel();
135  virtual int solveImpactsGaussSeidel();
136  virtual int solveConstraintsRootFinding();
137  virtual int solveImpactsRootFinding();
138  virtual void checkConstraintsForTermination();
139  virtual void checkImpactsForTermination();
140  /***************************************************/
141 
142  /* INHERITED INTERFACE OF OBJECTINTERFACE */
143  virtual void updateT();
144  virtual void updateh(int i=0);
145  virtual void updateM();
146  virtual void updateLLM();
147  virtual void updatezd();
148  /***************************************************/
149 
150  /* INHERITED INTERFACE OF LINKINTERFACE */
151 
162  virtual void updater(int j=0);
163  virtual void updaterdt();
164  virtual void updatewb();
165  virtual void updateg();
166  virtual void updategd();
167  virtual void updateW(int j=0);
168  virtual void updateV(int j=0);
169  virtual void updatebc();
170  virtual void updatebi();
171  virtual void updatela();
172  virtual void updateLa();
173  /***************************************************/
174 
175  /* INHERITED INTERFACE OF ELEMENT */
177  virtual std::string getType() const { return "DynamicSystemSolver"; }
178 
179  virtual void closePlot();
180  /***************************************************/
181 
182  /* GETTER / SETTER */
183 
184  void setConstraintSolver(Solver solver_) { contactSolver = solver_; }
185  void setImpactSolver(Solver solver_) { impactSolver = solver_; }
186  const Solver& getConstraintSolver() { return contactSolver; }
187  const Solver& getImpactSolver() { return impactSolver; }
188  void setTermination(bool term_) { term = term_; }
189  void setStrategy(Strategy strategy_) { strategy = strategy_; }
190  void setMaxIter(int iter) { maxIter = iter; }
191  void setHighIter(int iter) { highIter = iter; }
192  void setNumJacProj(bool numJac_) { numJac = numJac_; }
193  void setMaxDampingSteps(int maxDSteps) { maxDampingSteps = maxDSteps; }
194  void setLevenbergMarquardtParam(double lmParm_) { lmParm = lmParm_; }
195  void setLinAlg(LinAlg linAlg_) { linAlg = linAlg_; }
196 
197  void setUseOldla(bool flag) { useOldla = flag; }
198  void setDecreaseLevels(const fmatvec::VecInt &decreaseLevels_) { decreaseLevels = decreaseLevels_; }
199  void setCheckTermLevels(const fmatvec::VecInt &checkTermLevels_) { checkTermLevels = checkTermLevels_; }
200  void setCheckGSize(bool checkGSize_) { checkGSize = checkGSize_; }
201  void setLimitGSize(int limitGSize_) { limitGSize = limitGSize_; checkGSize = false; }
202 
203  double& getTime() { return t; }
204  const double& getTime() const { return t; }
205  void setTime(double t_) { t = t_; }
206 
207  double getStepSize() const { return dt; }
208  void setStepSize(double dt_) { dt = dt_; }
209 
210  fmatvec::Vec& getState() { return zParent; }
211  const fmatvec::Vec& getState() const { return zParent; }
212  void setState(const fmatvec::Vec &z) { zParent = z; }
213 
214  void setzd(const fmatvec::Vec &zd) { zdParent = zd; }
215 
216  const fmatvec::SqrMat& getG(bool check=true) const { assert((not check) or (not updG)); return G; }
217  const fmatvec::SparseMat& getGs(bool check=true) const { assert((not check) or (not updG)); return Gs; }
218  const fmatvec::Vec& getbc(bool check=true) const { assert((not check) or (not updbc)); return bc; }
219  const fmatvec::Vec& getbi(bool check=true) const { assert((not check) or (not updbi)); return bi; }
220  const fmatvec::SqrMat& getJprox() const { return Jprox; }
221  fmatvec::SqrMat& getG(bool check=true) { assert((not check) or (not updG)); return G; }
222  fmatvec::SparseMat& getGs(bool check=true) { assert((not check) or (not updG)); return Gs; }
223  fmatvec::Vec& getbc(bool check=true) { assert((not check) or (not updbc)); return bc; }
224  fmatvec::Vec& getbi(bool check=true) { assert((not check) or (not updbi)); return bi; }
225  fmatvec::SqrMat& getJprox() { return Jprox; }
226 
227  const fmatvec::Vec& evaldq() { if(upddq) updatedq(); return dq; }
228  const fmatvec::Vec& evaldu() { if(upddu) updatedu(); return du; }
229  const fmatvec::Vec& evaldx() { if(upddx) updatedx(); return dx; }
230  const fmatvec::Vec& evalzd();
231  const fmatvec::SqrMat& evalG() { if(updG) updateG(); return G; }
232  const fmatvec::SparseMat& evalGs() { if(updG) updateG(); return Gs; }
233  const fmatvec::Vec& evalbc() { if(updbc) updatebc(); return bc; }
234  const fmatvec::Vec& evalbi() { if(updbi) updatebi(); return bi; }
235  const fmatvec::Vec& evalsv();
236  const fmatvec::Vec& evalz0();
237  const fmatvec::Vec& evalla() { if(updla) updatela(); return la; }
238  const fmatvec::Vec& evalLa() { if(updLa) updateLa(); return La; }
239 
240  const fmatvec::Mat& getWParent(int i=0) const { return WParent[i]; }
241  const fmatvec::Mat& getVParent(int i=0) const { return VParent[i]; }
242  const fmatvec::Vec& getlaParent() const { return laParent; }
243  const fmatvec::Vec& getLaParent() const { return LaParent; }
244  const fmatvec::Vec& getgdParent() const { return gdParent; }
245  const fmatvec::Vec& getresParent() const { return resParent; }
246  const fmatvec::Vec& getrFactorParent() const { return rFactorParent; }
247 
248  DynamicSystemSolver* getDynamicSystemSolver() { return this; }
249  bool getIntegratorExitRequest() { return integratorExitRequest; }
250  int getMaxIter() {return maxIter;}
251  int getIterC() {return iterc;}
252  int getIterI() {return iteri;}
253  /***************************************************/
254 
259 
265 
272 
279 
287 
292  virtual void updateG();
293 
297  void decreaserFactors();
298 
305  virtual const fmatvec::Vec& shift();
306 
311  void getLinkStatus(fmatvec::VecInt &LinkStatusExt);
312 
317  void getLinkStatusReg(fmatvec::VecInt &LinkStatusRegExt);
318 
322  bool positionDriftCompensationNeeded(double gmax);
323 
327  bool velocityDriftCompensationNeeded(double gdmax);
328 
332  void projectGeneralizedPositions(int mode, bool fullUpdate=false);
333 
337  void projectGeneralizedVelocities(int mode);
338 
343  void savela();
344 
349  void initla();
350 
355  void saveLa();
356 
361  void initLa();
362 
366  double evalKineticEnergy() { return 0.5*u.T()*M*u; }
367 
372  double evalPotentialEnergy();
373 
378  void addElement(Element *element_);
379 
385  Element* getElement(const std::string &name);
386 
390  std::string getSolverInfo();
391 
396  void setStopIfNoConvergence(bool flag, bool dropInfo = false) { stopIfNoConvergence = flag; dropContactInfo=dropInfo; }
397 
401  void dropContactMatrices();
402 
403  // install the MBSim signal handler
404  static void installSignalHandler();
405 
409  static void sigInterruptHandler(int);
410 
414  static void sigAbortHandler(int);
415 
416  static void sigSegfaultHandler(int);
417 
418  void checkExitRequest();
419 
420  // TODO just for testing
421  void setPartialEventDrivenSolver(bool peds_) { peds = peds_; }
422 
428  void writez(std::string fileName, bool formatH5=true);
429 
434  void readz0(std::string fileName);
435 
436  virtual void initializeUsingXML(xercesc::DOMElement *element);
437 
438  static DynamicSystemSolver* readXMLFile(const std::string &filename);
439 
444  void setProjectionTolerance(double tol) { tolProj = tol; }
445 
449  void setgTol(double tol) {gTol = tol; Group::setgTol(tol);}
450 
454  void setgdTol(double tol) {gdTol = tol; Group::setgdTol(tol);}
455 
459  void setgddTol(double tol) {gddTol = tol; Group::setgddTol(tol);}
460 
464  void setlaTol(double tol) {laTol = tol; Group::setlaTol(tol);}
465 
469  void setLaTol(double tol) {LaTol = tol; Group::setLaTol(tol);}
470 
475  void updatezRef(const fmatvec::Vec &ext);
476 
481  void updatezdRef(const fmatvec::Vec &ext);
482 
487  void setFlushEvery(unsigned int every) {flushEvery = every;}
488 
489  void setAlwaysConsiderContact(bool alwaysConsiderContact_) {alwaysConsiderContact = alwaysConsiderContact_;}
490 
491  void setInverseKinetics(bool inverseKinetics_) {inverseKinetics = inverseKinetics_;}
492 
493  void setInitialProjection(bool initialProjection_) {initialProjection = initialProjection_;}
494 
495  void setUseConstraintSolverForPlot(bool useConstraintSolverForPlot_) {useConstraintSolverForPlot = useConstraintSolverForPlot_;}
496  bool getUseConstraintSolverForPlot() const { return useConstraintSolverForPlot; }
497 
498  fmatvec::Mat dhdq(int lb=0, int ub=0);
499  fmatvec::Mat dhdu(int lb=0, int ub=0);
500  fmatvec::Mat dhdx();
501  fmatvec::Vec dhdt();
502 
503  void setRootID(int ID) {rootID = ID;}
504  int getRootID() const {return rootID;}
505 
506  void resetUpToDate();
507 
508  bool getUpdateT() { return updT; }
509  bool getUpdateM() { return updM; }
510  bool getUpdateLLM() { return updLLM; }
511  bool getUpdateh(int j) { return updh[j]; }
512  bool getUpdater(int j) { return updr[j]; }
513  bool getUpdaterdt() { return updrdt; }
514  bool getUpdateW(int j) { return updW[j]; }
515  bool getUpdateV(int j) { return updV[j]; }
516  bool getUpdatewb() { return updwb; }
517  bool getUpdateg() { return updg; }
518  bool getUpdategd() { return updgd; }
519  bool getUpdatela() { return updla; }
520  bool getUpdateLa() { return updLa; }
521  bool getUpdatezd() { return updzd; }
522  bool getUpdatedq() { return upddq; }
523  bool getUpdatedu() { return upddu; }
524  bool getUpdatedx() { return upddx; }
525  void setUpdatela(bool updla_) { updla = updla_; }
526  void setUpdateLa(bool updLa_) { updLa = updLa_; }
527  void setUpdatebi(bool updbi_) { updbi = updbi_; }
528  void setUpdatebc(bool updbc_) { updbc = updbc_; }
529  void setUpdatezd(bool updzd_) { updzd = updzd_; }
530 
531  void resize_();
532 
537  void updategRef(const fmatvec::Vec &ref) { Group::updategRef(ref); updg = true; }
538 
543  void updategdRef(const fmatvec::Vec &ref) { Group::updategdRef(ref); updgd = true; }
544 
549  void updatewbRef(const fmatvec::Vec &ref) { Group::updatewbRef(ref); updwb = true; }
550 
556  void updateWRef(const fmatvec::Mat &ref, int i=0) { Group::updateWRef(ref,i); updW[i] = true; }
557 
563  void updateVRef(const fmatvec::Mat &ref, int i=0) { Group::updateVRef(ref,i); updV[i] = true; }
564 
568  virtual void updatelaInverseKinetics();
569 
570  virtual void updatedq();
571  virtual void updatedu();
572  virtual void updatedx();
573 
574  virtual void updateStopVector();
575 
576  void plot();
577 
578  protected:
582  double t;
583 
587  double dt;
588 
593 
598 
603 
608 
613 
618 
623 
628 
633 
638 
643 
648 
653 
658 
659  fmatvec::Vec dxParent, dqParent, duParent;
660 
665 
669  fmatvec::Vec rParent[2], rdtParent;
670 
675 
679  fmatvec::VecInt jsvParent;
680 
684  fmatvec::VecInt LinkStatusParent;
685 
690  fmatvec::VecInt LinkStatusRegParent;
691 
695  fmatvec::SparseMat Gs;
696 
701 
706 
711 
715  bool term;
716 
720  int maxIter, highIter, maxDampingSteps, iterc, iteri;
721 
725  double lmParm;
726 
730  Solver contactSolver, impactSolver;
731 
736 
741 
746 
751 
755  bool useOldla;
756 
760  bool numJac;
761 
765  fmatvec::VecInt decreaseLevels;
766 
770  fmatvec::VecInt checkTermLevels;
771 
776 
781 
786 
790  bool peds;
791 
796  unsigned int flushEvery;
797 
802  unsigned int flushCount;
803 
807  double tolProj;
808 
813  virtual void updaterFactors();
814 
818  fmatvec::Vec laInverseKineticsParent;
819  fmatvec::Mat bInverseKineticsParent;
820 
824  fmatvec::Mat WInverseKineticsParent;
825 
826  bool alwaysConsiderContact;
827  bool inverseKinetics;
828  bool initialProjection;
829  bool useConstraintSolverForPlot;
830 
831  fmatvec::Vec corrParent;
832 
833  int rootID;
834 
835  double gTol, gdTol, gddTol, laTol, LaTol;
836 
837  bool updT, updh[2], updr[2], updrdt, updM, updLLM, updW[2], updV[2], updwb, updg, updgd, updG, updbc, updbi, updsv, updzd, updla, updLa, upddq, upddu, upddx;
838 
839  bool solveDirectly;
840 
841  private:
845  void constructor();
846 
851 
855  static bool exitRequest;
856 
860  bool READZ0;
861 
870  void addToGraph(Graph* graph, fmatvec::SqrMat &A, int i, std::vector<Element*> &objList);
871 
872  bool truncateSimulationFiles;
873 
874  // Holds the dynamic systems before the "reorganize hierarchy" takes place.
875  // This is required since all elements of all other containers from DynamicSystem are readded to DynamicSystemSolver,
876  // except this container (which is a "pure" container = no calculation is done in DynamicSystem)
877  // However, we must hold this container until the dtor of DynamicSystemSolver is called to avoid the deallocation of other
878  // elements hold by DynamicSystem elements (especially (currently only) hdf5File)
879  std::vector<DynamicSystem*> dynamicsystemPreReorganize;
880 
881  double facSizeGs;
882  };
883 
884 }
885 
886 #endif /* _DYNAMIC_SYSTEM_SOLVER_H_ */
bool useOldla
flag if contac force parameter of last time step should be used
Definition: dynamic_system_solver.h:755
void updateVRef(const fmatvec::Mat &ref, int i=0)
references to condensed contact force direction matrix of dynamic system parent
Definition: dynamic_system_solver.h:563
void updatewbRef(const fmatvec::Vec &ref)
references to TODO of dynamic system parent
Definition: dynamic_system_solver.h:549
int solveConstraintsLinearEquations()
solution of contact equations with Cholesky decomposition
Definition: dynamic_system_solver.cc:885
solver-interface for dynamic systems
Definition: solver.h:38
virtual int solveConstraintsRootFinding()
solve contact equations with Newton scheme
Definition: dynamic_system_solver.cc:544
virtual void setlaTol(double tol)
Definition: dynamic_system.cc:1302
fmatvec::Vec gParent
relative distances
Definition: dynamic_system_solver.h:642
void updatezRef(const fmatvec::Vec &ext)
references to external state
Definition: dynamic_system_solver.cc:1347
double t
time
Definition: dynamic_system_solver.h:582
bool checkGSize
boolean if force action matrix should be resized in each step
Definition: dynamic_system_solver.h:775
int solveImpactsLinearEquations()
solution of contact equations with Cholesky decomposition on velocity level
Definition: dynamic_system_solver.cc:890
fmatvec::Vec rParent[2]
nonsmooth right hand side
Definition: dynamic_system_solver.h:669
void initLa()
load contact impulses for use as starting value
Definition: dynamic_system_solver.cc:1159
bool term
boolean to check for termination of contact equations solution
Definition: dynamic_system_solver.h:715
void readz0(std::string fileName)
reads state from a file
Definition: dynamic_system_solver.cc:1339
bool dropContactInfo
flag if contact matrices for debugging should be dropped in no-convergence case
Definition: dynamic_system_solver.h:750
double tolProj
Tolerance for projection of generalized position.
Definition: dynamic_system_solver.h:807
solver interface for modelling and simulation of dynamic systems
Definition: dynamic_system_solver.h:48
void initla()
load contact forces for use as starting value
Definition: dynamic_system_solver.cc:1149
fmatvec::SparseMat Gs
sparse mass action matrix
Definition: dynamic_system_solver.h:695
Strategy strategy
relaxarion strategy for solution of fixed-point scheme
Definition: dynamic_system_solver.h:735
void setgddTol(double tol)
Definition: dynamic_system_solver.h:459
fmatvec::VecInt checkTermLevels
TODO.
Definition: dynamic_system_solver.h:770
DynamicSystemSolver(const std::string &name="")
constructor
Definition: dynamic_system_solver.cc:62
Solver contactSolver
solver for contact equations and impact equations
Definition: dynamic_system_solver.h:730
fmatvec::Vec bc
TODO.
Definition: dynamic_system_solver.h:710
void updateWRef(const fmatvec::Mat &ref, int i=0)
references to contact force direction matrix of dynamic system parent
Definition: dynamic_system_solver.h:556
void setgdTol(double tol)
Definition: dynamic_system_solver.h:454
virtual void plot()
plots time dependent data
Definition: dynamic_system.cc:330
bool stopIfNoConvergence
flag if the contact equations should be stopped if there is no convergence
Definition: dynamic_system_solver.h:745
virtual const fmatvec::Vec & shift()
update for event driven integrator for event
Definition: dynamic_system_solver.cc:1497
virtual void updateWRef(const fmatvec::Mat &ref, int i=0)
references to contact force direction matrix of dynamic system parent
Definition: dynamic_system.cc:824
void updategRef(const fmatvec::Vec &ref)
references to relative distances of dynamic system parent
Definition: dynamic_system_solver.h:537
void setgTol(double tol)
Definition: dynamic_system_solver.h:449
fmatvec::Mat TParent
matrix of linear relation between differentiated positions and velocities
Definition: dynamic_system_solver.h:602
fmatvec::Vec sParent
TODO.
Definition: dynamic_system_solver.h:632
fmatvec::Vec u
velocities, differentiated velocities, initial velocities
Definition: dynamic_system.h:836
void plot()
plots time dependent data
Definition: dynamic_system_solver.cc:1662
fmatvec::VecInt LinkStatusRegParent
status vector of single valued links
Definition: dynamic_system_solver.h:690
void getLinkStatusReg(fmatvec::VecInt &LinkStatusRegExt)
Definition: dynamic_system_solver.cc:1019
fmatvec::Vec gdParent
relative velocities
Definition: dynamic_system_solver.h:647
void projectGeneralizedVelocities(int mode)
drift projection for positions
Definition: dynamic_system_solver.cc:1104
fmatvec::VecInt LinkStatusParent
status vector of set valued links with piecewise link equation (which piece is valid) ...
Definition: dynamic_system_solver.h:684
virtual void updateLLM()
compute Cholesky decomposition of mass matrix TODO necessary?
Definition: dynamic_system_solver.cc:814
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:64
virtual void closePlot()
closes plot file
Definition: dynamic_system_solver.cc:862
unsigned int flushCount
counts plot-calls until files to be flushed TODO
Definition: dynamic_system_solver.h:802
std::string getSolverInfo()
Definition: dynamic_system_solver.cc:1216
virtual void setgTol(double tol)
Definition: dynamic_system.cc:1281
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:725
void addElement(Element *element_)
Definition: dynamic_system_solver.cc:1177
virtual void setLaTol(double tol)
Definition: dynamic_system.cc:1309
fmatvec::Vec rFactorParent
relaxation parameters for contact equations
Definition: dynamic_system_solver.h:627
fmatvec::VecInt decreaseLevels
decreasing relaxation factors is done in levels containing the number of contact iterations as condit...
Definition: dynamic_system_solver.h:765
Strategy
relaxation strategies in solution of contact equations
Definition: dynamic_system_solver.h:59
virtual void updatelaInverseKinetics()
update inverse kinetics constraint forces
Definition: dynamic_system_solver.cc:1590
void setStopIfNoConvergence(bool flag, bool dropInfo=false)
Definition: dynamic_system_solver.h:396
virtual void updategdRef(const fmatvec::Vec &ref)
references to relative velocities of dynamic system parent
Definition: dynamic_system.cc:786
virtual void updater(int j=0)
update smooth link force law
Definition: dynamic_system_solver.cc:819
void dropContactMatrices()
Definition: dynamic_system_solver.cc:1255
void initialize()
Initialize the system.
Definition: dynamic_system_solver.cc:82
virtual void updateG()
updates mass action matrix
Definition: dynamic_system_solver.cc:895
bool numJac
flag if Jacobian for Newton scheme should be calculated numerically
Definition: dynamic_system_solver.h:760
int(DynamicSystemSolver::* solveImpacts_)()
function pointer for election of prox-solver for impact equations on velocity level ...
Definition: dynamic_system_solver.h:271
unsigned int flushEvery
flushes all hdf5-files every x-times the plot-routine is called TODO
Definition: dynamic_system_solver.h:796
fmatvec::Mat WParent[2]
contact force directions
Definition: dynamic_system_solver.h:607
int warnLevel
level for warning output (0-2)
Definition: dynamic_system_solver.h:785
virtual void updategRef(const fmatvec::Vec &ref)
references to relative distances of dynamic system parent
Definition: dynamic_system.cc:779
double evalPotentialEnergy()
compute potential energy of entire dynamic system change? TODO
Definition: dynamic_system_solver.cc:1164
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
void setProjectionTolerance(double tol)
set tolerance for projection of generalized position
Definition: dynamic_system_solver.h:444
int(DynamicSystemSolver::* solveConstraints_)()
Definition: dynamic_system_solver.h:264
Element * getElement(const std::string &name)
Definition: dynamic_system_solver.cc:1189
bool positionDriftCompensationNeeded(double gmax)
check if drift compensation on position level is needed
Definition: dynamic_system_solver.cc:1027
fmatvec::VecInt jsvParent
boolean evaluation of stopvector
Definition: dynamic_system_solver.h:679
void savela()
save contact forces for use as starting value in next time step
Definition: dynamic_system_solver.cc:1144
virtual int solveImpactsRootFinding()
solve impact equations with Newton scheme on velocity level
Definition: dynamic_system_solver.cc:620
LinAlg linAlg
linear system solver used for Newton scheme in contact equations
Definition: dynamic_system_solver.h:740
virtual void setgdTol(double tol)
Definition: dynamic_system.cc:1288
virtual void updateVRef(const fmatvec::Mat &ref, int i=0)
references to condensed contact force direction matrix of dynamic system parent
Definition: dynamic_system.cc:845
void projectGeneralizedPositions(int mode, bool fullUpdate=false)
drift projection for positions
Definition: dynamic_system_solver.cc:1053
fmatvec::SymMat LLMParent
Cholesky decomposition of mass matrix.
Definition: dynamic_system_solver.h:597
std::string name
name of element
Definition: element.h:298
virtual void checkImpactsForTermination()
validate force laws concerning given tolerances on velocity level
Definition: dynamic_system_solver.cc:706
int limitGSize
TODO.
Definition: dynamic_system_solver.h:780
virtual int solveConstraintsFixpointSingle()
solve contact equations with single step fixed point scheme
Definition: dynamic_system_solver.cc:440
virtual void updaterFactors()
update relaxation factors for contact equations
Definition: dynamic_system_solver.cc:1371
virtual int solveImpactsFixpointSingle()
solve impact equations with single step fixed point scheme on velocity level
Definition: dynamic_system_solver.cc:472
static bool exitRequest
boolean signal evaluation for end integration set by user
Definition: dynamic_system_solver.h:855
fmatvec::SymMat MParent
mass matrix
Definition: dynamic_system_solver.h:592
bool peds
TODO, flag for occuring impact and sticking in event driven solver.
Definition: dynamic_system_solver.h:790
fmatvec::Vec hParent[2]
smooth, smooth with respect to objects, smooth with respect to links right hand side ...
Definition: dynamic_system_solver.h:664
void constructor()
set plot feature default values
Definition: dynamic_system_solver.cc:1393
fmatvec::SqrMat G
mass action matrix
Definition: dynamic_system_solver.h:705
double dt
step size
Definition: dynamic_system_solver.h:587
void setlaTol(double tol)
Definition: dynamic_system_solver.h:464
void saveLa()
save contact impulses for use as starting value in next time step
Definition: dynamic_system_solver.cc:1154
bool integratorExitRequest
boolean signal evaluation for end integration set by program
Definition: dynamic_system_solver.h:850
fmatvec::Vec wbParent
TODO.
Definition: dynamic_system_solver.h:617
void addToGraph(Graph *graph, fmatvec::SqrMat &A, int i, std::vector< Element * > &objList)
adds list of objects to tree
Definition: dynamic_system_solver.cc:1484
virtual void checkConstraintsForTermination()
validate force laws concerning given tolerances
Definition: dynamic_system_solver.cc:695
void updategdRef(const fmatvec::Vec &ref)
references to relative velocities of dynamic system parent
Definition: dynamic_system_solver.h:543
void init(InitStage stage)
plots time series header
Definition: dynamic_system_solver.cc:101
fmatvec::Vec la
contact force parameters
Definition: dynamic_system.h:856
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:720
void getLinkStatus(fmatvec::VecInt &LinkStatusExt)
Definition: dynamic_system_solver.cc:1011
Solver
solver for contact equations
Definition: dynamic_system_solver.h:54
fmatvec::Vec zParent
state
Definition: dynamic_system_solver.h:652
virtual int solveImpactsGaussSeidel()
solve impact equations with Gauss-Seidel scheme on velocity level
Definition: dynamic_system_solver.cc:524
fmatvec::Vec laParent
contact force parameters
Definition: dynamic_system_solver.h:622
virtual std::string getType() const
Definition: dynamic_system_solver.h:177
fmatvec::Vec resParent
residuum of contact equations
Definition: dynamic_system_solver.h:637
void setLaTol(double tol)
Definition: dynamic_system_solver.h:469
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:487
virtual ~DynamicSystemSolver()
destructor
Definition: dynamic_system_solver.cc:72
bool velocityDriftCompensationNeeded(double gdmax)
check if drift compensation on velocity level is needed
Definition: dynamic_system_solver.cc:1040
fmatvec::Vec svParent
stopvector (rootfunctions for event driven integration
Definition: dynamic_system_solver.h:674
fmatvec::SqrMat Jprox
JACOBIAN of contact equations for Newton scheme.
Definition: dynamic_system_solver.h:700
fmatvec::Mat VParent[2]
condensed contact force directions
Definition: dynamic_system_solver.h:612
void writez(std::string fileName, bool formatH5=true)
writes state to a file
Definition: dynamic_system_solver.cc:1319
virtual int solveConstraintsGaussSeidel()
solve contact equations with Gauss-Seidel scheme
Definition: dynamic_system_solver.cc:504
static void sigAbortHandler(int)
handler for abort signals
Definition: dynamic_system_solver.cc:1290
fmatvec::SymMat M
mass matrix
Definition: dynamic_system.h:821
void decreaserFactors()
decrease relaxation factors if mass action matrix is not diagonal dominant
Definition: dynamic_system_solver.cc:990
void computeInitialCondition()
compute initial condition for links for event driven integrator
Definition: dynamic_system_solver.cc:868
void updatezdRef(const fmatvec::Vec &ext)
references to differentiated external state
Definition: dynamic_system_solver.cc:1358
bool READZ0
is a state read from a file
Definition: dynamic_system_solver.h:860
virtual void updatewbRef(const fmatvec::Vec &ref)
references to TODO of dynamic system parent
Definition: dynamic_system.cc:817
virtual void setgddTol(double tol)
Definition: dynamic_system.cc:1295
fmatvec::Vec zdParent
differentiated state
Definition: dynamic_system_solver.h:657
double evalKineticEnergy()
compute kinetic energy of entire dynamic system
Definition: dynamic_system_solver.h:366

Impressum / Disclaimer / Datenschutz Generated by doxygen 1.8.5 Valid HTML