1 #ifndef OPENSIM_MANAGER_H_ 2 #define OPENSIM_MANAGER_H_ 3 /* -------------------------------------------------------------------------- * 4 * OpenSim: Manager.h * 5 * -------------------------------------------------------------------------- * 6 * The OpenSim API is a toolkit for musculoskeletal modeling and simulation. * 7 * See http://opensim.stanford.edu and the NOTICE file for more information. * 8 * OpenSim is developed at Stanford University and supported by the US * 9 * National Institutes of Health (U54 GM072970, R24 HD065690) and by DARPA * 10 * through the Warrior Web program. * 11 * * 12 * Copyright (c) 2005-2017 Stanford University and the Authors * 13 * Author(s): Frank C. Anderson * 14 * * 15 * Licensed under the Apache License, Version 2.0 (the "License"); you may * 16 * not use this file except in compliance with the License. You may obtain a * 17 * copy of the License at http://www.apache.org/licenses/LICENSE-2.0. * 18 * * 19 * Unless required by applicable law or agreed to in writing, software * 20 * distributed under the License is distributed on an "AS IS" BASIS, * 21 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. * 22 * See the License for the specific language governing permissions and * 23 * limitations under the License. * 24 * -------------------------------------------------------------------------- */ 25 26 /* Note: This code was originally developed by Realistic Dynamics Inc. 27 * Author: Frank C. Anderson 28 */ 29 30 31 // INCLUDES 32 #include <OpenSim/Common/Array.h> 33 #include "OpenSim/Common/TimeSeriesTable.h" 34 #include <OpenSim/Simulation/osimSimulationDLL.h> 35 #include <SimTKcommon/internal/ReferencePtr.h> 36 37 namespace SimTK { 38 class Integrator; 39 class State; 40 class System; 41 class TimeStepper; 42 } 43 44 namespace OpenSim { 45 46 class Model; 47 class Storage; 48 class ControllerSet; 49 50 //============================================================================= 51 //============================================================================= 52 /** 53 * A class that manages the execution of a simulation. This class uses a 54 * SimTK::Integrator and SimTK::TimeStepper to perform the simulation. By 55 * default, a Runge-Kutta-Merson integrator is used, but can be changed by 56 * using setIntegratorMethod(). 57 * 58 * In order to prevent an inconsistency between the Integrator and TimeStepper, 59 * we only create a TimeStepper once, specifically when we call 60 * initialize(). To ensure this, the Manager will throw 61 * an exception if initialize() is called more than once. Note 62 * that editing the SimTK::State after calling initialize() 63 * will not affect the simulation. 64 * 65 * Note that this interface means that you cannot "reinitialize" a Manager. 66 * If you make changes to the SimTK::State, a new Manager must be created 67 * before integrating again. 68 */ 69 class OSIMSIMULATION_API Manager 70 { 71 72 //============================================================================= 73 // DATA 74 //============================================================================= 75 private: 76 /** Simulation session name. */ 77 std::string _sessionName; 78 /** Model for which the simulation is performed. */ 79 SimTK::ReferencePtr<Model> _model; 80 81 /** Integrator. */ 82 std::unique_ptr<SimTK::Integrator> _integ; 83 84 /** TimeStepper */ 85 std::unique_ptr<SimTK::TimeStepper> _timeStepper; 86 87 /** Storage for the states. */ 88 std::unique_ptr<Storage> _stateStore; 89 90 /** Flag for signaling a desired halt. */ 91 bool _halt; 92 93 /** Flag to indicate whether or not specified integration time steps 94 should be used. The specified integration time steps are held in _tVec. 95 If _tVec does not contain time steps appropriate for the integration, 96 an exception is thrown. */ 97 bool _specifiedDT; 98 /** Flag to indicate whether or not constant (fixed) integration time 99 steps should be used. The constant integration time step is set using 100 setDT(). */ 101 bool _constantDT; 102 /** Constant integration time step. */ 103 double _dt; 104 /** Vector of integration time steps. */ 105 Array<double> _tArray; 106 /** Vector of integration time step deltas. */ 107 Array<double> _dtArray; 108 109 /** Name to be shown by the UI */ 110 static std::string _displayName; 111 112 /** flag indicating if manager should call Analyses after each step */ 113 bool _performAnalyses; 114 115 /** flag indicating if manager should write to storage each step */ 116 bool _writeToStorage; 117 118 /** controllerSet used for the integration */ 119 ControllerSet* _controllerSet; 120 121 122 //============================================================================= 123 // METHODS 124 //============================================================================= 125 public: 126 /** Constructor that takes a model only and internally uses a 127 * SimTK::RungeKuttaMersonIntegrator with default settings (accuracy, 128 * constraint tolerance, etc.). */ 129 Manager(Model& model); 130 131 /** Convenience constructor for creating and initializing a Manager (by 132 * calling initialize()). 133 * Do not use this constructor if you intend to change integrator settings; 134 * changes to the integrator may not take effect after initializing. */ 135 Manager(Model& model, const SimTK::State& state); 136 137 /** <b>(Deprecated)</b> A Constructor that does not take a model or 138 * controllerSet. This constructor also does not set an integrator; you 139 * must call setIntegrator() on your own. You should use one of the other 140 * constructors. */ 141 DEPRECATED_14("There will be no replacement for this constructor.") 142 Manager(); 143 144 // This class would not behave properly if copied (we would need to write a 145 // complex custom copy constructor, etc.), so don't allow copies. 146 Manager(const Manager&) = delete; 147 void operator=(const Manager&) = delete; 148 149 private: 150 void setNull(); 151 bool constructStorage(); 152 //-------------------------------------------------------------------------- 153 // GET AND SET 154 //-------------------------------------------------------------------------- 155 public: 156 void setSessionName(const std::string &name); 157 void setModel(Model& model); 158 const std::string& getSessionName() const; 159 const std::string& toString() const; 160 setPerformAnalyses(bool performAnalyses)161 void setPerformAnalyses(bool performAnalyses) 162 { _performAnalyses = performAnalyses; } setWriteToStorage(bool writeToStorage)163 void setWriteToStorage(bool writeToStorage) 164 { _writeToStorage = writeToStorage; } 165 166 /** @name Configure the Integrator 167 * @note Call these functions before calling `Manager::initialize()`. 168 * @{ */ 169 170 // Integrator 171 /** Supported integrator methods. For MATLAB, int's must be used rather 172 than enum's (see example in setIntegratorMethod(IntegratorMethod)). */ 173 enum class IntegratorMethod { 174 ExplicitEuler = 0, ///< 0 : For details, see SimTK::ExplicitEulerIntegrator. 175 RungeKutta2 = 1, ///< 1 : For details, see SimTK::RungeKutta2Integrator. 176 RungeKutta3 = 2, ///< 2 : For details, see SimTK::RungeKutta3Integrator. 177 RungeKuttaFeldberg = 3, ///< 3 : For details, see SimTK::RungeKuttaFeldbergIntegrator. 178 RungeKuttaMerson = 4, ///< 4 : For details, see SimTK::RungeKuttaMersonIntegrator. 179 SemiExplicitEuler2 = 5, ///< 5 : For details, see SimTK::SemiExplicitEuler2Integrator. 180 Verlet = 6 ///< 6 : For details, see SimTK::VerletIntegrator. 181 182 // Not included 183 //CPodes, stochastic segfaults when destructed via unique_ptr::reset() 184 //SemiExplicitEuler, no error ctrl, requires fixed stepSize arg on construction 185 }; 186 187 /** Sets the integrator method used via IntegratorMethod enum. The 188 * integrator will be set to its default options, even if the caller 189 * requests the same integrator method. Note that this function must 190 * be called before `Manager::initialize()`. 191 192 <b>C++ example</b> 193 \code{.cpp} 194 auto manager = Manager(model); 195 manager.setIntegratorMethod(Manager::IntegratorMethod::SemiExplicitEuler2); 196 \endcode 197 198 <b>Python example</b> 199 \code{.py} 200 import opensim 201 manager = opensim.Manager(model) 202 manager.setIntegratorMethod(opensim.Manager.IntegratorMethod_SemiExplicitEuler2) 203 \endcode 204 205 <b>MATLAB example</b> 206 \code{.m} 207 import org.opensim.modeling.* 208 manager = Manager(model); 209 manager.setIntegratorMethod(5); 210 \endcode 211 */ 212 void setIntegratorMethod(IntegratorMethod integMethod); 213 214 SimTK::Integrator& getIntegrator() const; 215 216 /** Sets the accuracy of the integrator. 217 * For more details, see `SimTK::Integrator::setAccuracy(SimTK::Real)`. */ 218 void setIntegratorAccuracy(double accuracy); 219 220 /** Sets the minimum step size of the integrator. 221 * For more details, see `SimTK::Integrator::setMinimumStepSize(SimTK::Real)`. */ 222 void setIntegratorMinimumStepSize(double hmin); 223 224 /** Sets the maximum step size of the integrator. 225 * For more details, see `SimTK::Integrator::setMaximumStepSize(SimTK::Real)`. */ 226 void setIntegratorMaximumStepSize(double hmax); 227 228 /** Sets the limit of steps the integrator can take per call of `stepTo()`. 229 * Note that Manager::integrate() calls `stepTo()` for each interval when a fixed 230 * step size is used. 231 * For more details, see SimTK::Integrator::setInternalStepLimit(int). */ 232 void setIntegratorInternalStepLimit(int nSteps); 233 234 //void setIntegratorFixedStepSize(double stepSize); 235 236 /** @} */ 237 238 // SPECIFIED TIME STEP 239 void setUseSpecifiedDT(bool aTrueFalse); 240 bool getUseSpecifiedDT() const; 241 // CONSTANT TIME STEP 242 void setUseConstantDT(bool aTrueFalse); 243 bool getUseConstantDT() const; 244 // DT VECTOR 245 const Array<double>& getDTArray(); 246 void setDTArray(const SimTK::Vector_<double>& aDT, double aTI = 0.0); 247 double getDTArrayDT(int aStep); 248 void printDTArray(const char *aFileName=NULL); 249 // TIME VECTOR 250 const Array<double>& getTimeArray(); 251 double getTimeArrayTime(int aStep); 252 int getTimeArrayStep(double aTime); 253 void printTimeArray(const char *aFileName=NULL); 254 void resetTimeAndDTArrays(double aTime); 255 double getNextTimeArrayTime(double aTime); 256 257 258 259 //-------------------------------------------------------------------------- 260 // EXECUTION 261 //-------------------------------------------------------------------------- 262 /** 263 * Initializes the Manager by creating and initializing the underlying 264 * SimTK::TimeStepper. This must be called before calling 265 * Manager::integrate() Subsequent changes to the State object passed in 266 * here will not affect the simulation. Calling this function multiple 267 * times with the same Manager will trigger an Exception. 268 * 269 * Changes to the integrator (e.g., setIntegratorAccuracy()) after calling 270 * initialize() may not have any effect. 271 */ 272 void initialize(const SimTK::State& s); 273 274 /** 275 * Integrate the equations of motion for the specified model, given the current 276 * state (at which the integration will start) and a finalTime. You must call 277 * Manager::initialize() before calling this function. 278 * 279 * If you must update states or controls in a loop, you must recreate the 280 * manager within the loop (such discontinuous changes are considered "events" 281 * and cannot be handled during integration of the otherwise continuous system). 282 * The proper way to handle this situation is to create a SimTK::EventHandler. 283 * 284 * Example: Integrating from time = 1s to time = 2s 285 * @code 286 * SimTK::State state = model.initSystem(); 287 * Manager manager(model); 288 * state.setTime(1.0); 289 * manager.initialize(state); 290 * state = manager.integrate(2.0); 291 * @endcode 292 * 293 * Example: Integrating from time = 1s to time = 2s using the 294 * convenience constructor 295 * @code 296 * SimTK::State state = model.initSystem(); 297 * state.setTime(1.0); 298 * Manager manager(model, state); 299 * state = manager.integrate(2.0); 300 * @endcode 301 * 302 * Example: Integrate from time = 0s to time = 10s, in 2s increments 303 * @code 304 * dTime = 2.0; 305 * finalTime = 10.0; 306 * int n = int(round(finalTime/dTime)); 307 * state.setTime(0.0); 308 * manager.initialize(state); 309 * for (int i = 1; i <= n; ++i) { 310 * state = manager.integrate(i*dTime); 311 * } 312 * @endcode 313 * 314 * Example: Integrate from time = 0s to time = 10s, updating the state 315 * (e.g., the model's first coordinate value) every 2s 316 * @code 317 * dTime = 2.0; 318 * finalTime = 10.0; 319 * int n = int(round(finalTime/dTime)); 320 * state.setTime(0.0); 321 * for (int i = 0; i < n; ++i) { 322 * model.getCoordinateSet().get(0).setValue(state, 0.1*i); 323 * Manager manager(model); 324 * state.setTime(i*dTime); 325 * manager.initialize(state); 326 * state = manager.integrate((i+1)*dTime); 327 * } 328 * @endcode 329 * 330 */ 331 const SimTK::State& integrate(double finalTime); 332 333 /** Get the current State from the Integrator associated with this 334 * Manager. */ 335 const SimTK::State& getState() const; 336 337 double getFixedStepSize(int tArrayStep) const; 338 339 // STATE STORAGE 340 bool hasStateStorage() const; 341 /** Set the Storage object to be used for storing states. The Manager takes 342 ownership of the passed-in Storage. */ 343 void setStateStorage(Storage& aStorage); 344 Storage& getStateStorage() const; 345 TimeSeriesTable getStatesTable() const; 346 347 //-------------------------------------------------------------------------- 348 // INTERRUPT 349 //-------------------------------------------------------------------------- 350 void halt(); 351 void clearHalt(); 352 bool checkHalt(); 353 354 private: 355 356 // Handles common tasks of some of the other constructors. 357 Manager(Model& model, bool dummyVar); 358 359 // Helper functions during initialization of integration 360 void initializeStorageAndAnalyses(const SimTK::State& s); 361 362 // Helper to record state and analysis values at integration steps. 363 // step = 0 is the beginning, step = -1 used to denote the end/final step 364 void record(const SimTK::State& s, const int& step); 365 366 //============================================================================= 367 }; // END of class Manager 368 369 }; //namespace 370 //============================================================================= 371 //============================================================================= 372 373 #endif // OPENSIM_MANAGER_H_ 374 375