Description
Chrono Simulation System.
This class is the main simulation object used to collect:
- the simulation model itself, through the ChAssembly object;
- the solver and the timestepper;
- the collision system;
- system-wise parameters (time, gravitational acceleration, etc.)
- counters and timers; as well as other secondary features.
This class is in charge of managing dynamic, kinematic and static simulations.
Consult the Simulation system manual page for additional details.
#include <ChSystem.h>
Classes | |
class | CustomCollisionCallback |
Class to be used as a callback interface for user defined actions performed at each collision detection step. More... | |
Public Member Functions | |
ChSystem () | |
Create a physical system. | |
ChSystem (const ChSystem &other) | |
Copy constructor. | |
virtual | ~ChSystem () |
Destructor. | |
virtual ChSystem * | Clone () const =0 |
"Virtual" copy constructor. More... | |
void | SetTimestepperType (ChTimestepper::Type type) |
Set the method for time integration (time stepper type). More... | |
ChTimestepper::Type | GetTimestepperType () const |
Get the current method for time integration (time stepper type). | |
void | SetTimestepper (std::shared_ptr< ChTimestepper > stepper) |
Set the timestepper object to be used for time integration. | |
std::shared_ptr< ChTimestepper > | GetTimestepper () const |
Get the timestepper currently used for time integration. | |
virtual void | SetCollisionSystemType (ChCollisionSystem::Type type) |
Set the collision detection system used by this Chrono system to the specified type. | |
virtual void | SetCollisionSystem (std::shared_ptr< ChCollisionSystem > coll_system) |
Set the collision detection system used by this Chrono system. | |
std::shared_ptr< ChCollisionSystem > | GetCollisionSystem () const |
Access the underlying collision system. More... | |
virtual void | SetMaterialCompositionStrategy (std::unique_ptr< ChContactMaterialCompositionStrategy > &&strategy) |
Change the default composition laws for contact surface materials (coefficient of friction, cohesion, compliance, etc.) | |
const ChContactMaterialCompositionStrategy & | GetMaterialCompositionStrategy () const |
Accessor for the current composition laws for contact surface material. | |
void | SetMaxPenetrationRecoverySpeed (double value) |
Set the speed limit of exiting from penetration situations (default: 0.6). More... | |
virtual void | SetSolver (std::shared_ptr< ChSolver > newsolver) |
Attach a solver (derived from ChSolver) for use by this system. More... | |
virtual std::shared_ptr< ChSolver > | GetSolver () |
Access the solver currently associated with this system. | |
void | SetSolverType (ChSolver::Type type) |
Choose the solver type, to be used for the simultaneous solution of the constraints in dynamical simulations (as well as in kinematics, statics, etc.) Notes: More... | |
ChSolver::Type | GetSolverType () const |
Gets the current solver type. | |
void | SetSystemDescriptor (std::shared_ptr< ChSystemDescriptor > newdescriptor) |
Instead of using the default 'system descriptor', you can create your own custom descriptor (inherited from ChSystemDescriptor) and plug it into the system using this function. | |
std::shared_ptr< ChSystemDescriptor > | GetSystemDescriptor () |
Access directly the 'system descriptor'. | |
void | SetGravitationalAcceleration (const ChVector3d &gacc) |
Set the gravitational acceleration vector. | |
const ChVector3d & | GetGravitationalAcceleration () const |
Get the gravitatoinal acceleration vector. | |
double | GetChTime () const |
Get the simulation time of this system. | |
void | SetChTime (double time) |
Set (overwrite) the simulation time of this system. | |
double | GetStep () const |
Gets the current time step used for integration (dynamic and kinematic simulation). More... | |
size_t | GetNumSteps () const |
Return the total number of time steps taken so far. | |
void | ResetNumSteps () |
Reset to 0 the total number of time steps. | |
int | DoStepDynamics (double step_size) |
Advance the dynamics simulation by a single time step of given length. More... | |
bool | DoFrameDynamics (double frame_time, double step_size) |
Advance the dynamics simulation until the specified frame end time is reached. More... | |
AssemblyAnalysis::ExitFlag | DoStepKinematics (double step_size) |
Advance the kinematics simulation for a single step of given length. | |
AssemblyAnalysis::ExitFlag | DoFrameKinematics (double frame_time, double step_size) |
Advance the kinematics simulation to the specified frame end time. More... | |
bool | DoStaticAnalysis (ChStaticAnalysis &analysis) |
Perform a generic static analysis. | |
bool | DoStaticLinear () |
Solve the position of static equilibrium (and the reactions). More... | |
bool | DoStaticNonlinear (int nsteps=10, bool verbose=false) |
Solve the position of static equilibrium (and the reactions). More... | |
bool | DoStaticNonlinearRheonomic (int max_num_iterations=10, bool verbose=false, std::shared_ptr< ChStaticNonLinearRheonomicAnalysis::IterationCallback > callback=nullptr) |
Solve the position of static equilibrium (and the reactions). More... | |
bool | DoStaticRelaxing (double step_size, int num_iterations=10) |
Find the static equilibrium configuration (and the reactions) starting from the current position. More... | |
unsigned int | GetSolverSolveCount () const |
Return the number of calls to the solver's Solve() function. More... | |
unsigned int | GetSolverSetupCount () const |
Return the number of calls to the solver's Setup() function. More... | |
AssemblyAnalysis::ExitFlag | DoAssembly (int action, int max_num_iterationsNR=6, double abstol_residualNR=1e-10, double reltol_updateNR=1e-6, double abstol_updateNR=1e-6) |
Assemble the system. More... | |
unsigned int | RemoveRedundantConstraints (bool remove_links=false, double qr_tol=1e-6, bool verbose=false) |
Remove redundant constraints through QR decomposition of the constraints Jacobian matrix. More... | |
virtual void | SetNumThreads (int num_threads_chrono, int num_threads_collision=0, int num_threads_eigen=0) |
Set the number of OpenMP threads used by Chrono itself, Eigen, and the collision detection system. More... | |
unsigned int | GetNumThreadsChrono () const |
unsigned int | GetNumThreadsCollision () const |
unsigned int | GetNumThreadsEigen () const |
const ChAssembly & | GetAssembly () const |
Get the underlying assembly containing all physics items. | |
virtual void | AddBody (std::shared_ptr< ChBody > body) |
Attach a body to the underlying assembly. | |
virtual void | AddShaft (std::shared_ptr< ChShaft > shaft) |
Attach a shaft to the underlying assembly. | |
virtual void | AddLink (std::shared_ptr< ChLinkBase > link) |
Attach a link to the underlying assembly. | |
virtual void | AddMesh (std::shared_ptr< fea::ChMesh > mesh) |
Attach a mesh to the underlying assembly. | |
virtual void | AddOtherPhysicsItem (std::shared_ptr< ChPhysicsItem > item) |
Attach a ChPhysicsItem object that is not a body, link, or mesh. | |
void | Add (std::shared_ptr< ChPhysicsItem > item) |
Attach an arbitrary ChPhysicsItem (e.g. More... | |
void | AddBatch (std::shared_ptr< ChPhysicsItem > item) |
Items added in this way are added like in the Add() method, but not instantly, they are simply queued in a batch of 'to add' items, that are added automatically at the first Setup() call. More... | |
void | FlushBatch () |
If some items are queued for addition in the assembly, using AddBatch(), this will effectively add them and clean the batch. More... | |
virtual void | RemoveBody (std::shared_ptr< ChBody > body) |
Remove a body from this assembly. | |
virtual void | RemoveShaft (std::shared_ptr< ChShaft > shaft) |
Remove a shaft from this assembly. | |
virtual void | RemoveLink (std::shared_ptr< ChLinkBase > link) |
Remove a link from this assembly. | |
virtual void | RemoveMesh (std::shared_ptr< fea::ChMesh > mesh) |
Remove a mesh from the assembly. | |
virtual void | RemoveOtherPhysicsItem (std::shared_ptr< ChPhysicsItem > item) |
Remove a ChPhysicsItem object that is not a body or a link. | |
void | Remove (std::shared_ptr< ChPhysicsItem > item) |
Remove arbitrary ChPhysicsItem that was added to the underlying assembly. | |
void | RemoveAllBodies () |
Remove all bodies from the underlying assembly. | |
void | RemoveAllShafts () |
Remove all shafts from the underlying assembly. | |
void | RemoveAllLinks () |
Remove all links from the underlying assembly. | |
void | RemoveAllMeshes () |
Remove all meshes from the underlying assembly. | |
void | RemoveAllOtherPhysicsItems () |
Remove all physics items not in the body, link, or mesh lists. | |
const std::vector< std::shared_ptr< ChBody > > & | GetBodies () const |
Get the list of bodies. | |
const std::vector< std::shared_ptr< ChShaft > > & | GetShafts () const |
Get the list of shafts. | |
const std::vector< std::shared_ptr< ChLinkBase > > & | GetLinks () const |
Get the list of links. | |
const std::vector< std::shared_ptr< fea::ChMesh > > & | GetMeshes () const |
Get the list of meshes. | |
const std::vector< std::shared_ptr< ChPhysicsItem > > & | GetOtherPhysicsItems () const |
Get the list of physics items that are not in the body or link lists. | |
std::shared_ptr< ChBody > | SearchBody (const std::string &name) const |
Search a body by its name. | |
std::shared_ptr< ChBody > | SearchBodyID (int id) const |
Search a body by its ID. | |
std::shared_ptr< ChShaft > | SearchShaft (const std::string &name) const |
Search a shaft by its name. | |
std::shared_ptr< ChLinkBase > | SearchLink (const std::string &name) const |
Search a link by its name. | |
std::shared_ptr< fea::ChMesh > | SearchMesh (const std::string &name) const |
Search a mesh by its name. | |
std::shared_ptr< ChPhysicsItem > | SearchOtherPhysicsItem (const std::string &name) const |
Search from other ChPhysics items (not bodies, links, or meshes) by name. | |
std::shared_ptr< ChMarker > | SearchMarker (const std::string &name) const |
Search a marker by its name. | |
std::shared_ptr< ChMarker > | SearchMarker (int id) const |
Search a marker by its unique ID. | |
std::shared_ptr< ChPhysicsItem > | Search (const std::string &name) const |
Search an item (body, link or other ChPhysics items) by name. | |
virtual unsigned int | GetNumBodies () const |
Get the total number of bodies added to the system, including fixed and sleeping bodies. | |
virtual unsigned int | GetNumBodiesActive () const |
Get the number of active bodies, excluding sleeping or fixed. | |
virtual unsigned int | GetNumBodiesSleeping () const |
Get the number of sleeping bodies. | |
virtual unsigned int | GetNumBodiesFixed () const |
Get the number of bodies fixed to ground. | |
virtual unsigned int | GetNumShafts () const |
Get the number of shafts. | |
virtual unsigned int | GetNumShaftsSleeping () const |
Get the number of shafts that are in sleeping mode (excluding fixed shafts). | |
virtual unsigned int | GetNumShaftsFixed () const |
Get the number of shafts that are fixed to ground. | |
virtual unsigned int | GetNumShaftsTotal () const |
Get the total number of shafts added to the assembly, including the grounded and sleeping shafts. | |
virtual unsigned int | GetNumLinks () const |
Get the number of links (including non active). | |
virtual unsigned int | GetNumLinksActive () const |
Get the number of active links. | |
virtual unsigned int | GetNumMeshes () const |
Get the number of meshes. | |
virtual unsigned int | GetNumOtherPhysicsItems () const |
Get the number of other physics items (including non active). | |
virtual unsigned int | GetNumOtherPhysicsItemsActive () const |
Get the number of other active physics items. | |
void | ShowHierarchy (std::ostream &m_file, int level=0) const |
Write the hierarchy of contained bodies, markers, etc. More... | |
void | Clear () |
Removes all bodies/marker/forces/links/contacts, also resets timers and events. | |
virtual ChContactMethod | GetContactMethod () const =0 |
Return the contact method supported by this system. More... | |
virtual void | CustomEndOfStep () |
Executes custom processing at the end of step. More... | |
double | ComputeCollisions () |
Perform the collision detection. More... | |
void | RegisterCustomCollisionCallback (std::shared_ptr< CustomCollisionCallback > callback) |
Specify a callback object to be invoked at each collision detection step. More... | |
void | UnregisterCustomCollisionCallback (std::shared_ptr< CustomCollisionCallback > callback) |
Remove the given collision callback from this system. | |
virtual void | SetContactContainer (std::shared_ptr< ChContactContainer > container) |
Change the underlying contact container. More... | |
std::shared_ptr< ChContactContainer > | GetContactContainer () const |
Access the underlying contact container. More... | |
void | SetSleepingAllowed (bool ms) |
Turn on this feature to let the system put to sleep the bodies whose motion has almost come to a rest. More... | |
bool | IsSleepingAllowed () const |
Tell if the system will put to sleep the bodies whose motion has almost come to a rest. | |
ChVisualSystem * | GetVisualSystem () const |
Get the visual system to which this ChSystem is attached (if any). | |
virtual unsigned int | GetNumContacts () |
Gets the number of contacts. | |
virtual double | GetTimerStep () const |
Return the time (in seconds) spent for computing the time step. | |
virtual double | GetTimerAdvance () const |
Return the time (in seconds) for time integration, within the time step. | |
virtual double | GetTimerLSsolve () const |
Return the time (in seconds) for the solver, within the time step. More... | |
virtual double | GetTimerLSsetup () const |
Return the time (in seconds) for the solver Setup phase, within the time step. | |
virtual double | GetTimerJacobian () const |
Return the time (in seconds) for calculating/loading Jacobian information, within the time step. | |
virtual double | GetTimerCollision () const |
Return the time (in seconds) for runnning the collision detection step, within the time step. | |
virtual double | GetTimerSetup () const |
Return the time (in seconds) for system setup, within the time step. | |
virtual double | GetTimerUpdate () const |
Return the time (in seconds) for updating auxiliary data, within the time step. | |
double | GetTimerCollisionBroad () const |
Return the time (in seconds) for broadphase collision detection, within the time step. | |
double | GetTimerCollisionNarrow () const |
Return the time (in seconds) for narrowphase collision detection, within the time step. | |
double | GetRTF () const |
Get current estimated RTF (real time factor). More... | |
void | SetRTF (double rtf) |
Set (overwrite) the RTF value for this system (if calculated externally). | |
void | ResetTimers () |
Resets the timers. | |
void | EnableSolverMatrixWrite (bool val, const std::string &out_dir=".") |
DEBUGGING. More... | |
bool | IsSolverMatrixWriteEnabled () const |
Return true if solver matrices are being written to disk. | |
void | WriteSystemMatrices (bool save_M, bool save_K, bool save_R, bool save_Cq, const std::string &path, bool one_indexed=true) |
Write the mass (M), damping (R), stiffness (K), and constraint Jacobian (Cq) matrices at current configuration. More... | |
void | GetMassMatrix (ChSparseMatrix &M) |
Compute the system-level mass matrix and load in the provided sparse matrix. | |
void | GetStiffnessMatrix (ChSparseMatrix &K) |
Compute the system-level stiffness matrix and load in the provided sparse matrix. More... | |
void | GetDampingMatrix (ChSparseMatrix &R) |
Compute the system-level damping matrix and load in the provided sparse matrix. More... | |
void | GetConstraintJacobianMatrix (ChSparseMatrix &Cq) |
Compute the system-level constraint Jacobian matrix and load in the provided sparse matrix. More... | |
virtual void | ArchiveOut (ChArchiveOut &archive_out) |
Method to allow serialization of transient data to archives. | |
virtual void | ArchiveIn (ChArchiveIn &archive_in) |
Method to allow deserialization of transient data from archives. | |
virtual void | Setup () |
Counts the number of bodies and links. More... | |
void | Update (double mytime, bool update_assets=true) |
Updates all the auxiliary data and children of bodies, forces, links, given their current state. | |
void | Update (bool update_assets=true) |
Updates all the auxiliary data and children of bodies, forces, links, given their current state. | |
void | ForceUpdate () |
In normal usage, no system update is necessary at the beginning of a new dynamics step (since an update is performed at the end of a step). More... | |
void | IntToDescriptor (const unsigned int off_v, const ChStateDelta &v, const ChVectorDynamic<> &R, const unsigned int off_L, const ChVectorDynamic<> &L, const ChVectorDynamic<> &Qc) |
void | IntFromDescriptor (const unsigned int off_v, ChStateDelta &v, const unsigned int off_L, ChVectorDynamic<> &L) |
void | InjectVariables (ChSystemDescriptor &sys_descriptor) |
Register with the given system descriptor all ChVariable objects associated with items in the system. | |
void | InjectConstraints (ChSystemDescriptor &sys_descriptor) |
Register with the given system descriptor any ChConstraint objects associated with items in the system. | |
void | LoadConstraintJacobians () |
Compute and load current Jacobians in encapsulated ChConstraint objects. | |
void | InjectKRMMatrices (ChSystemDescriptor &sys_descriptor) |
Register with the given system descriptor any ChKRMBlock objects associated with items in the system. | |
void | LoadKRMMatrices (double Kfactor, double Rfactor, double Mfactor) |
Compute and load current stiffnes (K), damping (R), and mass (M) matrices in encapsulated ChKRMBlock objects. More... | |
void | VariablesFbReset () |
void | VariablesFbLoadForces (double factor=1) |
void | VariablesQbLoadSpeed () |
void | VariablesFbIncrementMq () |
void | VariablesQbSetSpeed (double step=0) |
void | VariablesQbIncrementPosition (double step) |
void | ConstraintsBiReset () |
void | ConstraintsBiLoad_C (double factor=1, double recovery_clamp=0.1, bool do_clamp=false) |
void | ConstraintsBiLoad_Ct (double factor=1) |
void | ConstraintsBiLoad_Qc (double factor=1) |
void | ConstraintsFbLoadForces (double factor=1) |
void | ConstraintsFetch_react (double factor=1) |
virtual unsigned int | GetNumCoordsPosLevel () override |
Get the number of coordinates at the position level. More... | |
virtual unsigned int | GetNumCoordsVelLevel () override |
Get the number of coordinates at the velocity level. More... | |
virtual unsigned int | GetNumConstraints () override |
Get the number of scalar constraints in the system. | |
virtual unsigned int | GetNumConstraintsBilateral () |
Get the number of bilateral scalar constraints. | |
virtual unsigned int | GetNumConstraintsUnilateral () |
Get the number of unilateral scalar constraints. | |
virtual void | StateGather (ChState &x, ChStateDelta &v, double &T) override |
From system to state y={x,v}. | |
virtual void | StateScatter (const ChState &x, const ChStateDelta &v, const double T, bool full_update) override |
From state Y={x,v} to system. This also triggers an update operation. | |
virtual void | StateGatherAcceleration (ChStateDelta &a) override |
From system to state derivative (acceleration), some timesteppers might need last computed accel. | |
virtual void | StateScatterAcceleration (const ChStateDelta &a) override |
From state derivative (acceleration) to system, sometimes might be needed. | |
virtual void | StateGatherReactions (ChVectorDynamic<> &L) override |
From system to reaction forces (last computed) - some timestepper might need this. | |
virtual void | StateScatterReactions (const ChVectorDynamic<> &L) override |
From reaction forces to system, ex. store last computed reactions in ChLink objects for plotting etc. | |
virtual void | StateIncrementX (ChState &x_new, const ChState &x, const ChStateDelta &Dx) override |
Perform x_new = x + dx, for x in Y = {x, dx/dt}. More... | |
virtual bool | StateSolveCorrection (ChStateDelta &Dv, ChVectorDynamic<> &DL, const ChVectorDynamic<> &R, const ChVectorDynamic<> &Qc, const double c_a, const double c_v, const double c_x, const ChState &x, const ChStateDelta &v, const double T, bool force_state_scatter, bool full_update, bool force_setup) override |
Assuming a DAE of the form. More... | |
virtual void | LoadResidual_F (ChVectorDynamic<> &R, const double c) override |
Increment a vector R with the term c*F: R += c*F. More... | |
virtual void | LoadResidual_Mv (ChVectorDynamic<> &R, const ChVectorDynamic<> &w, const double c) override |
Increment a vector R with a term that has M multiplied a given vector w: R += c*M*w. More... | |
virtual void | LoadLumpedMass_Md (ChVectorDynamic<> &Md, double &err, const double c) override |
Adds the lumped mass to a Md vector, representing a mass diagonal matrix. More... | |
virtual void | LoadResidual_CqL (ChVectorDynamic<> &R, const ChVectorDynamic<> &L, const double c) override |
Increment a vectorR with the term Cq'*L: R += c*Cq'*L. More... | |
virtual void | LoadConstraint_C (ChVectorDynamic<> &Qc, const double c, const bool do_clamp=false, const double clamp=1e30) override |
Increment a vector Qc with the term C: Qc += c*C. More... | |
virtual void | LoadConstraint_Ct (ChVectorDynamic<> &Qc, const double c) override |
Increment a vector Qc with the term Ct = partial derivative dC/dt: Qc += c*Ct. More... | |
Public Member Functions inherited from chrono::ChIntegrableIIorder | |
virtual void | StateSetup (ChState &x, ChStateDelta &v, ChStateDelta &a) |
Set up the system state with separate II order components x, v, a for y = {x, v} and dy/dt={v, a}. | |
virtual bool | StateSolveA (ChStateDelta &Dvdt, ChVectorDynamic<> &L, const ChState &x, const ChStateDelta &v, const double T, const double dt, bool force_state_scatter, bool full_update, ChLumpingParms *lumping=nullptr) |
Solve for accelerations: a = f(x,v,t) Given current state y={x,v} , computes acceleration a in the state derivative dy/dt={v,a} and lagrangian multipliers L (if any). More... | |
virtual void | StateGather (ChState &y, double &T) override |
Gather system state in specified array. More... | |
virtual void | StateScatter (const ChState &y, const double T, bool full_update) override |
Scatter the states from the provided array to the system. More... | |
virtual void | StateGatherDerivative (ChStateDelta &Dydt) override |
Gather from the system the state derivatives in specified array. More... | |
virtual void | StateScatterDerivative (const ChStateDelta &Dydt) override |
Scatter the state derivatives from the provided array to the system. More... | |
virtual void | StateIncrement (ChState &y_new, const ChState &y, const ChStateDelta &Dy) override |
Increment state array: y_new = y + Dy. More... | |
virtual bool | StateSolve (ChStateDelta &dydt, ChVectorDynamic<> &L, const ChState &y, const double T, const double dt, bool force_state_scatter, bool full_update, ChLumpingParms *lumping=nullptr) override |
Solve for state derivatives: dy/dt = f(y,t). More... | |
virtual bool | StateSolveCorrection (ChStateDelta &Dy, ChVectorDynamic<> &L, const ChVectorDynamic<> &R, const ChVectorDynamic<> &Qc, const double a, const double b, const ChState &y, const double T, const double dt, bool force_state_scatter, bool full_update, bool force_setup) override final |
Override of method for Ist order implicit integrators. More... | |
Public Member Functions inherited from chrono::ChIntegrable | |
virtual unsigned int | GetNumCoordsAccLevel () |
Return the number of coordinates at the acceleration level. | |
virtual void | StateSetup (ChState &y, ChStateDelta &dy) |
Set up the system state. | |
virtual void | LoadResidual_Hv (ChVectorDynamic<> &R, const ChVectorDynamic<> &v, const double c) |
Increment a vector R (usually the residual in a Newton Raphson iteration for solving an implicit integration step) with a term that has H multiplied a given vector w: R += c*H*w. More... | |
Protected Member Functions | |
virtual void | DescriptorPrepareInject (ChSystemDescriptor &sys_descriptor) |
Pushes all ChConstraints and ChVariables contained in links, bodies, etc. into the system descriptor. | |
void | Initialize () |
Initial system setup before analysis. More... | |
virtual ChVector3d | GetBodyAppliedForce (ChBody *body) |
Return the resultant applied force on the specified body. More... | |
virtual ChVector3d | GetBodyAppliedTorque (ChBody *body) |
Return the resultant applied torque on the specified body. More... | |
bool | ManageSleepingBodies () |
Put bodies to sleep if possible. More... | |
virtual bool | AdvanceDynamics () |
Performs a single dynamics simulation step, advancing the system state by the current step size. | |
Protected Attributes | |
ChAssembly | assembly |
underlying mechanical assembly | |
std::shared_ptr< ChContactContainer > | contact_container |
the container of contacts | |
ChVector3d | G_acc |
gravitational acceleration | |
bool | is_initialized |
if false, an initial setup is required (i.e. a call to Initialize) | |
bool | is_updated |
if false, a new update is required (i.e. a call to Update) | |
unsigned int | m_num_coords_pos |
num of scalar coordinates at position level for all active bodies | |
unsigned int | m_num_coords_vel |
num of scalar coordinates at velocity level for all active bodies | |
unsigned int | m_num_constr |
num of scalar constraints (at velocity level) for all active constraints | |
unsigned int | m_num_constr_bil |
num of scalar bilateral active constraints (velocity level) | |
unsigned int | m_num_constr_uni |
num of scalar unilateral active constraints (velocity level) | |
double | ch_time |
simulation time of the system | |
double | step |
time step | |
bool | use_sleeping |
if true, put to sleep objects that come to rest | |
std::shared_ptr< ChSystemDescriptor > | descriptor |
system descriptor | |
std::shared_ptr< ChSolver > | solver |
solver for DVI or DAE problem | |
double | max_penetration_recovery_speed |
limit for speed of penetration recovery (positive) | |
size_t | stepcount |
internal counter for steps | |
unsigned int | setupcount |
number of calls to the solver's Setup() | |
unsigned int | solvecount |
number of StateSolveCorrection (reset to 0 at each timestep of static analysis) | |
bool | write_matrix |
write current system matrix to file(s); for debugging | |
std::string | output_dir |
output directory for writing system matrices | |
unsigned int | ncontacts |
total number of contacts | |
std::shared_ptr< ChCollisionSystem > | collision_system |
collision engine | |
std::vector< std::shared_ptr< CustomCollisionCallback > > | collision_callbacks |
user-defined collision callbacks | |
std::unique_ptr< ChContactMaterialCompositionStrategy > | composition_strategy |
ChVisualSystem * | visual_system |
material composition strategy More... | |
int | nthreads_chrono |
int | nthreads_eigen |
int | nthreads_collision |
ChTimer | timer_step |
timer for integration step | |
ChTimer | timer_advance |
timer for time integration | |
ChTimer | timer_ls_solve |
timer for solver (excluding setup phase) | |
ChTimer | timer_ls_setup |
timer for solver setup | |
ChTimer | timer_jacobian |
timer for computing/loading Jacobian information | |
ChTimer | timer_collision |
timer for collision detection | |
ChTimer | timer_setup |
timer for system setup | |
ChTimer | timer_update |
timer for system update | |
double | m_RTF |
real-time factor (simulation time / simulated time) | |
std::shared_ptr< ChTimestepper > | timestepper |
time-stepper object | |
ChVectorDynamic | applied_forces |
system-wide vector of applied forces (lazy evaluation) | |
bool | applied_forces_current |
indicates if system-wide vector of forces is up-to-date | |
Friends | |
class | ChAssembly |
class | ChBody |
class | fea::ChMesh |
class | ChContactContainerNSC |
class | ChContactContainerSMC |
class | ChVisualSystem |
class | ChCollisionSystem |
class | modal::ChModalAssembly |
Member Function Documentation
◆ Add()
void chrono::ChSystem::Add | ( | std::shared_ptr< ChPhysicsItem > | item | ) |
Attach an arbitrary ChPhysicsItem (e.g.
ChBody, ChParticles, ChLink, etc.) to the assembly. It will take care of adding it to the proper list of bodies, links, meshes, or generic physic item. (i.e. it calls AddBody, AddShaft(), AddLink(), AddMesh(), or AddOtherPhysicsItem()). Note, you cannot call Add() during an Update (i.e. items like particle generators that are already inserted in the assembly cannot call this) because not thread safe; instead, use AddBatch().
◆ AddBatch()
|
inline |
◆ Clone()
|
pure virtual |
"Virtual" copy constructor.
Concrete derived classes must implement this.
Implemented in chrono::ChSystemMulticoreSMC, chrono::ChSystemMulticoreNSC, chrono::ChSystemSMC, and chrono::ChSystemNSC.
◆ ComputeCollisions()
double chrono::ChSystem::ComputeCollisions | ( | ) |
Perform the collision detection.
New contacts are inserted in the ChContactContainer object(s), and old ones are removed. This is mostly called automatically by time integration.
◆ CustomEndOfStep()
|
inlinevirtual |
Executes custom processing at the end of step.
By default it does nothing, but if you inherit a special ChSystem you can implement this.
◆ DoAssembly()
AssemblyAnalysis::ExitFlag chrono::ChSystem::DoAssembly | ( | int | action, |
int | max_num_iterationsNR = 6 , |
||
double | abstol_residualNR = 1e-10 , |
||
double | reltol_updateNR = 1e-6 , |
||
double | abstol_updateNR = 1e-6 |
||
) |
Assemble the system.
The assembling is performed by satisfying constraints at position, velocity, and acceleration levels. Position-level assembling requires Newton-Raphson iterations. Velocity-level assembling is performed by taking a small integration step. Acceleration-level assembling is obtained through finite differentation. Argument 'action' can be one of AssemblyLevel enum values (POSITION, VELOCITY, ACCELERATION, or FULL). These values can also be combined using bit operations. Returns true if the assembling converged, false otherwise (impossible assembly?) The maximum number of iterations and the tolerance refer to the Newton-Raphson iteration for position-level. Iterations and tolerance for the inner linear solver must be set on the solver itself.
◆ DoFrameDynamics()
bool chrono::ChSystem::DoFrameDynamics | ( | double | frame_time, |
double | step_size | ||
) |
Advance the dynamics simulation until the specified frame end time is reached.
Integration proceeds with the specified time step size which may be adjusted to exactly reach the frame time.
◆ DoFrameKinematics()
AssemblyAnalysis::ExitFlag chrono::ChSystem::DoFrameKinematics | ( | double | frame_time, |
double | step_size | ||
) |
Advance the kinematics simulation to the specified frame end time.
A system assembly analysis (inverse kinematics) is performed at step_size intervals. The last step may be adjusted to exactly reach the frame end time.
◆ DoStaticLinear()
bool chrono::ChSystem::DoStaticLinear | ( | ) |
Solve the position of static equilibrium (and the reactions).
This is a one-step only approach that solves the linear equilibrium. Appropriate mostly for FEM problems with small deformations.
◆ DoStaticNonlinear()
bool chrono::ChSystem::DoStaticNonlinear | ( | int | nsteps = 10 , |
bool | verbose = false |
||
) |
Solve the position of static equilibrium (and the reactions).
This function solves the equilibrium for the nonlinear problem (large displacements). This version uses a nonlinear static analysis solver with default parameters.
◆ DoStaticNonlinearRheonomic()
bool chrono::ChSystem::DoStaticNonlinearRheonomic | ( | int | max_num_iterations = 10 , |
bool | verbose = false , |
||
std::shared_ptr< ChStaticNonLinearRheonomicAnalysis::IterationCallback > | callback = nullptr |
||
) |
Solve the position of static equilibrium (and the reactions).
This function solves the equilibrium for the nonlinear problem (large displacements), but differently from DoStaticNonlinear it considers rheonomic constraints (ex. ChLinkMotorRotationSpeed) that can impose steady-state speeds&accelerations to the mechanism, ex. to generate centrifugal forces in turbine blades. This version uses a nonlinear static analysis solver with default parameters.
◆ DoStaticRelaxing()
bool chrono::ChSystem::DoStaticRelaxing | ( | double | step_size, |
int | num_iterations = 10 |
||
) |
Find the static equilibrium configuration (and the reactions) starting from the current position.
Since a truncated iterative method is used, you may need to call this method multiple times in case of large nonlinearities before coming to the precise static solution.
◆ DoStepDynamics()
int chrono::ChSystem::DoStepDynamics | ( | double | step_size | ) |
Advance the dynamics simulation by a single time step of given length.
This function is typically called many times in a loop in order to simulate up to a desired end time.
◆ EnableSolverMatrixWrite()
void chrono::ChSystem::EnableSolverMatrixWrite | ( | bool | val, |
const std::string & | out_dir = "." |
||
) |
DEBUGGING.
Enable/disable debug output of system matrices. Set this to "true" to enable automatic saving of solver matrices at each time step, for debugging purposes. Matrices will be saved in 'out_dir' (default to the the working directory of the executable). The name pattern is solve_numstep_numsubstep_objectname.dat where 'objectname' can be either:
- Z: the assembled optimization matrix [H, -Cq'; Cq, -E]
- rhs: the assembled right-hand side vector [f; -b]
- x_pre: the position vector before the integration step
- v_pre: the velocity vector before the integration step
- F_pre: unscaled loads before the integration step
- Dv: the state variation
- Dl: the lagrangian multipliers variation if the solver is direct then also the following are output:
- f, b: subparts of 'rhs'
- H, Cq, E: the submatrices of Z as passed to the solver in the problem
| H Cq'|*| q|-| f|=|0| | Cq E | |-l| |-b| |c|
where l \(\in Y, c \in Ny\), normal cone to Y
◆ FlushBatch()
|
inline |
If some items are queued for addition in the assembly, using AddBatch(), this will effectively add them and clean the batch.
Called automatically at each Setup().
◆ ForceUpdate()
void chrono::ChSystem::ForceUpdate | ( | ) |
In normal usage, no system update is necessary at the beginning of a new dynamics step (since an update is performed at the end of a step).
However, this is not the case if external changes to the system are made. Most such changes are discovered automatically (addition/removal of items, input of mesh loads). For special cases, this function allows the user to trigger a system update at the beginning of the step immediately following this call.
◆ GetBodyAppliedForce()
|
protectedvirtual |
Return the resultant applied force on the specified body.
This resultant force includes all external applied loads acting on the body (from gravity, loads, springs, etc). However, this does not include any constraint forces. In particular, contact forces are not included if using the NSC formulation, but are included when using the SMC formulation.
Reimplemented in chrono::ChSystemMulticore.
◆ GetBodyAppliedTorque()
|
protectedvirtual |
Return the resultant applied torque on the specified body.
This resultant torque includes all external applied loads acting on the body (from gravity, loads, springs, etc). However, this does not include any constraint forces. In particular, contact torques are not included if using the NSC formulation, but are included when using the SMC formulation.
Reimplemented in chrono::ChSystemMulticore.
◆ GetCollisionSystem()
|
inline |
Access the underlying collision system.
Usually this is not needed, as the collision system is automatically handled by the ChSystem.
◆ GetConstraintJacobianMatrix()
void chrono::ChSystem::GetConstraintJacobianMatrix | ( | ChSparseMatrix & | Cq | ) |
Compute the system-level constraint Jacobian matrix and load in the provided sparse matrix.
This is the Jacobian Cq=-dC/dq, where C are constraints (the lower left part of the KKT matrix).
◆ GetContactContainer()
|
inline |
Access the underlying contact container.
Usually this is not needed, as the contact container is automatically handled by the ChSystem.
◆ GetContactMethod()
|
pure virtual |
Return the contact method supported by this system.
Contactables (bodies, FEA nodes, FEA traiangles, etc.) added to this system must be compatible.
Implemented in chrono::ChSystemNSC, chrono::ChSystemMulticoreSMC, chrono::ChSystemMulticoreNSC, and chrono::ChSystemSMC.
◆ GetDampingMatrix()
void chrono::ChSystem::GetDampingMatrix | ( | ChSparseMatrix & | R | ) |
Compute the system-level damping matrix and load in the provided sparse matrix.
This is the Jacobian -dF/dv, where F are stiff loads. Note that not all loads provide a Jacobian, as this is optional in their implementation.
◆ GetNumCoordsPosLevel()
|
inlineoverridevirtual |
Get the number of coordinates at the position level.
Might differ from coordinates at the velocity level if quaternions are used for rotations.
Implements chrono::ChIntegrable.
◆ GetNumCoordsVelLevel()
|
inlineoverridevirtual |
Get the number of coordinates at the velocity level.
Might differ from coordinates at the position level if quaternions are used for rotations.
Implements chrono::ChIntegrable.
◆ GetRTF()
|
inline |
Get current estimated RTF (real time factor).
This represents the real time factor for advancing the dynamic state of the system only and as such does not take into account any other operations performed during a step (e.g., run-time visualization). During each call to DoStepDynamics(), this value is calculated as T/step_size, where T includes the time spent in system setup, collision detection, and integration.
◆ GetSolverSetupCount()
|
inline |
Return the number of calls to the solver's Setup() function.
This counter is reset at each timestep.
◆ GetSolverSolveCount()
|
inline |
Return the number of calls to the solver's Solve() function.
This counter is reset at each timestep.
◆ GetStep()
|
inline |
Gets the current time step used for integration (dynamic and kinematic simulation).
The value is set automatically when a dynamic or kinematic simulation is run.
◆ GetStiffnessMatrix()
void chrono::ChSystem::GetStiffnessMatrix | ( | ChSparseMatrix & | K | ) |
Compute the system-level stiffness matrix and load in the provided sparse matrix.
This is the Jacobian -dF/dq, where F are stiff loads. Note that not all loads provide a jacobian, as this is optional in their implementation.
◆ GetTimerLSsolve()
|
inlinevirtual |
Return the time (in seconds) for the solver, within the time step.
Note that this time excludes any calls to the solver's Setup function.
Reimplemented in chrono::ChSystemMulticore.
◆ Initialize()
|
protected |
Initial system setup before analysis.
This function performs an initial system setup, once system construction is completed and before an analysis. This function also initializes the collision system (if any), as well as any visualization system to which this Chrono system was attached. The initialization function is called automatically before starting any type of analysis.
◆ LoadConstraint_C()
|
overridevirtual |
Increment a vector Qc with the term C: Qc += c*C.
- Parameters
-
Qc result: the Qc residual, Qc += c*C c a scaling factor do_clamp enable optional clamping of Qc clamp clamping value
Reimplemented from chrono::ChIntegrableIIorder.
◆ LoadConstraint_Ct()
|
overridevirtual |
Increment a vector Qc with the term Ct = partial derivative dC/dt: Qc += c*Ct.
- Parameters
-
Qc result: the Qc residual, Qc += c*Ct c a scaling factor
Reimplemented from chrono::ChIntegrableIIorder.
◆ LoadKRMMatrices()
void chrono::ChSystem::LoadKRMMatrices | ( | double | Kfactor, |
double | Rfactor, | ||
double | Mfactor | ||
) |
Compute and load current stiffnes (K), damping (R), and mass (M) matrices in encapsulated ChKRMBlock objects.
The resulting KRM blocks represent linear combinations of the K, R, and M matrices, with the specified coefficients Kfactor, Rfactor,and Mfactor, respectively.
◆ LoadLumpedMass_Md()
|
overridevirtual |
Adds the lumped mass to a Md vector, representing a mass diagonal matrix.
Used by lumped explicit integrators. If mass lumping is impossible or approximate, adds scalar error to "error" parameter. Md += c*diag(M) or Md += c*HRZ(M)
- Parameters
-
Md result: Md vector, diagonal of the lumped mass matrix err result: not touched if lumping does not introduce errors c a scaling factor
Reimplemented from chrono::ChIntegrableIIorder.
◆ LoadResidual_CqL()
|
overridevirtual |
Increment a vectorR with the term Cq'*L: R += c*Cq'*L.
- Parameters
-
R result: the R residual, R += c*Cq'*L L the L vector c a scaling factor
Reimplemented from chrono::ChIntegrableIIorder.
◆ LoadResidual_F()
|
overridevirtual |
Increment a vector R with the term c*F: R += c*F.
- Parameters
-
R result: the R residual, R += c*F c a scaling factor
Reimplemented from chrono::ChIntegrableIIorder.
◆ LoadResidual_Mv()
|
overridevirtual |
Increment a vector R with a term that has M multiplied a given vector w: R += c*M*w.
- Parameters
-
R result: the R residual, R += c*M*v w the w vector c a scaling factor
Reimplemented from chrono::ChIntegrableIIorder.
◆ ManageSleepingBodies()
|
protected |
Put bodies to sleep if possible.
Also awakens sleeping bodies, if needed. Returns true if some body changed from sleep to no sleep or viceversa, returns false if nothing changed. In the former case also performs Setup() since the system changed.
If some body still must change from no sleep-> sleep, do it
◆ RegisterCustomCollisionCallback()
void chrono::ChSystem::RegisterCustomCollisionCallback | ( | std::shared_ptr< CustomCollisionCallback > | callback | ) |
Specify a callback object to be invoked at each collision detection step.
Multiple such callback objects can be registered with a system. If present, their OnCustomCollision() method is invoked. Use this if you want that some specific callback function is executed at each collision detection step (ex. all the times that ComputeCollisions() is automatically called by the integration method). For example some other collision engine could add further contacts using this callback.
◆ RemoveRedundantConstraints()
unsigned int chrono::ChSystem::RemoveRedundantConstraints | ( | bool | remove_links = false , |
double | qr_tol = 1e-6 , |
||
bool | verbose = false |
||
) |
Remove redundant constraints through QR decomposition of the constraints Jacobian matrix.
Remove redundant constraints present in ChSystem through QR decomposition of constraints Jacobian matrix.
This function can be used to improve the stability and performance of the system by removing redundant constraints. The function returns the number of redundant constraints that were deactivated/removed. Please consider that some constraints might falsely appear as redundant in a specific configuration, while they might be not in general.
- Parameters
-
remove_links false: redundant links are just deactivated; true: redundant links get removed qr_tol tolerance in QR decomposition to identify linearly dependent constraints verbose set verbose output
◆ SetContactContainer()
|
virtual |
Change the underlying contact container.
The contact container collects information from the underlying collision detection system required for contact force generation. Usually this is not needed, as the contact container is automatically handled by the ChSystem. Make sure the provided contact container is compatible with both the collision detection system and the contact force formulation (NSC or SMC).
Reimplemented in chrono::ChSystemMulticoreSMC, chrono::ChSystemMulticoreNSC, chrono::ChSystemSMC, and chrono::ChSystemNSC.
◆ SetMaxPenetrationRecoverySpeed()
|
inline |
Set the speed limit of exiting from penetration situations (default: 0.6).
Usually set a positive value, (about 0.1...2 m/s, as exiting speed). Used for unilateral constraints with the EULER_IMPLICIT_LINEARIZED time stepper.
◆ SetNumThreads()
|
virtual |
Set the number of OpenMP threads used by Chrono itself, Eigen, and the collision detection system.
num_threads_chrono - used in FEA (parallel evaluation of internal forces and Jacobians) and in SCM deformable terrain calculations. num_threads_collision - used in parallelization of collision detection (if applicable). If passing 0, then num_threads_collision = num_threads_chrono. num_threads_eigen - used in the Eigen sparse direct solvers and a few linear algebra operations. Note that Eigen enables multi-threaded execution only under certain size conditions. See the Eigen documentation. If passing 0, then num_threads_eigen = num_threads_chrono. By default (if this function is not called), the following values are used: num_threads_chrono = 1 num_threads_collision = 1 num_threads_eigen = 1
Note that a derived class may ignore some or all of these settings.
Reimplemented in chrono::ChSystemMulticore.
◆ SetSleepingAllowed()
|
inline |
Turn on this feature to let the system put to sleep the bodies whose motion has almost come to a rest.
This feature will allow faster simulation of large scenarios for real-time purposes, but it will affect the precision! This functionality can be turned off selectively for specific ChBodies.
◆ SetSolver()
|
virtual |
Attach a solver (derived from ChSolver) for use by this system.
- Suggested solver for speed, but lower precision: PSOR
- Suggested solver for higher precision: BARZILAIBORWEIN or APGD
- For problems that involve a stiffness matrix: GMRES, MINRES
◆ SetSolverType()
void chrono::ChSystem::SetSolverType | ( | ChSolver::Type | type | ) |
Choose the solver type, to be used for the simultaneous solution of the constraints in dynamical simulations (as well as in kinematics, statics, etc.) Notes:
- This function is a shortcut, internally equivalent to a call to SetSolver().
- Only a subset of available Chrono solvers can be set through this mechanism.
- Prefer explicitly creating a solver, setting solver parameters, and then attaching the solver with SetSolver.
- Deprecated:
- This function does not support all available Chrono solvers. Prefer using SetSolver.
◆ SetTimestepperType()
void chrono::ChSystem::SetTimestepperType | ( | ChTimestepper::Type | type | ) |
Set the method for time integration (time stepper type).
- Suggested for fast dynamics with hard (NSC) contacts: EULER_IMPLICIT_LINEARIZED
- Suggested for fast dynamics with hard (NSC) contacts and low inter-penetration: EULER_IMPLICIT_PROJECTED
- Suggested for finite element smooth dynamics: HHT, EULER_IMPLICIT_LINEARIZED For full access to a time stepper's settings, use SetTimestepper()
◆ Setup()
|
virtual |
Counts the number of bodies and links.
Computes the offsets of object states in the global state. Assumes that offset_x, offset_w, and offset_L are already set as starting point for offsetting all the contained sub objects.
Reimplemented in chrono::ChSystemMulticoreSMC, and chrono::ChSystemMulticore.
◆ ShowHierarchy()
|
inline |
Write the hierarchy of contained bodies, markers, etc.
in ASCII readable form, mostly for debugging purposes. Level is the tab spacing at the left.
◆ StateIncrementX()
|
overridevirtual |
Perform x_new = x + dx, for x in Y = {x, dx/dt}.
It takes care of the fact that x has quaternions, dx has angular vel etc. NOTE: the system is not updated automatically after the state increment, so one might need to call StateScatter() if needed.
- Parameters
-
x_new resulting x_new = x + Dx x initial state x Dx state increment Dx
Reimplemented from chrono::ChIntegrableIIorder.
◆ StateSolveCorrection()
|
overridevirtual |
Assuming a DAE of the form.
M*a = F(x,v,t) + Cq'*L C(x,t) = 0
this function computes the solution of the change Du (in a or v or x) for a Newton iteration within an implicit integration scheme.
| Du| = [ G Cq' ]^-1 * | R | |-DL| [ Cq 0 ] |-Qc|
for residual R and G = [ c_a*M + c_v*dF/dv + c_x*dF/dx ].
This function returns true if successful and false otherwise.
- Parameters
-
Dv result: computed Dv DL result: computed lagrangian multipliers. Note the sign in system above. R the R residual Qc the Qc residual. Note the sign in system above. c_a the factor in c_a*M c_v the factor in c_v*dF/dv c_x the factor in c_x*dF/dv x current state, x part v current state, v part T current time T force_state_scatter if false, x and v are not scattered to the system full_update if true, perform a full update during scatter force_setup if true, call the solver's Setup() function
Reimplemented from chrono::ChIntegrableIIorder.
◆ WriteSystemMatrices()
void chrono::ChSystem::WriteSystemMatrices | ( | bool | save_M, |
bool | save_K, | ||
bool | save_R, | ||
bool | save_Cq, | ||
const std::string & | path, | ||
bool | one_indexed = true |
||
) |
Write the mass (M), damping (R), stiffness (K), and constraint Jacobian (Cq) matrices at current configuration.
These can be used for linearized motion, modal analysis, buckling analysis, etc. The sparse matrices are saved in COO format in 'path' folder according to the naming: solve_M.dat solve_K.dat solve_R.dat, and solve_Cq.dat. By default, uses 1-based indices (as in Matlab).
Member Data Documentation
◆ visual_system
|
protected |
material composition strategy
run-time visualization engine
The documentation for this class was generated from the following files:
- /builds/uwsbel/chrono/src/chrono/physics/ChSystem.h
- /builds/uwsbel/chrono/src/chrono/physics/ChSystem.cpp