ArtiSynth supports an inverse simulation capability that allows the computation of the excitation signals (for muscles and other actuators) needed to track target trajectories for a set of source components within the model (typically points or frames). This is achieved using a specialized controller component known as a TrackingController, defined in the package
artisynth.core.inverse
An application model creates the tracking controller and configures it with the source components, the exciter components needed to effect the tracking, and probes or controllers to specify the target trajectories. Then at each simulation time step the controller computes the excitations (i.e., signals to the exciter components) that allow the trajectories to be followed as closely as possible, in a manner conceptually similar to computed muscle control in OpenSim [7].
It is currently recommended to run inverse simulation with hybrid solving disabled. Hybrid solving combines direct and iterative solution techniques to speed up sparse matrix solves within ArtiSynth’s integrators. However, it can occasionally cause stability issues when used with inverse control. Hybrid solves can be disabled in several ways:
- 1.
By setting hybridSolvesEnabled to false under “Settings > Simulation ...” in the GUI; see “Simulation” under “Settings and Preferences” in the ArtiSynth User Interface Guide. This change can be made to persist across ArtiSynth restarts.
- 2.
By calling the static method setHybridSolvesEnabled() of MechSystemSolver in the model’s build() method:
MechSystemSolver.setHybridSolvesEnabled(false);
Let , and be composite vectors describing the positions, velocities and forces of all the components in the application model, and let be a vector of excitation signals for all the excitation components available to the controller. can be decomposed into
(10.1) |
where are the passive forces corresponding to zero excitation, and are the active forces. During simulation, the controller computes values for to best track the target trajectories.
In the case of motion tracking, let and be the subvectors of and containing the positions and velocities of the source components, and let and be the corresponding target trajectories. At the beginning of each time step, the controller compares the current source velocity state with the desired target state and uses this to determine a desired source velocity for the next time step; this computation is done by the motion tracker (Figure 10.1, Section 10.1.2). An excitation generator then computes in order to try and make match (Section 10.1.3).
The excitation generator accepts a desired velocity instead of a desired acceleration because ArtiSynth’s physics engine computes velocity changes directly.
As mentioned above, the motion tracker computes a desired source velocity for the next time step based on the current source state and desired target state . This can be done in two ways: chase control (the default), and PD control.
Chase control simply sets to , plus an additional velocity that tries to reduce the position error over a specified time interval , called the chase time:
(10.2) |
In general, should be greater than or equal to the simulation step size . If it is greater, then will tend to lag behind , but this will also reduce the potential for overshoot due to system nonlinearities. Conversely, if if less than , then is much more likely to overshoot. The default value of is 0.01.
PD control computes a desired source acceleration based on
(10.3) |
and then integrates this to determine :
(10.4) |
where is the simulation step size. PD control offers a greater ability to adjust the tracking behavior than chase control, but it is often necessary to tune the gain parameters and . One rule of thumb is to set their initial values such that (10.2) and (10.4) are equal, which leads to
The default values for and are and , corresponding to and both equaling . Lowering the value of will reduce overshoot but increase tracking lag.
PD control is enabled and adjusted by setting the usePDControl, Kp and Kd properties of the motion target term to true (Section 10.3.3).
Given , the excitation generator computes to try and ensure that the velocity at the end of the subsequent time step, , satisfies
(10.5) |
This is accomplished using a quadratic program.
First, for a broad class of problems, is linearly related to , so that (10.1) simplifies to
(10.6) |
where is an excitation response matrix. Given (10.6), it is possible to show (section “Inverse modeling” in [11]) that
where is the velocity with zero excitations and is a motion excitation response matrix. To best follow the trajectory, we can compute to minimize the quadratic cost function
(10.7) |
If we add in the constraint that excitations lie in the range , the problem takes the form of a quadratic program (QP)
(10.8) |
In order to prioritize some target terms over others, can be modified to include weights, according to
(10.9) |
where is a diagonal weighting matrix. To handle excitation redundancy, where the size of exceeds the size of , we can add a regularization term , where is a diagonal weight matrix that can be used to adjust the relative importance of specific excitations:
(10.10) |
The matrix described here is the inverse of the matrix presented in [11].
Other cost terms can be added to the tracking controller. For instance, we can request a force trajectory for a selected set of force effectors. As with motions, the forces of these effectors at the end of step are linearly related to via
where is the force with zero excitations and is the force excitation response matrix. The force trajectory can then be tracking by minimizing
(10.11) |
where is a diagonal weighting matrix. When tracking force targets, the target force trajectory is fed directly into the excitation generator (Figure 10.2); there is no equivalent of the motion tracker.
To balance the effect of different cost terms, each is associated with a weight (e.g., , , and for the motion, force and regularization terms), so that the QP takes a form such as
(10.12) |
In some cases, the system forces are not linear with respect to the excitations (i.e., equation (10.6) is not valid). One situation where this occurs is when equilibrium muscle models are used (Section 4.5.3).
When forces are not linear in , the force relationship can be linearized and the quadratic program can be reformulated in terms of excitation changes ; for example,
(10.13) |
where denotes the excitation values at the beginning of the time step. Excitations are then updated according to
Incremental computation can be enabled, at the expense of slower computation, by setting the tracking controller property computeIncrementally to true (Section 10.3.2).
Applications will generally set up a tracking controller using the following steps inside the application’s build() method:
Create an instance of TrackingController and add it to the application using the root model’s addController() method.
Configure the controller with the available excitation components using its addExciter(ex) and addExciter(weight,ex) methods.
Configure the controller to track position or force trajectories for selected source components, which may include Point, Frame, or ForceEffector components. These are specified to the controller using methods such as addPointTarget(point), addFrameTarget(frame), and addForceEffectorTarget(forceEffector). These methods return a target object, which is allocated and maintained by the controller, and which is used for specifying the target positions or forces.
Create input probes or controllers to specify the desired trajectories for the sources added in Step 3. Trajectories are specified by using probes or controllers to set the appropriate target properties of the target objects returned by the addXXXTarget() methods. Likewise, output probes or monitors can be created to record the excitations and the tracked positions or forces of the sources.
Add other cost terms as needed. These may include an L2 regularization term, added using the controller’s addL2RegularizationTerm() method. Regularization attempts to minimize the norm of the excitation values and is needed to resolve redundancy if the excitation set has more degrees of freedom that the target space.
Set controller configuration parameters.
An simple example illustrating the steps of Section 10.1.6 is given by
artisynth.demos.tutorial.InverseParticle
which uses a set of point-to-point muscles to drive the position of a single particle. The model consists of a single dynamic particle, initially located at the origin, attached to a set of 16 point-to-point muscles arranged radially around it. A tracking controller is created to control the excitations of these muscles so as to enable the particle to follow a trajectory specified by a probe. The code for the model, without the package and include directives, is listed below:
The build() method begins by creating a MechModel in the usual way and disabling gravity (lines 35-37). The particle to be controlled is then created and placed at the origin and with a damping factor of 0.1 (lines 40-42). Next, the particle is attached to a set of point-to-point muscles that are arranged radially around it (lines 45-48). These are created using the method createMuscle() (lines 13-31), which attaches each to a fixed non-dynamic particle located a distance of dist from the origin, and sets its material to a SimpleAxialMuscle (Section 4.5.1.1) with a passive stiffness, damping and maximum force (at excitation 1) defined by muscleStiffness, muscleDamping, and muscleFmax (lines 4-6). Each muscle’s excitationColor and maxColoredExcitation property is set so that its color transitions to red as its excitation value varies from to (lines 28-29, Section 10.2.1.1)), and its rest length is initialized to its current length (line 30).
After the model components have been built, a TrackingController is created and added to the root model’s controller set (lines 51-52). All of the muscles are then added to the controller as exciters to be controlled (lines 54-58). The particle is added to the controller as a motion source using the addPointTarget() method (line 60), which returns a target object in the form of a TargetPoint. Since the number of exciters exceeds the degrees-of-freedom of the target space, an L2-regularization term is added (line 63) to resolve the redundancy.
Rendering properties are set at lines 67-69: points in the model are rendered as gray spheres, except for the dynamic particle which is colored white, and the muscles are drawn as blue cylinders. (The target particle, which is contained within the controller, is displayed using its default color of cyan.)
Probes are created to specify the target trajectory and record the computed excitations. The first, targetprobe, is an input probe attached to the position property of the target component, running from time 0 to 1 (lines 72-79). Its data is specified in code using addData(), which specifies three knot points with a time step of 0.5. Interpolation is set to cubic (line 78) for greater smoothness. The second, exprobe, is a probe that records all the excitations computed by the controller (lines 82-85). It is created by the utility method createOutputProbe() supplied by the InverseManager class (Section 10.4.1).
To run this example, select All demos > tutorial > InverseParticle from the Models menu. The model should load and initially appear as in Figure 10.3 (left). When run, the controller computes excitations to move the particle along the trajectory specified by the probe, which is upward and to the right (Figure 10.3, right) and back again.
The target point created by the controller is rendered separately from the source point as a cyan-colored sphere (see Section 10.2.2.2). As the simulation proceeds and certain muscles are excited, their coloring changes from blue to red, in proportion to their excitation, in accordance with the properties excitationColor and maxColoredExcitation. Recorded data for both the trajectory and the computed excitation probes are shown in Figure 10.4.