mbsim  4.0.0
MBSim Kernel
dynamic_system_solver.h
1/* Copyright (C) 2004-2018 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#include "mbsim/functions/function.h"
26#include "mbsim/environment.h"
27
28#include <atomic>
29
30namespace MBSim {
31
32 class Graph;
33 class MBSimEnvironment;
34 class MultiDimNewtonMethod;
35 class ConstraintResiduum;
36 class ConstraintJacobian;
37
38 struct StateTable {
39 std::string name;
40 char label;
41 int number;
42 StateTable() = default;
43 StateTable(const std::string &name_, char label_, int number_) : name(name_), label(label_), number(number_) { }
44 };
45
61 class DynamicSystemSolver : public Group {
62
63 class Residuum : public MBSim::Function<fmatvec::Vec(fmatvec::Vec)> {
64 public:
65 Residuum(MBSim::DynamicSystemSolver *sys_) : sys(sys_) { }
66 fmatvec::Vec operator()(const fmatvec::Vec &z);
67 private:
69 };
70
71 public:
72
76 enum Solver { fixedpoint, GaussSeidel, direct, rootfinding, unknownSolver, directNonlinear };
77
82 DynamicSystemSolver(const std::string &name="");
83
87 ~DynamicSystemSolver() override;
88
127 void initialize();
128
129 /*
130 * If true (the default) then
131 * the simulation output files (h5 files) are deleted/truncated/regenerated.
132 * If false, then these files left are untouched from a previous run.
133 * This is useful to regenerate the
134 * e.g. OpenMBV XML files without doing a new integration.
135 */
136 void setTruncateSimulationFiles(bool trunc) { truncateSimulationFiles=trunc; }
137
138 /* INHERITED INTERFACE OF GROUP */
139 void init(InitStage stage, const InitConfigSet &config) override;
140 using Group::plot;
141 /***************************************************/
142
143 /* INHERITED INTERFACE OF DYNAMICSYSTEM */
144 int solveConstraintsFixpointSingle() override;
145 int solveImpactsFixpointSingle() override;
146 int solveConstraintsGaussSeidel() override;
147 int solveImpactsGaussSeidel() override;
148 int solveConstraintsRootFinding() override;
149 int solveImpactsRootFinding() override;
150 void checkConstraintsForTermination() override;
151 void checkImpactsForTermination() override;
152 /***************************************************/
153
154 /* INHERITED INTERFACE OF OBJECTINTERFACE */
155 void updateT() override;
156 void updateh(int i=0) override;
157 void updateM() override;
158 void updateLLM() override;
159 void updatezd() override;
160 /***************************************************/
161
162 /* INHERITED INTERFACE OF LINKINTERFACE */
163
174 void updater(int j=0) override;
175 void updateJrla(int j=0) override;
176 virtual void updaterdt();
177 void updatewb() override;
178 void updateg() override;
179 void updategd() override;
180 void updateW(int j=0) override;
181 void updateV(int j=0) override;
182 virtual void updatebc();
183 virtual void updatebi();
184 virtual void updatela();
185 virtual void updateLa();
186 /***************************************************/
187
188 /* INHERITED INTERFACE OF ELEMENT */
189 /***************************************************/
190
191 /* GETTER / SETTER */
192
193 void setSmoothSolver(Solver solver_) { smoothSolver = solver_; }
194 void setConstraintSolver(Solver solver_) { contactSolver = solver_; }
195 void setImpactSolver(Solver solver_) { impactSolver = solver_; }
196 const Solver& getSmoothSolver() { return smoothSolver; }
197 const Solver& getConstraintSolver() { return contactSolver; }
198 const Solver& getImpactSolver() { return impactSolver; }
199 void setTermination(bool term_) { term = term_; }
200 void setMaximumNumberOfIterations(int iter) { maxIter = iter; }
201 void setHighNumberOfIterations(int iter) { highIter = iter; }
202 void setNumericalJacobian(bool numJac_) { numJac = numJac_; }
203 void setMaximumDampingSteps(int maxDSteps) { maxDampingSteps = maxDSteps; }
204 void setLevenbergMarquardtParamater(double lmParm_) { lmParm = lmParm_; }
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 double& getTime() { return t; }
213 const double& getTime() const { return t; }
214 void setTime(double t_) { t = t_; }
215
216 double getStepSize() const { return dt; }
217 void setStepSize(double dt_) { dt = dt_; }
218
219 int getzSize() const { return zSize; }
220
221 fmatvec::Vec& getState() { return z; }
222 const fmatvec::Vec& getState() const { return z; }
223 void setState(const fmatvec::Vec &z_) { z = z_; }
224
225 const fmatvec::Vec& getzd(bool check=true) const { assert((not check) or (not updzd)); return zd; }
226 void setzd(const fmatvec::Vec &zd_) { zd = zd_; }
227
228 const fmatvec::SqrMat& getG(bool check=true) const { assert((not check) or (not updG)); return G; }
229 const fmatvec::SparseMat& getGs(bool check=true) const { assert((not check) or (not updG)); return Gs; }
230 const fmatvec::Vec& getbc(bool check=true) const { assert((not check) or (not updbc)); return bc; }
231 const fmatvec::Vec& getbi(bool check=true) const { assert((not check) or (not updbi)); return bi; }
232 const fmatvec::SqrMat& getJprox() const { return Jprox; }
233 fmatvec::SqrMat& getG(bool check=true) { assert((not check) or (not updG)); return G; }
234 fmatvec::SparseMat& getGs(bool check=true) { assert((not check) or (not updG)); return Gs; }
235 fmatvec::Vec& getbc(bool check=true) { assert((not check) or (not updbc)); return bc; }
236 fmatvec::Vec& getbi(bool check=true) { assert((not check) or (not updbi)); return bi; }
237 fmatvec::SqrMat& getJprox() { return Jprox; }
238 const fmatvec::Vec& getr(int j=0, bool check=true) const { assert((not check) or (not updr[j])); return r[j]; }
239 fmatvec::Vec& getr(int j=0, bool check=true) { assert((not check) or (not updr[j])); return r[j]; }
240
241 // Jacobian of dr/dla
242 const fmatvec::Mat& getJrla(int j=0, bool check=true) const { assert((not check) or (not updJrla[j])); return Jrla[j]; }
243 fmatvec::Mat& getJrla(int j=0, bool check=true) { assert((not check) or (not updJrla[j])); return Jrla[j]; }
244
245 const fmatvec::Vec& evaldq() { if(upddq) updatedq(); return dq; }
246 const fmatvec::Vec& evaldu() { if(upddu) updatedu(); return du; }
247 const fmatvec::Vec& evaldx() { if(upddx) updatedx(); return dx; }
248 const fmatvec::Vec& evalzd();
249 const fmatvec::SqrMat& evalG() { if(updG) updateG(); return G; }
250 const fmatvec::SparseMat& evalGs() { if(updG) updateG(); return Gs; }
251 const fmatvec::Vec& evalbc() { if(updbc) updatebc(); return bc; }
252 const fmatvec::Vec& evalbi() { if(updbi) updatebi(); return bi; }
253 const fmatvec::Vec& evalsv();
254 const fmatvec::Vec& evalz0();
255 const fmatvec::Vec& evalla() { if(updla) updatela(); return la; }
256 const fmatvec::Vec& evalLa() { if(updLa) updateLa(); return La; }
257
258 fmatvec::Vec& getzParent() { return zParent; }
259 fmatvec::Vec& getzdParent() { return zdParent; }
260 fmatvec::Vec& getlaParent() { return laParent; }
261 fmatvec::Vec& getLaParent() { return LaParent; }
262// const fmatvec::Vec& getzParent() const { return zParent; }
263// const fmatvec::Vec& getzdParent() const { return zdParent; }
264// const fmatvec::Mat& getWParent(int i=0) const { return WParent[i]; }
265 fmatvec::Mat& getWParent(int i=0) { return WParent[i]; }
266 fmatvec::Mat& getVParent(int i=0) { return VParent[i]; }
267// fmatvec::Vec& getlaParent() { return laParent; }
268// fmatvec::Vec& getLaParent() { return LaParent; }
269// const fmatvec::Vec& getgParent() const { return gParent; }
270 fmatvec::Vec& getgParent() { return gParent; }
271// const fmatvec::Vec& getgdParent() const { return gdParent; }
272 fmatvec::Vec& getgdParent() { return gdParent; }
273 fmatvec::Vec& getresParent() { return resParent; }
274 fmatvec::Vec& getrFactorParent() { return rFactorParent; }
275
276 void resizezParent(int nz) { zParent.resize(nz); }
277 void resizezdParent(int nz) { zdParent.resize(nz); }
278
279 DynamicSystemSolver* getDynamicSystemSolver() { return this; }
280 int getMaxIter() {return maxIter;}
281 int getHighIter() {return highIter;}
282 int getIterC() {return iterc;}
283 int getIterI() {return iteri;}
284 /***************************************************/
285
290
291 int (DynamicSystemSolver::*solveSmooth_)();
292
298
305
312
319
327
328 int solveImpactsNonlinearEquations();
329
334 virtual void updateG();
335
339 void decreaserFactors();
340
347 virtual const fmatvec::Vec& shift(std::optional<std::reference_wrapper<bool>> &&velProjWasCalled={},
348 std::optional<std::reference_wrapper<bool>> &&posProjWasCalled={});
349
354 void getLinkStatus(fmatvec::VecInt &LinkStatusExt);
355
360 void getLinkStatusReg(fmatvec::VecInt &LinkStatusRegExt);
361
365 bool positionDriftCompensationNeeded(double gmax);
366
370 bool velocityDriftCompensationNeeded(double gdmax);
371
375 void projectGeneralizedPositions(int mode, bool fullUpdate=false);
376
380 void projectGeneralizedVelocities(int mode);
381
386 void savela();
387
392 void initla();
393
398 void saveLa();
399
404 void initLa();
405
410 void addElement(Element *element_);
411
417 Element* getElement(const std::string &name);
418
422 std::string getSolverInfo();
423
428 void setStopIfNoConvergence(bool flag, bool dropInfo = false) { stopIfNoConvergence = flag; dropContactInfo=dropInfo; }
429
430 bool getStopIfNoConvergence() { return stopIfNoConvergence; }
431
435 void dropContactMatrices();
436
437 // MBSim signal handler
438#ifndef SWIG
440 public:
441 // This resets the DSS flags exitRequest and silentExit to false
444 private:
445 using SigHandle = void (*)(int);
446 #ifndef _WIN32
447 SigHandle oldSigHup;
448 #endif
449 SigHandle oldSigInt;
450 SigHandle oldSigTerm;
451 };
452#endif
453
454 // This function should be called repeatetly in mbsim, at least in heavy work parts of the code.
455 // If the DSS exitRequest flag is set a error message is printed and a std::runtime_error exception is thrown to interrupt.
456 // If silent is true or the DSS silentExit flag is set then not message is printed and a SilentError is thrown.
457 static void throwIfExitRequested(bool silent=false) {
458 if(exitRequest) {
459 if(!silent && !silentExit) {
460 msgStatic(fmatvec::Atom::Error)<<"User requested a exit (throw exception now)."<<std::endl;
461 throw std::runtime_error("Exception due to user requested exit.");
462 }
463 msgStatic(fmatvec::Atom::Info)<<"User requested a silent exit (throw exception now)."<<std::endl;
464 throw SilentError();
465 }
466 }
467
468 // TODO just for testing
469 void setPartialEventDrivenSolver(bool peds_) { peds = peds_; }
470
476 void writez(std::string fileName, bool formatH5=true);
477 using DynamicSystem::writez;
478
484 void writeStateTable(std::string fileName);
485
490 void readz0(std::string fileName);
491 using DynamicSystem::readz0;
492
493 void initializeUsingXML(xercesc::DOMElement *element) override;
494
499 void setProjectionTolerance(double tol) { tolProj = tol; }
500
505 void setLocalSolverTolerance(double tol) { tolLocalSolver = tol; }
506 double getLocalSolverTolerance() { return tolLocalSolver; }
507
512 void setDynamicSystemSolverTolerance(double tol) { tolDSS = tol; }
513 double getDynamicSystemSolverTolerance() { return tolDSS; }
514
519 void updatezRef(fmatvec::Vec &ext);
520
525 void updatezdRef(fmatvec::Vec &ext);
526
527 void setAlwaysConsiderContact(bool alwaysConsiderContact_) { alwaysConsiderContact = alwaysConsiderContact_; }
528
529 void setInverseKinetics(bool inverseKinetics_) { inverseKinetics = inverseKinetics_; }
530 bool getInverseKinetics() const { return inverseKinetics; }
531
532 void setInitialProjection(bool initialProjection_) { initialProjection = initialProjection_; }
533 bool getInitialProjection() const { return initialProjection; }
534
535 void setDetermineEquilibriumState(bool determineEquilibriumState_) { determineEquilibriumState = determineEquilibriumState_; }
536 bool getDetermineEquilibriumState() const { return determineEquilibriumState; }
537
538 void setUseConstraintSolverForPlot(bool useConstraintSolverForPlot_) { useConstraintSolverForPlot = useConstraintSolverForPlot_; }
539 bool getUseConstraintSolverForPlot() const { return useConstraintSolverForPlot; }
540
541 fmatvec::Mat dhdq(int lb=0, int ub=0);
542 fmatvec::Mat dhdu(int lb=0, int ub=0);
543 fmatvec::Mat dhdx();
544 fmatvec::Vec dhdt();
545
546 void setRootID(int ID) {rootID = ID;}
547 int getRootID() const {return rootID;}
548
549 void resetUpToDate() override;
550
551 bool getUpdateT() { return updT; }
552 bool getUpdateM() { return updM; }
553 bool getUpdateLLM() { return updLLM; }
554 bool getUpdateh(int j) { return updh[j]; }
555 bool getUpdater(int j) { return updr[j]; }
556 bool getUpdateJrla(int j) { return updJrla[j]; }
557 bool getUpdaterdt() { return updrdt; }
558 bool getUpdateW(int j) { return updW[j]; }
559 bool getUpdateV(int j) { return updV[j]; }
560 bool getUpdatewb() { return updwb; }
561 bool getUpdateg() { return updg; }
562 bool getUpdategd() { return updgd; }
563 bool getUpdatela() { return updla; }
564 bool getUpdateLa() { return updLa; }
565 bool getUpdatezd() { return updzd; }
566 bool getUpdatedq() { return upddq; }
567 bool getUpdatedu() { return upddu; }
568 bool getUpdatedx() { return upddx; }
569 void setUpdatela(bool updla_) { updla = updla_; }
570 void setUpdateLa(bool updLa_) { updLa = updLa_; }
571 void setUpdateG(bool updG_) { updG = updG_; }
572 void setUpdatebi(bool updbi_) { updbi = updbi_; }
573 void setUpdatebc(bool updbc_) { updbc = updbc_; }
574 void setUpdatezd(bool updzd_) { updzd = updzd_; }
575 void setUpdateW(bool updW_, int i=0) { updW[i] = updW_; }
576
581 void updategRef(fmatvec::Vec &ref) override { Group::updategRef(ref); updg = true; }
582
587 void updategdRef(fmatvec::Vec &ref) override { Group::updategdRef(ref); updgd = true; }
588
593 void updatewbRef(fmatvec::Vec &ref) override { Group::updatewbRef(ref); updwb = true; }
594
600 void updateWRef(fmatvec::Mat &ref, int i=0) override { Group::updateWRef(ref,i); updW[i] = true; }
601
607 void updateVRef(fmatvec::Mat &ref, int i=0) override { Group::updateVRef(ref,i); updV[i] = true; }
608
612 virtual void updatelaInverseKinetics();
613
614 void updatedq() override;
615 void updatedu() override;
616 void updatedx() override;
617
618 void updateStopVector() override;
619
620 void updateInternalState();
621
622 void plot() override;
623
624 void addEnvironment(Environment* env);
625
630 template<class Env>
631 Env* getEnvironment();
632
636
637 std::vector<StateTable>& getStateTable() { return tabz; }
638
639 void setCompressionLevel(int level) { compressionLevel=level; }
640 void setChunkSize(int size) { chunkSize=size; }
641 void setCacheSize(int size) { cacheSize=size; }
642
643 void setqdequ(bool qdequ_) { qdequ = qdequ_; }
644 bool getqdequ() { return qdequ; }
645
647 static void interrupt(bool silent=false);
648
649 protected:
653 double t;
654
658 double dt;
659
663 int zSize;
664
668 fmatvec::Vec z;
669
673 fmatvec::Vec zd;
674
678 fmatvec::SymMat MParent;
679
683 fmatvec::SymMat LLMParent;
684
688 fmatvec::Mat TParent;
689
693 fmatvec::Mat WParent[2];
694
698 fmatvec::Mat VParent[2];
699
703 fmatvec::Vec wbParent;
704
708 fmatvec::Vec laParent, LaParent;
709
713 fmatvec::Vec rFactorParent;
714
718 fmatvec::Vec sParent;
719
723 fmatvec::Vec resParent;
724
725 fmatvec::Vec curisParent, nextisParent;
726
730 fmatvec::Vec gParent;
731
735 fmatvec::Vec gdParent;
736
740 fmatvec::Vec zParent;
741
745 fmatvec::Vec zdParent;
746
747 fmatvec::Vec dxParent, dqParent, duParent;
748
752 fmatvec::Vec hParent[2];
753
757 fmatvec::Vec rParent[2], rdtParent;
758
759 fmatvec::Mat JrlaParent[2];
760
764 fmatvec::Vec svParent;
765
769 fmatvec::VecInt jsvParent;
770
774 fmatvec::VecInt LinkStatusParent;
775
780 fmatvec::VecInt LinkStatusRegParent;
781
785 fmatvec::SparseMat Gs;
786
790 fmatvec::SqrMat Jprox;
791
795 fmatvec::SqrMat G;
796
800 fmatvec::Vec bc, bi;
801
805 bool term;
806
810 int maxIter, highIter, maxDampingSteps, iterc, iteri;
811
815 double lmParm;
816
820 Solver smoothSolver, contactSolver, impactSolver;
821
826
831
836
840 bool numJac;
841
845 fmatvec::VecInt decreaseLevels;
846
850 fmatvec::VecInt checkTermLevels;
851
856
861
865 bool peds;
866
870 double tolProj;
871
875 double tolLocalSolver { 1e-10 };
876
880 double tolDSS { 1e-9 };
881
885 fmatvec::Vec laInverseKineticsParent;
886 fmatvec::Mat bInverseKineticsParent;
887
891 fmatvec::Mat WInverseKineticsParent;
892
893 bool alwaysConsiderContact;
894 bool inverseKinetics;
895 bool initialProjection;
896 bool determineEquilibriumState;
897 bool useConstraintSolverForPlot;
898
899 fmatvec::Vec corrParent;
900
901 int rootID;
902
903 double gTol, gdTol, gddTol, laTol, LaTol;
904
905 bool updT, updh[2], updr[2], updJrla[2], updrdt, updM, updLLM, updW[2], updV[2], updwb, updg, updgd, updG, updbc, updbi, updsv, updzd, updla, updLa, upddq, upddu, upddx;
906
907 bool useSmoothSolver;
908
909 bool qdequ;
910
911 std::vector<StateTable> tabz;
912
913 int compressionLevel { H5::File::getDefaultCompression() };
914 int chunkSize { H5::File::getDefaultChunkSize() };
915 int cacheSize { H5::File::getDefaultCacheSize() };
916
917 private:
922 static void sigInterruptHandler(int);
923
927 void constructor();
928
932 static std::atomic<bool> exitRequest;
933 static_assert(decltype(exitRequest)::is_always_lock_free);
934 // flag to indicate if a user requested interruption should be silent.
935 static std::atomic<bool> silentExit;
936
940 bool READZ0;
941
950 void addToGraph(Graph* graph, fmatvec::SqrMat &A, int i, std::vector<Element*> &objList);
951
952 bool truncateSimulationFiles;
953
954 // Holds the dynamic systems before the "reorganize hierarchy" takes place.
955 // This is required since all elements of all other containers from DynamicSystem are readded to DynamicSystemSolver,
956 // except this container (which is a "pure" container = no calculation is done in DynamicSystem)
957 // However, we must hold this container until the dtor of DynamicSystemSolver is called to avoid the deallocation of other
958 // elements hold by DynamicSystem elements (especially (currently only) hdf5File)
959 std::vector<DynamicSystem*> dynamicsystemPreReorganize;
960
961 double facSizeGs;
962
963 std::vector<std::unique_ptr<Environment>> environments;
964 MBSimEnvironment *mbsimEnvironment = nullptr;
965
966 bool firstPlot { true };
967
968 std::unique_ptr<MultiDimNewtonMethod> nonlinearConstraintNewtonSolver;
969 std::unique_ptr<ConstraintResiduum> constraintResiduum;
970 std::unique_ptr<ConstraintJacobian> constraintJacobian;
971 };
972
973 template<class Env>
975 // get the Environment of type Env
976 auto &reqType=typeid(Env);
977 for(auto &e : environments) {
978 auto &e_=*e;
979 if(reqType==typeid(e_))
980 return static_cast<Env*>(e.get());
981 }
982 auto newEnv=new Env;
983 addEnvironment(newEnv);
984 return newEnv;
985 }
986
987}
988
989#endif /* _DYNAMIC_SYSTEM_SOLVER_H_ */
Definition: dynamic_system_solver.h:63
Definition: dynamic_system_solver.h:439
solver interface for modelling and simulation of dynamic systems
Definition: dynamic_system_solver.h:61
void checkConstraintsForTermination() override
validate force laws concerning given tolerances
Definition: dynamic_system_solver.cc:726
fmatvec::Vec sParent
TODO.
Definition: dynamic_system_solver.h:718
int(DynamicSystemSolver::* solveImpacts_)()
function pointer for election of prox-solver for impact equations on velocity level
Definition: dynamic_system_solver.h:304
void plot() override
plots time dependent data
Definition: dynamic_system_solver.cc:1761
bool term
boolean to check for termination of contact equations solution
Definition: dynamic_system_solver.h:805
int solveConstraintsFixpointSingle() override
solve contact equations with single step fixed point scheme
Definition: dynamic_system_solver.cc:494
void computeInitialCondition()
compute initial condition for links for event driven integrator
Definition: dynamic_system_solver.cc:901
fmatvec::Vec hParent[2]
smooth, smooth with respect to objects, smooth with respect to links right hand side
Definition: dynamic_system_solver.h:752
void setProjectionTolerance(double tol)
set tolerance for projection of generalized position
Definition: dynamic_system_solver.h:499
fmatvec::Vec zd
derivative of state vector
Definition: dynamic_system_solver.h:673
bool useOldla
flag if contac force parameter of last time step should be used
Definition: dynamic_system_solver.h:835
fmatvec::Mat WParent[2]
contact force directions
Definition: dynamic_system_solver.h:693
fmatvec::VecInt LinkStatusParent
status vector of set valued links with piecewise link equation (which piece is valid)
Definition: dynamic_system_solver.h:774
fmatvec::Vec bc
TODO.
Definition: dynamic_system_solver.h:800
int solveConstraintsLinearEquations()
solution of contact equations with direct linear solver using minimum of least squares
Definition: dynamic_system_solver.cc:951
void setLocalSolverTolerance(double tol)
set tolerance for local none-linear solver (solvers on element level), like the Newton-Solver in Join...
Definition: dynamic_system_solver.h:505
fmatvec::SqrMat Jprox
JACOBIAN of contact equations for Newton scheme.
Definition: dynamic_system_solver.h:790
void initialize()
Initialize the system.
Definition: dynamic_system_solver.cc:105
void updateWRef(fmatvec::Mat &ref, int i=0) override
references to contact force direction matrix of dynamic system parent
Definition: dynamic_system_solver.h:600
fmatvec::Mat VParent[2]
condensed contact force directions
Definition: dynamic_system_solver.h:698
double lmParm
Levenberg-Marquard parameter.
Definition: dynamic_system_solver.h:815
fmatvec::Vec gParent
relative distances
Definition: dynamic_system_solver.h:730
virtual void updatelaInverseKinetics()
update inverse kinetics constraint forces
Definition: dynamic_system_solver.cc:1681
fmatvec::Vec z
state vector
Definition: dynamic_system_solver.h:668
bool peds
TODO, flag for occuring impact and sticking in event driven solver.
Definition: dynamic_system_solver.h:865
static void sigInterruptHandler(int)
handler for user interrupt signal Set the DSS flag exitRequest to true.
Definition: dynamic_system_solver.cc:1369
void writeStateTable(std::string fileName)
writes state table to a file
Definition: dynamic_system_solver.cc:1400
bool numJac
flag if Jacobian for Newton scheme should be calculated numerically
Definition: dynamic_system_solver.h:840
void saveLa()
save contact impulses for use as starting value in next time step
Definition: dynamic_system_solver.cc:1264
fmatvec::Vec wbParent
TODO.
Definition: dynamic_system_solver.h:703
double t
time
Definition: dynamic_system_solver.h:653
Solver smoothSolver
solver for contact equations and impact equations
Definition: dynamic_system_solver.h:820
int solveConstraintsRootFinding() override
solve contact equations with Newton scheme
Definition: dynamic_system_solver.cc:594
double tolProj
Tolerance for projection of generalized position.
Definition: dynamic_system_solver.h:870
void updater(int j=0) override
update smooth link force law
Definition: dynamic_system_solver.cc:850
fmatvec::Vec resParent
residuum of contact equations
Definition: dynamic_system_solver.h:723
int(DynamicSystemSolver::* solveConstraints_)()
Definition: dynamic_system_solver.h:297
bool dropContactInfo
flag if contact matrices for debugging should be dropped in no-convergence case
Definition: dynamic_system_solver.h:830
fmatvec::VecInt decreaseLevels
decreasing relaxation factors is done in levels containing the number of contact iterations as condit...
Definition: dynamic_system_solver.h:845
fmatvec::VecInt checkTermLevels
TODO.
Definition: dynamic_system_solver.h:850
bool READZ0
is a state read from a file
Definition: dynamic_system_solver.h:940
int solveConstraintsGaussSeidel() override
solve contact equations with Gauss-Seidel scheme
Definition: dynamic_system_solver.cc:554
fmatvec::SqrMat G
mass action matrix
Definition: dynamic_system_solver.h:795
std::string getSolverInfo()
Definition: dynamic_system_solver.cc:1313
int solveConstraintsNonlinearEquations()
solution of contact equations with direct nonlinear newton solver using minimum of least squares in e...
Definition: dynamic_system_solver.cc:970
fmatvec::VecInt LinkStatusRegParent
status vector of single valued links
Definition: dynamic_system_solver.h:780
void updategRef(fmatvec::Vec &ref) override
references to relative distances of dynamic system parent
Definition: dynamic_system_solver.h:581
void updatezRef(fmatvec::Vec &ext)
references to external state
Definition: dynamic_system_solver.cc:1419
void setStopIfNoConvergence(bool flag, bool dropInfo=false)
Definition: dynamic_system_solver.h:428
int solveImpactsLinearEquations()
solution of contact equations with Cholesky decomposition on velocity level
Definition: dynamic_system_solver.cc:980
void writez(std::string fileName, bool formatH5=true)
writes state to a file
Definition: dynamic_system_solver.cc:1378
void addToGraph(Graph *graph, fmatvec::SqrMat &A, int i, std::vector< Element * > &objList)
adds list of objects to tree
Definition: dynamic_system_solver.cc:1564
int solveImpactsRootFinding() override
solve impact equations with Newton scheme on velocity level
Definition: dynamic_system_solver.cc:661
static std::atomic< bool > exitRequest
boolean flags to store a user requested interruption
Definition: dynamic_system_solver.h:932
~DynamicSystemSolver() override
destructor
Definition: dynamic_system_solver.cc:97
double tolLocalSolver
Tolerance for local none-linear solvers (solvers on element level)
Definition: dynamic_system_solver.h:875
void updatezdRef(fmatvec::Vec &ext)
references to differentiated external state
Definition: dynamic_system_solver.cc:1431
void addElement(Element *element_)
Definition: dynamic_system_solver.cc:1274
void checkImpactsForTermination() override
validate force laws concerning given tolerances on velocity level
Definition: dynamic_system_solver.cc:737
bool positionDriftCompensationNeeded(double gmax)
check if drift compensation on position level is needed
Definition: dynamic_system_solver.cc:1112
fmatvec::Vec gdParent
relative velocities
Definition: dynamic_system_solver.h:735
DynamicSystemSolver(const std::string &name="")
constructor
Definition: dynamic_system_solver.cc:86
MBSimEnvironment * getMBSimEnvironment()
Definition: dynamic_system_solver.cc:1804
Env * getEnvironment()
Definition: dynamic_system_solver.h:974
fmatvec::Vec zParent
state
Definition: dynamic_system_solver.h:740
bool stopIfNoConvergence
flag if the contact equations should be stopped if there is no convergence
Definition: dynamic_system_solver.h:825
fmatvec::Vec laParent
contact force parameters
Definition: dynamic_system_solver.h:708
void projectGeneralizedPositions(int mode, bool fullUpdate=false)
drift projection for positions
Definition: dynamic_system_solver.cc:1138
fmatvec::SymMat LLMParent
Cholesky decomposition of mass matrix.
Definition: dynamic_system_solver.h:683
fmatvec::Vec rFactorParent
relaxation parameters for contact equations
Definition: dynamic_system_solver.h:713
void init(InitStage stage, const InitConfigSet &config) override
plots time series header
Definition: dynamic_system_solver.cc:122
Element * getElement(const std::string &name)
Definition: dynamic_system_solver.cc:1286
fmatvec::Vec svParent
stopvector (rootfunctions for event driven integration
Definition: dynamic_system_solver.h:764
void constructor()
set plot feature default values
Definition: dynamic_system_solver.cc:1445
double dt
step size
Definition: dynamic_system_solver.h:658
int zSize
size of state vector
Definition: dynamic_system_solver.h:663
fmatvec::Mat TParent
matrix of linear relation between differentiated positions and velocities
Definition: dynamic_system_solver.h:688
void updategdRef(fmatvec::Vec &ref) override
references to relative velocities of dynamic system parent
Definition: dynamic_system_solver.h:587
virtual const fmatvec::Vec & shift(std::optional< std::reference_wrapper< bool > > &&velProjWasCalled={}, std::optional< std::reference_wrapper< bool > > &&posProjWasCalled={})
update for event driven integrator for event
Definition: dynamic_system_solver.cc:1577
virtual void updateG()
updates mass action matrix
Definition: dynamic_system_solver.cc:989
void updateVRef(fmatvec::Mat &ref, int i=0) override
references to condensed contact force direction matrix of dynamic system parent
Definition: dynamic_system_solver.h:607
fmatvec::Vec zdParent
differentiated state
Definition: dynamic_system_solver.h:745
void setDynamicSystemSolverTolerance(double tol)
set tolerance for global none-linear solver (solvers on DynamicSystemSolver level),...
Definition: dynamic_system_solver.h:512
void initla()
load contact forces for use as starting value
Definition: dynamic_system_solver.cc:1259
void readz0(std::string fileName)
reads state from a file
Definition: dynamic_system_solver.cc:1411
fmatvec::SparseMat Gs
sparse mass action matrix
Definition: dynamic_system_solver.h:785
void updatewbRef(fmatvec::Vec &ref) override
references to TODO of dynamic system parent
Definition: dynamic_system_solver.h:593
fmatvec::SymMat MParent
mass matrix
Definition: dynamic_system_solver.h:678
void dropContactMatrices()
Definition: dynamic_system_solver.cc:1330
static void interrupt(bool silent=false)
Definition: dynamic_system_solver.cc:1373
bool velocityDriftCompensationNeeded(double gdmax)
check if drift compensation on velocity level is needed
Definition: dynamic_system_solver.cc:1125
fmatvec::VecInt jsvParent
boolean evaluation of stopvector
Definition: dynamic_system_solver.h:769
int solveImpactsFixpointSingle() override
solve impact equations with single step fixed point scheme on velocity level
Definition: dynamic_system_solver.cc:524
fmatvec::Vec rParent[2]
nonsmooth right hand side
Definition: dynamic_system_solver.h:757
Solver
solver for contact equations
Definition: dynamic_system_solver.h:76
double tolDSS
Tolerance for global none-linear solvers (solvers on DynamicSystemSolver level)
Definition: dynamic_system_solver.h:880
int maxIter
maximum number of contact iterations, high number of contact iterations for warnings,...
Definition: dynamic_system_solver.h:810
int solveImpactsGaussSeidel() override
solve impact equations with Gauss-Seidel scheme on velocity level
Definition: dynamic_system_solver.cc:574
void projectGeneralizedVelocities(int mode)
drift projection for positions
Definition: dynamic_system_solver.cc:1205
bool checkGSize
boolean if force action matrix should be resized in each step
Definition: dynamic_system_solver.h:855
void initLa()
load contact impulses for use as starting value
Definition: dynamic_system_solver.cc:1269
void savela()
save contact forces for use as starting value in next time step
Definition: dynamic_system_solver.cc:1254
int limitGSize
TODO.
Definition: dynamic_system_solver.h:860
void decreaserFactors()
decrease relaxation factors if mass action matrix is not diagonal dominant
Definition: dynamic_system_solver.cc:1090
virtual void updatewbRef(fmatvec::Vec &wbParent)
references to TODO of dynamic system parent
Definition: dynamic_system.cc:854
void plot() override
plots time dependent data
Definition: dynamic_system.cc:369
fmatvec::Mat Jrla[2]
Jacobian dr/dla.
Definition: dynamic_system.h:868
fmatvec::Vec la
contact force parameters
Definition: dynamic_system.h:878
virtual void updateVRef(fmatvec::Mat &VParent, int j=0)
references to condensed contact force direction matrix of dynamic system parent
Definition: dynamic_system.cc:882
virtual void updateWRef(fmatvec::Mat &WParent, int j=0)
references to contact force direction matrix of dynamic system parent
Definition: dynamic_system.cc:861
virtual void updategRef(fmatvec::Vec &gParent)
references to relative distances of dynamic system parent
Definition: dynamic_system.cc:819
virtual void updategdRef(fmatvec::Vec &gdParent)
references to relative velocities of dynamic system parent
Definition: dynamic_system.cc:826
Element(const std::string &name)
constructor
Definition: element.cc:63
InitStage
The stages of the initialization.
Definition: element.h:63
std::string name
name of element
Definition: element.h:261
Definition: environment.h:37
Definition: function.h:53
class for tree-structured mechanical systems with recursive and flat memory mechanism
Definition: graph.h:36
group ingredients do not depend on each other
Definition: group.h:35
Definition: environment.h:59
Definition: mbsim_event.h:79
namespace MBSim
Definition: bilateral_constraint.cc:30
Definition: dynamic_system_solver.h:38