This page explains some of the specific code related to this project. Some of the code is similar to the example given in the Creating an Optimization example that is distributed with OpenSim. Please consult that page for support on generally setting up an optimization. This page covers the following topics:
The global variables are modified from the optimization example.
A forceThreshold is added for the event handler. This is the total amount of vertical force for the contact points when the event handler will trigger. A settleTime is added so that the feet can come off initially without the event handler terminating the simulation early. Four files are also added to help with output.
The main function handles reading in the model file, adding all of the controllers, and calling on the optimizer. The first part of the code reads in the model file and initializes the controller parameters.
The line to setSerializeAllDefaults(true) ensures that all parts of the model are printed out, but is not necessary after OpenSim 3.1. The vector of controllerParameters represents the initial time of the "bang" of the control and the duration of the "bang". Thus, the even indices correspond to initial time for each actuator, and the odd indices for the durations. The controllerParameters are initialized in this block of code and can be changed by the user in this .cpp file.
The next part of the code adds controllers to each of the actuators, enforcing symmetry.
The order of the actuators are assumed to be in the same order for the left and right sides, with all of the actuators on one side listed first, then the other side. Symmetry is ensured by assigning the same PiecewiseConstantFunction for the controller added to each side.
The last part is to set the settings of the Optimizer and OptimizerSystem.
The lower_bounds for the durations were set to be small, but not 0, in order to decrease the search space of arbitrarily small excitation patterns that will not affect the force output in the end.
An event handler is used to stop the simulation once the model leaves the ground (i.e. the vertical force at the contact points reaches forceThreshold). At this point, the final height is known by the ballistic trajectory, so we can save time by not simulating time in the air.
The class TerminateSimulation inherits properties from TriggeredEventHandler. The constructor needs to know what stage it needs to reach to interpret the witness function, which is dynamics for this case since we are interested in forces. It also needs to know that it will trigger when the sign change is falling (i.e. goes from positive to negative). The witness function returns a value at each integration time step. We see that before settleTime is reached, a positive number is always returned. After we reach the settleTime, we retrieve the forces on both feet and return back the sum of the vertical ground reaction force minus the forceThreshold. The last part is to fill in handleEvent, which tells the event handler what to do when the event is triggered.
The objective function first updates the control values with the newControllerParameters.
A new PiecewiseConstantFunction is defined using the newControllerParameters passed in by the optimizer. These are then used to replace the prescribed controller for each of the controllers in the model, again enforcing symmetry.
The next step is to initialize the system. We see here, though, that initSystem() cannot be used because the event handler must be added between the steps of buildSystem() and initializeState().
Similar steps as the optimization example are needed to initialize activations, equilibrate the muscles, and setup and run the forward integration.
The last key part is evaluating the objective function to pass to the Optimizer.
We must first realize the system to velocity since we need position and velocity to evaluate the objection function. We get the position and velocity of the system's center of mass (the ground has a mass of 0 in the model) by going to the MultibodySystem. Gravity is found by calling on the model. As with the Arm26 Optimization Example, the negative of the value is taken since optimizers minimize objective functions.